/[escript]/trunk/finley/src/CPPAdapter/MeshAdapter.h
ViewVC logotype

Contents of /trunk/finley/src/CPPAdapter/MeshAdapter.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3991 - (show annotations)
Tue Sep 25 23:52:22 2012 UTC (7 years ago) by caltinay
File MIME type: text/plain
File size: 21920 byte(s)
Updated doxygen cfg file and made a few first changes to doco.

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2012 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development since 2012 by School of Earth Sciences
13 *
14 *****************************************************************************/
15
16
17 #if !defined finley_MeshAdapter_20040526_H
18 #define finley_MeshAdapter_20040526_H
19 #include "system_dep.h"
20
21 extern "C" {
22 #include "finley/Mesh.h"
23 #include "finley/Finley.h"
24 #include "finley/Assemble.h"
25 #include "esysUtils/Esys_MPI.h"
26 }
27
28 #include "FinleyError.h"
29 #include "FinleyAdapterException.h"
30
31 #include <pasowrap/SystemMatrixAdapter.h>
32 #include <pasowrap/TransportProblemAdapter.h>
33 #include "escript/AbstractContinuousDomain.h"
34 #include "escript/FunctionSpace.h"
35 #include "escript/FunctionSpaceFactory.h"
36
37 #include <boost/shared_ptr.hpp>
38 #include <boost/python/dict.hpp>
39 #include <boost/python/extract.hpp>
40
41 #include <map>
42 #include <vector>
43 #include <string>
44 #include <sstream>
45
46 namespace finley {
47
48 // These are friends implemented in MeshAdapterFactory.cpp
49 // They are only fwd declared here so that vis.studio will accept the freind decls
50 FINLEY_DLL_API
51 escript::Domain_ptr brick(int n0,int n1,int n2,int order,
52 double l0,double l1,double l2,
53 int periodic0,int periodic1,
54 int periodic2,
55 int integrationOrder,
56 int reducedIntegrationOrder,
57 int useElementsOnFace,
58 int useFullElementOrder,
59 int optimize,
60 const std::vector<double>& points,
61 const std::vector<int>& tags,
62 const std::map<std::string, int>& tagnamestonums
63 );
64
65 FINLEY_DLL_API
66 escript::Domain_ptr rectangle(int n0,int n1,int order,
67 double l0, double l1,
68 int periodic0,int periodic1,
69 int integrationOrder,
70 int reducedIntegrationOrder,
71 int useElementsOnFace,
72 int useFullElementOrder,
73 int optimize,
74 const std::vector<double>& points,
75 const std::vector<int>& tags,
76 const std::map<std::string, int>& tagnamestonums);
77
78 }
79
80
81
82 namespace finley {
83
84 const boost::python::list EmptyPythonList = boost::python::list();
85 struct null_deleter
86 {
87 void operator()(void const *ptr) const
88 {
89 }
90 };
91
92
93 /**
94 \brief implements the AbstractContinuousDomain interface for the Finley
95 library.
96 */
97
98 class MeshAdapter : public escript::AbstractContinuousDomain {
99
100 public:
101
102 //
103 // Codes for function space types supported
104 static const int DegreesOfFreedom;
105 static const int ReducedDegreesOfFreedom;
106 static const int Nodes;
107 static const int ReducedNodes;
108 static const int Elements;
109 static const int ReducedElements;
110 static const int FaceElements;
111 static const int ReducedFaceElements;
112 static const int Points;
113 static const int ContactElementsZero;
114 static const int ReducedContactElementsZero;
115 static const int ContactElementsOne;
116 static const int ReducedContactElementsOne;
117
118 /**
119 \brief
120 Constructor for MeshAdapter
121
122 Description:
123 Constructor for MeshAdapter. The pointer passed to MeshAdapter
124 is deleted using a call to Finley_Mesh_free in the
125 MeshAdapter destructor.
126
127 Throws:
128 May throw an exception derived from EsysException
129
130 \param finleyMesh Input - A pointer to the externally constructed
131 finley mesh.The pointer passed to MeshAdapter
132 is deleted using a call to
133 Finley_Mesh_free in the MeshAdapter
134 destructor.
135 */
136 FINLEY_DLL_API
137 MeshAdapter(Finley_Mesh* finleyMesh=0);
138
139 /**
140 \brief
141 Copy constructor.
142 */
143 FINLEY_DLL_API
144 MeshAdapter(const MeshAdapter& in);
145
146 /**
147 \brief
148 Destructor for MeshAdapter. As specified in the constructor
149 this calls Finley_Mesh_free for the pointer given to the
150 constructor.
151 */
152 FINLEY_DLL_API
153 ~MeshAdapter();
154
155 /**
156 \brief
157 return the number of processors used for this domain
158 */
159 FINLEY_DLL_API
160 virtual int getMPISize() const;
161 /**
162 \brief
163 return the number MPI rank of this processor
164 */
165
166 FINLEY_DLL_API
167 virtual int getMPIRank() const;
168
169 /**
170 \brief
171 If compiled for MPI then execute an MPI_Barrier, else do nothing
172 */
173
174 FINLEY_DLL_API
175 virtual void MPIBarrier() const;
176
177 /**
178 \brief
179 Return true if on MPI processor 0, else false
180 */
181
182 FINLEY_DLL_API
183 virtual bool onMasterProcessor() const;
184
185 FINLEY_DLL_API
186 #ifdef ESYS_MPI
187 MPI_Comm
188 #else
189 unsigned int
190 #endif
191 getMPIComm() const;
192
193 /**
194 \brief
195 Write the current mesh to a file with the given name.
196 \param fileName Input - The name of the file to write to.
197 */
198 FINLEY_DLL_API
199 void write(const std::string& fileName) const;
200
201 /**
202 \brief
203 \param full
204 */
205 FINLEY_DLL_API
206 void Print_Mesh_Info(const bool full=false) const;
207
208 /**
209 \brief
210 dumps the mesh to a file with the given name.
211 \param fileName Input - The name of the file
212 */
213 FINLEY_DLL_API
214 void dump(const std::string& fileName) const;
215
216 /**
217 \brief
218 return the pointer to the underlying finley mesh structure
219 */
220 FINLEY_DLL_API
221 Finley_Mesh* getFinley_Mesh() const;
222
223 /**
224 \brief
225 Return the tag key for the given sample number.
226 \param functionSpaceType Input - The function space type.
227 \param sampleNo Input - The sample number.
228 */
229 FINLEY_DLL_API
230 int getTagFromSampleNo(int functionSpaceType, int sampleNo) const;
231
232 /**
233 \brief
234 Return the reference number of the given sample number.
235 \param functionSpaceType Input - The function space type.
236 */
237 FINLEY_DLL_API
238 const int* borrowSampleReferenceIDs(int functionSpaceType) const;
239
240 /**
241 \brief
242 Returns true if the given integer is a valid function space type
243 for this domain.
244 */
245 FINLEY_DLL_API
246 virtual bool isValidFunctionSpaceType(int functionSpaceType) const;
247
248 /**
249 \brief
250 Return a description for this domain
251 */
252 FINLEY_DLL_API
253 virtual std::string getDescription() const;
254
255 /**
256 \brief
257 Return a description for the given function space type code
258 */
259 FINLEY_DLL_API
260 virtual std::string functionSpaceTypeAsString(int functionSpaceType) const;
261
262 /**
263 \brief
264 Build the table of function space type names
265 */
266 FINLEY_DLL_API
267 void setFunctionSpaceTypeNames();
268
269 /**
270 \brief
271 Return a continuous FunctionSpace code
272 */
273 FINLEY_DLL_API
274 virtual int getContinuousFunctionCode() const;
275
276 /**
277 \brief
278 Return a continuous on reduced order nodes FunctionSpace code
279 */
280 FINLEY_DLL_API
281 virtual int getReducedContinuousFunctionCode() const;
282
283 /**
284 \brief
285 Return a function FunctionSpace code
286 */
287 FINLEY_DLL_API
288 virtual int getFunctionCode() const;
289
290 /**
291 \brief
292 Return a function with reduced integration order FunctionSpace code
293 */
294 FINLEY_DLL_API
295 virtual int getReducedFunctionCode() const;
296
297 /**
298 \brief
299 Return a function on boundary FunctionSpace code
300 */
301 FINLEY_DLL_API
302 virtual int getFunctionOnBoundaryCode() const;
303
304 /**
305 \brief
306 Return a function on boundary with reduced integration order FunctionSpace code
307 */
308 FINLEY_DLL_API
309 virtual int getReducedFunctionOnBoundaryCode() const;
310
311 /**
312 \brief
313 Return a FunctionOnContactZero code
314 */
315 FINLEY_DLL_API
316 virtual int getFunctionOnContactZeroCode() const;
317
318 /**
319 \brief
320 Return a FunctionOnContactZero code with reduced integration order
321 */
322 FINLEY_DLL_API
323 virtual int getReducedFunctionOnContactZeroCode() const;
324
325 /**
326 \brief
327 Return a FunctionOnContactOne code
328 */
329 FINLEY_DLL_API
330 virtual int getFunctionOnContactOneCode() const;
331
332 /**
333 \brief
334 Return a FunctionOnContactOne code with reduced integration order
335 */
336 FINLEY_DLL_API
337 virtual int getReducedFunctionOnContactOneCode() const;
338
339 /**
340 \brief
341 Return a Solution code
342 */
343 FINLEY_DLL_API
344 virtual int getSolutionCode() const;
345
346 /**
347 \brief
348 Return a ReducedSolution code
349 */
350 FINLEY_DLL_API
351 virtual int getReducedSolutionCode() const;
352
353 /**
354 \brief
355 Return a DiracDeltaFunctions code
356 */
357 FINLEY_DLL_API
358 virtual int getDiracDeltaFunctionsCode() const;
359
360 /**
361 5B
362 \brief
363 */
364 typedef std::map<int, std::string> FunctionSpaceNamesMapType;
365
366 /**
367 \brief
368 */
369 FINLEY_DLL_API
370 virtual int getDim() const;
371
372 /**
373 \brief
374 Returns a status indicator of the domain. The status identifier should be unique over
375 the live time if the object but may be updated if changes to the domain happen, e.g.
376 modifications to its geometry.
377
378 This has to be implemented by the actual Domain adapter.
379 */
380 FINLEY_DLL_API
381 virtual StatusType getStatus() const;
382
383
384 /**
385 \brief
386 Return the number of data points summed across all MPI processes
387 */
388 FINLEY_DLL_API
389 virtual int getNumDataPointsGlobal() const;
390
391 /**
392 \brief
393 Return the number of data points per sample, and the number of samples as a pair.
394 \param functionSpaceCode Input -
395 */
396 FINLEY_DLL_API
397 virtual std::pair<int,int> getDataShape(int functionSpaceCode) const;
398
399 /**
400 \brief
401 copies the location of data points into arg. The domain of arg has to match this.
402 has to be implemented by the actual Domain adapter.
403 */
404 FINLEY_DLL_API
405 virtual void setToX(escript::Data& arg) const;
406
407 /**
408 \brief
409 sets a map from a clear tag name to a tag key
410 \param name Input - tag name.
411 \param tag Input - tag key.
412 */
413 FINLEY_DLL_API
414 virtual void setTagMap(const std::string& name, int tag);
415
416 /**
417 \brief
418 Return the tag key for tag name.
419 \param name Input - tag name
420 */
421 FINLEY_DLL_API
422 virtual int getTag(const std::string& name) const;
423
424 /**
425 \brief
426 Returns true if name is a defined tage name.
427 \param name Input - tag name to be checked.
428 */
429 FINLEY_DLL_API
430 virtual bool isValidTagName(const std::string& name) const;
431
432 /**
433 \brief
434 Returns all tag names in a single string sperated by commas
435 */
436 FINLEY_DLL_API
437 virtual std::string showTagNames() const;
438
439 /**
440 \brief
441 assigns new location to the domain
442 */
443 FINLEY_DLL_API
444 virtual void setNewX(const escript::Data& arg);
445
446 /**
447 \brief
448 interpolates data given on source onto target where source and target have to be given on the same domain.
449 */
450 FINLEY_DLL_API
451 virtual void interpolateOnDomain(escript::Data& target,const escript::Data& source) const;
452
453
454 FINLEY_DLL_API
455 virtual bool probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const;
456
457 /**
458 \brief given a vector of FunctionSpace typecodes, pass back a code which then can all be interpolated to.
459 \return true is result is valid, false if not
460 */
461 FINLEY_DLL_API
462 bool
463 commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const;
464
465 /**
466 \brief
467 interpolates data given on source onto target where source and target are given on different domains.
468 has to be implemented by the actual Domain adapter.
469 */
470 FINLEY_DLL_API
471 virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
472
473 /**
474 \brief determines whether interpolation from source to target is possible.
475 Must be implemented by the actual Domain adapter
476 */
477 FINLEY_DLL_API
478 virtual bool probeInterpolationACross(int functionSpaceType_source,const escript::AbstractDomain& targetDomain, int functionSpaceType_target) const;
479
480 /**
481 \brief
482 copies the surface normals at data points into out. The actual function space to be considered
483 is defined by out. out has to be defined on this.
484 */
485 FINLEY_DLL_API
486 virtual void setToNormal(escript::Data& out) const;
487
488 /**
489 \brief
490 copies the size of samples into out. The actual function space to be considered
491 is defined by out. out has to be defined on this.
492 */
493 FINLEY_DLL_API
494 virtual void setToSize(escript::Data& out) const;
495
496 /**
497 \brief
498 copies the gradient of arg into grad. The actual function space to be considered
499 for the gradient is defined by grad. arg and grad have to be defined on this.
500 */
501 FINLEY_DLL_API
502 virtual void setToGradient(escript::Data& grad,const escript::Data& arg) const;
503
504 /**
505 \brief
506 copies the integrals of the function defined by arg into integrals.
507 arg has to be defined on this.
508 */
509 FINLEY_DLL_API
510 virtual void setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const;
511
512 /**
513 \brief
514 return the identifier of the matrix type to be used for the global stiffness matrix when a particular solver, package, perconditioner,
515 and symmetric matrix is used.
516 \param solver
517 \param preconditioner
518 \param package
519 \param symmetry
520 */
521 FINLEY_DLL_API
522 virtual int getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
523
524 /**
525 \brief
526 return the identifier of the transport problem type to be used when a particular solver, perconditioner, package
527 and symmetric matrix is used.
528 \param solver
529 \param preconditioner
530 \param package
531 \param symmetry
532 */
533 FINLEY_DLL_API
534 virtual int getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
535
536 /**
537 \brief
538 returns true if data on this domain and a function space of type functionSpaceCode has to
539 considered as cell centered data.
540 */
541 FINLEY_DLL_API
542 virtual bool isCellOriented(int functionSpaceCode) const;
543
544 /**
545 \brief
546 Saves a dictonary of Data objects to an OpenDX input file. The keywords are used as identifier
547
548 This has to be implemented by the actual Domain adapter.
549 */
550 FINLEY_DLL_API
551 virtual void saveDX(const std::string& filename,const boost::python::dict& arg) const;
552
553
554 /**
555 \brief
556 Saves a dictonary of Data objects to an VTK XML input file. The keywords are used as identifier
557
558 This has to be implemented by the actual Domain adapter.
559 */
560 FINLEY_DLL_API
561 virtual void saveVTK(const std::string& filename,const boost::python::dict& arg, const std::string& metadata, const std::string& metadata_schema) const;
562
563 FINLEY_DLL_API
564 virtual bool ownSample(int fs_code, index_t id) const;
565
566 /**
567 \brief
568 returns the function space representation of the type functionSpaceCode on this domain
569 as a vtkObject.
570 */
571 // vtkObject createVtkObject(int functionSpaceCode) const;
572
573 /**
574 \brief
575 adds a PDE onto the stiffness matrix mat and a rhs
576 */
577 FINLEY_DLL_API
578 virtual void addPDEToSystem(
579 escript::AbstractSystemMatrix& mat, escript::Data& rhs,
580 const escript::Data& A, const escript::Data& B, const escript::Data& C,
581 const escript::Data& D, const escript::Data& X, const escript::Data& Y,
582 const escript::Data& d, const escript::Data& y,
583 const escript::Data& d_contact, const escript::Data& y_contact,
584 const escript::Data& d_dirac, const escript::Data& y_dirac) const;
585 /**
586 \brief
587 adds a PDE onto the lumped stiffness matrix matrix
588 */
589 FINLEY_DLL_API
590 virtual void addPDEToLumpedSystem(
591 escript::Data& mat,
592 const escript::Data& D,
593 const escript::Data& d,
594 const escript::Data& d_dirac,
595 const bool useHRZ) const;
596
597 /**
598 \brief
599 adds a PDE onto the stiffness matrix mat and a rhs
600 */
601 FINLEY_DLL_API
602 virtual void addPDEToRHS(escript::Data& rhs,
603 const escript::Data& X, const escript::Data& Y,
604 const escript::Data& y, const escript::Data& y_contact, const escript::Data& y_dirac) const;
605 /**
606 \brief
607 adds a PDE onto a transport problem
608 */
609
610 FINLEY_DLL_API
611 virtual void addPDEToTransportProblem(
612 escript::AbstractTransportProblem& tp, escript::Data& source,
613 const escript::Data& M,
614 const escript::Data& A, const escript::Data& B, const escript::Data& C,const escript::Data& D,
615 const escript::Data& X,const escript::Data& Y,
616 const escript::Data& d, const escript::Data& y,
617 const escript::Data& d_contact,const escript::Data& y_contact, const escript::Data& d_dirac,const escript::Data& y_dirac) const;
618
619
620 /**
621 \brief
622 creates a SystemMatrixAdapter stiffness matrix and initializes it with zeros:
623 */
624 FINLEY_DLL_API
625 escript::ASM_ptr newSystemMatrix(
626 const int row_blocksize,
627 const escript::FunctionSpace& row_functionspace,
628 const int column_blocksize,
629 const escript::FunctionSpace& column_functionspace,
630 const int type) const;
631 /**
632 \brief
633 creates a TransportProblemAdapter
634
635 */
636
637 FINLEY_DLL_API
638 escript::ATP_ptr newTransportProblem(
639 const int blocksize,
640 const escript::FunctionSpace& functionspace,
641 const int type) const;
642
643 /**
644 \brief returns locations in the FEM nodes
645 */
646 FINLEY_DLL_API
647 virtual escript::Data getX() const;
648
649 /**
650 \brief return boundary normals at the quadrature point on the face elements
651 */
652 FINLEY_DLL_API
653 virtual escript::Data getNormal() const;
654
655 /**
656 \brief returns the element size
657 */
658 FINLEY_DLL_API
659 virtual escript::Data getSize() const;
660
661 /**
662 \brief comparison operators
663 */
664 FINLEY_DLL_API
665 virtual bool operator==(const escript::AbstractDomain& other) const;
666 FINLEY_DLL_API
667 virtual bool operator!=(const escript::AbstractDomain& other) const;
668
669 /**
670 \brief assigns new tag newTag to all samples of functionspace with a positive
671 value of mask for any its sample point.
672
673 */
674 FINLEY_DLL_API
675 virtual void setTags(const int functionSpaceType, const int newTag, const escript::Data& mask) const;
676
677 /**
678 \brief
679 return the number of tags in use and a pointer to an array with the number of tags in use
680 */
681 FINLEY_DLL_API
682 virtual int getNumberOfTagsInUse(int functionSpaceCode) const;
683
684 FINLEY_DLL_API
685 virtual const int* borrowListOfTagsInUse(int functionSpaceCode) const;
686
687
688 /**
689 \brief Checks if this domain allows tags for the specified functionSpaceCode.
690 */
691 FINLEY_DLL_API
692 virtual
693 bool canTag(int functionSpaceCode) const;
694
695 /**
696 \brief returns the approximation order used for a function space functionSpaceCode
697 */
698
699 FINLEY_DLL_API
700 virtual
701 int getApproximationOrder(const int functionSpaceCode) const;
702
703 FINLEY_DLL_API
704 bool supportsContactElements() const;
705
706
707 private:
708
709 /**
710 \brief adds points to support more Dirac delta function.
711
712 Do NOT call these at any time other than construction!
713 Using them later creates consistancy problems
714 */
715 FINLEY_DLL_API
716 void addDiracPoints( const std::vector<double>& points, const std::vector<int>& tags) const;
717 // FINLEY_DLL_API
718 // void addDiracPoint( const boost::python::list& points, const int tag=-1) const;
719 // FINLEY_DLL_API
720 // void addDiracPointWithTagName( const boost::python::list& points, const std::string& tag) const;
721
722 protected:
723
724 private:
725 void extractArgsFromDict(const boost::python::dict& arg, int& numData,
726 char**& names, escriptDataC*& data,
727 escriptDataC**& dataPtr) const;
728
729 //
730 // pointer to the externally created finley mesh
731 boost::shared_ptr<Finley_Mesh> m_finleyMesh;
732
733 // This is only provided so that the friends below can add tags during construction
734 // do not use for any other purpose
735 boost::shared_ptr<Finley_Mesh> getMesh()
736 {
737 return m_finleyMesh;
738 }
739
740 static FunctionSpaceNamesMapType m_functionSpaceTypeNames;
741
742 friend escript::Domain_ptr finley::brick(int n0,int n1,int n2,int order,
743 double l0,double l1,double l2,
744 int periodic0,int periodic1,
745 int periodic2,
746 int integrationOrder,
747 int reducedIntegrationOrder,
748 int useElementsOnFace,
749 int useFullElementOrder,
750 int optimize,
751 const std::vector<double>& points,
752 const std::vector<int>& tags,
753 const std::map<std::string, int>& tagnamestonums
754 );
755
756
757 friend escript::Domain_ptr finley::rectangle(int n0,int n1,int order,
758 double l0, double l1,
759 int periodic0,int periodic1,
760 int integrationOrder,
761 int reducedIntegrationOrder,
762 int useElementsOnFace,
763 int useFullElementOrder,
764 int optimize,
765 const std::vector<double>& points,
766 const std::vector<int>& tags,
767 const std::map<std::string, int>& tagnamestonums);
768 };
769
770 // Do not use this class. It is a convenience wrapper for the dataexporter.
771 class FINLEY_DLL_API ReferenceElementSetWrapper {
772 public:
773 ReferenceElementSetWrapper(Finley_ElementTypeId id, index_t order,
774 index_t reducedOrder);
775 ~ReferenceElementSetWrapper();
776
777 Finley_ReferenceElementSet* getElementSet() const { return m_refSet; }
778
779 private:
780 Finley_ReferenceElementSet* m_refSet;
781 };
782
783
784 } // end of namespace
785
786 #endif

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26