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

Contents of /branches/dirac/finley/src/CPPAdapter/MeshAdapter.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3555 - (show annotations)
Thu Aug 25 08:54:10 2011 UTC (7 years, 9 months ago) by jfenwick
File MIME type: text/plain
File size: 20584 byte(s)
Committing for home transfer

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26