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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26