escript  Revision_
finley/src/FinleyDomain.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2016 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16 
17 #ifndef __FINLEY_DOMAIN_H__
18 #define __FINLEY_DOMAIN_H__
19 
20 /****************************************************************************
21 
22  Finley: Domain
23 
24  A mesh is built from nodes and elements which describe the domain, surface,
25  and point sources (the latter are needed to establish links with other
26  codes, in particular particle codes). The nodes are stored in a NodeFile
27  and elements in ElementFiles. Finley domains have four ElementFiles
28  containing the elements, surface, contact and point sources, respectively.
29  Notice that the surface elements do not necessarily cover the entire
30  surface of the domain.
31 
32  The element type is fixed by the reference element, see ReferenceElement.h.
33  The numbering of the nodes starts with 0.
34 
35  Important: it is assumed that every node appears in at least one element or
36  surface element and that any node used in an element, surface element or as
37  a point is specified in the NodeFile, see also resolveNodeIds.
38 
39  In some cases it is useful to refer to a mesh entirely built from
40  order 1 (=linear) elements. The linear version of the mesh can be
41  accessed by referring to the first few nodes of each element
42  (thanks to the way the nodes are ordered). As the numbering of
43  these nodes is not continuous a relabeling vector is introduced
44  in the NodeFile. This feature is not fully implemented yet.
45 
46  All nodes and elements are tagged. The tag allows to group nodes and
47  elements. A typical application is to mark surface elements on a
48  certain portion of the domain with the same tag. All these surface
49  elements can then be assigned the same value e.g. for the pressure.
50 
51  The spatial dimensionality is determined by the type of elements
52  used and can be queried using getDim(). Notice that the element type
53  also determines the type of surface elements to be used.
54 
55 *****************************************************************************/
56 
57 #include <finley/Finley.h>
58 #include <finley/ElementFile.h>
59 #include <finley/NodeFile.h>
60 #include <finley/Util.h>
61 
62 #include <escript/AbstractContinuousDomain.h>
63 #include <escript/FunctionSpace.h>
64 #include <escript/FunctionSpaceFactory.h>
65 
66 #ifdef ESYS_HAVE_PASO
67 #include <paso/SystemMatrixPattern.h>
68 #endif
69 #ifdef ESYS_HAVE_TRILINOS
70 #include <trilinoswrap/types.h>
71 #endif
72 
73 #include <map>
74 #include <string>
75 #include <vector>
76 
77 namespace finley {
78 
79 typedef std::map<std::string, int> TagMap;
80 
82  SMT_PASO = 1<<8,
83  SMT_TRILINOS = 1<<10,
84  SMT_COMPLEX = 1<<16,
85  SMT_UNROLL = 1<<17
86 };
87 
94 {
95 public:
101  static escript::Domain_ptr load(const std::string& filename);
102 
115  static escript::Domain_ptr read(escript::JMPI mpiInfo,
116  const std::string& fileName,
117  int integrationOrder = -1,
118  int reducedIntegrationOrder = -1,
119  bool optimize = false);
120 
135  const std::string& filename,
136  int numDim, int integrationOrder = -1,
137  int reducedIntegrationOrder = -1,
138  bool optimize = false,
139  bool useMacroElements = false);
140 
158  static escript::Domain_ptr createRec4(dim_t NE0, dim_t NE1,
159  double L0, double L1,
160  bool periodic0, bool periodic1, int order,
161  int reducedOrder, bool useElementsOnFace,
162  bool optimize, escript::JMPI jmpi);
163 
184  static escript::Domain_ptr createRec8(dim_t NE0, dim_t NE1,
185  double l0, double l1,
186  bool periodic0, bool periodic1, int order,
187  int reducedOrder, bool useElementsOnFace,
188  bool useFullElementOrder,
189  bool useMacroElements, bool optimize,
190  escript::JMPI jmpi);
191 
212  static escript::Domain_ptr createHex8(dim_t NE0, dim_t NE1, dim_t NE2,
213  double l0, double l1, double l2,
214  bool periodic0, bool periodic1, bool periodic2,
215  int order, int reducedOrder,
216  bool useElementsOnFace,
217  bool optimize, escript::JMPI jmpi);
218 
241  static escript::Domain_ptr createHex20(dim_t NE0, dim_t NE1, dim_t NE2,
242  double l0, double l1, double l2,
243  bool periodic0, bool periodic1, bool periodic2,
244  int order, int reducedOrder,
245  bool useElementsOnFace,
246  bool useFullElementOrder,
247  bool useMacroElements, bool optimize,
248  escript::JMPI jmpi);
249 
258  FinleyDomain(const std::string& name, int numDim, escript::JMPI jmpi);
259 
264  FinleyDomain(const FinleyDomain& in);
265 
270  ~FinleyDomain();
271 
277  void addDiracPoints(const std::vector<double>& points,
278  const std::vector<int>& tags);
279 
284  NodeFile* getNodes() const { return m_nodes; }
285 
290  void setElements(ElementFile* elements);
291 
296  ElementFile* getElements() const { return m_elements; }
297 
302  void setFaceElements(ElementFile* elements);
303 
309 
314  void setContactElements(ElementFile* elements);
315 
321 
326  void setPoints(ElementFile* elements);
327 
332  ElementFile* getPoints() const { return m_points; }
333 
338  virtual escript::JMPI getMPI() const { return m_mpiInfo; }
339 
344  virtual int getMPISize() const { return m_mpiInfo->size; }
345 
350  virtual int getMPIRank() const { return m_mpiInfo->rank; }
351 
356  virtual void MPIBarrier() const;
357 
362  virtual bool onMasterProcessor() const { return getMPIRank() == 0; }
363 
364  MPI_Comm getMPIComm() const { return m_mpiInfo->comm; }
365 
372  void write(const std::string& fileName) const;
373 
378  void Print_Mesh_Info(bool full=false) const;
379 
385  void dump(const std::string& fileName) const;
386 
393  int getTagFromSampleNo(int functionSpaceType, index_t sampleNo) const;
394 
400  const index_t* borrowSampleReferenceIDs(int functionSpaceType) const;
401 
407  virtual bool isValidFunctionSpaceType(int functionSpaceType) const;
408 
413  virtual std::string getDescription() const;
414 
419  virtual std::string functionSpaceTypeAsString(int functionSpaceType) const;
420 
426 
431  virtual int getContinuousFunctionCode() const;
432 
437  virtual int getReducedContinuousFunctionCode() const;
438 
443  virtual int getFunctionCode() const;
444 
449  virtual int getReducedFunctionCode() const;
450 
455  virtual int getFunctionOnBoundaryCode() const;
456 
461  virtual int getReducedFunctionOnBoundaryCode() const;
462 
467  virtual int getFunctionOnContactZeroCode() const;
468 
473  virtual int getReducedFunctionOnContactZeroCode() const;
474 
479  virtual int getFunctionOnContactOneCode() const;
480 
485  virtual int getReducedFunctionOnContactOneCode() const;
486 
491  virtual int getSolutionCode() const;
492 
497  virtual int getReducedSolutionCode() const;
498 
503  virtual int getDiracDeltaFunctionsCode() const;
504 
508  typedef std::map<int, std::string> FunctionSpaceNamesMapType;
509 
513  virtual int getDim() const { return m_nodes->numDim; }
514 
521  virtual StatusType getStatus() const;
522 
527  virtual dim_t getNumDataPointsGlobal() const;
528 
534  virtual std::pair<int,dim_t> getDataShape(int functionSpaceCode) const;
535 
541  virtual void setToX(escript::Data& arg) const;
542 
549  virtual void setTagMap(const std::string& name, int tag);
550 
556  virtual int getTag(const std::string& name) const;
557 
563  virtual bool isValidTagName(const std::string& name) const;
564 
569  virtual std::string showTagNames() const;
570 
575  virtual void setNewX(const escript::Data& arg);
576 
581  virtual void interpolateOnDomain(escript::Data& target,
582  const escript::Data& source) const;
583 
584  virtual bool probeInterpolationOnDomain(int functionSpaceType_source,
585  int functionSpaceType_target) const;
586 
587  virtual signed char preferredInterpolationOnDomain(int functionSpaceType_source, int functionSpaceType_target) const;
588 
593  bool commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const;
594 
599  virtual void interpolateAcross(escript::Data& target, const escript::Data& source) const;
600 
604  virtual bool probeInterpolationAcross(int functionSpaceType_source,
605  const escript::AbstractDomain& targetDomain,
606  int functionSpaceType_target) const;
607 
613  virtual void setToNormal(escript::Data& out) const;
614 
620  virtual void setToSize(escript::Data& out) const;
621 
627  virtual void setToGradient(escript::Data& grad, const escript::Data& arg) const;
628 
634  virtual void setToIntegrals(std::vector<double>& integrals, const escript::Data& arg) const;
635 
644  virtual int getSystemMatrixTypeId(const boost::python::object& options) const;
645 
655  virtual int getTransportTypeId(int solver, int preconditioner, int package,
656  bool symmetry) const;
657 
663  virtual bool isCellOriented(int functionSpaceCode) const;
664 
665  virtual bool ownSample(int fsCode, index_t id) const;
666 
671  virtual void addPDEToSystem(
673  const escript::Data& A, const escript::Data& B,
674  const escript::Data& C, const escript::Data& D,
675  const escript::Data& X, const escript::Data& Y,
676  const escript::Data& d, const escript::Data& y,
677  const escript::Data& d_contact,
678  const escript::Data& y_contact,
679  const escript::Data& d_dirac,
680  const escript::Data& y_dirac) const;
681 
686  virtual void addPDEToLumpedSystem(escript::Data& mat,
687  const escript::Data& D,
688  const escript::Data& d,
689  const escript::Data& d_dirac,
690  bool useHRZ) const;
691 
696  virtual void addPDEToRHS(escript::Data& rhs, const escript::Data& X,
697  const escript::Data& Y, const escript::Data& y,
698  const escript::Data& y_contact,
699  const escript::Data& y_dirac) const;
700 
705  virtual void addPDEToTransportProblem(
707  escript::Data& source, const escript::Data& M,
708  const escript::Data& A, const escript::Data& B,
709  const escript::Data& C, const escript::Data& D,
710  const escript::Data& X, const escript::Data& Y,
711  const escript::Data& d, const escript::Data& y,
712  const escript::Data& d_contact,
713  const escript::Data& y_contact,
714  const escript::Data& d_dirac,
715  const escript::Data& y_dirac) const;
716 
722  int row_blocksize,
723  const escript::FunctionSpace& row_functionspace,
724  int column_blocksize,
725  const escript::FunctionSpace& column_functionspace,
726  int type) const;
727 
733  const escript::FunctionSpace& functionspace,
734  int type) const;
735 
739  virtual escript::Data getX() const;
740 
745  virtual escript::Data getNormal() const;
746 
750  virtual escript::Data getSize() const;
751 
755  virtual bool operator==(const escript::AbstractDomain& other) const;
756  virtual bool operator!=(const escript::AbstractDomain& other) const;
757 
762  virtual void setTags(int functionSpaceType, int newTag,
763  const escript::Data& mask) const;
764 
770  virtual int getNumberOfTagsInUse(int functionSpaceCode) const;
771 
772  virtual const int* borrowListOfTagsInUse(int functionSpaceCode) const;
773 
778  virtual bool canTag(int functionSpaceCode) const;
779 
783  virtual int getApproximationOrder(int functionSpaceCode) const;
784 
785  virtual bool supportsContactElements() const { return true; }
786 
788  const escript::FunctionSpace& what, long seed,
789  const boost::python::tuple& filter) const;
790 
795  const TagMap& getTagMap() const { return m_tagMap; }
796 
797  void createMappings(const IndexVector& dofDistribution,
798  const IndexVector& nodeDistribution);
799 
800 #ifdef ESYS_HAVE_PASO
801  paso::SystemMatrixPattern_ptr getPasoPattern(bool reducedRowOrder,
803  bool reducedColOrder) const;
804 #endif
805 
806 #ifdef ESYS_HAVE_TRILINOS
807  esys_trilinos::const_TrilinosGraph_ptr getTrilinosGraph(bool reducedOrder) const;
809 #endif
810 
811  void glueFaces(double safetyFactor, double tolerance, bool optimize);
812 
813  void joinFaces(double safetyFactor, double tolerance, bool optimize);
814 
817  static FinleyDomain* merge(const std::vector<const FinleyDomain*>& meshes);
818 
819 private:
820  void prepare(bool optimize);
821 
822  void setOrders();
823 
831  void resolveNodeIds();
832 
835  void relabelElementNodes(const IndexVector& newNode, index_t offset);
836 
837 #ifdef ESYS_HAVE_PASO
838  paso::SystemMatrixPattern_ptr makePasoPattern(bool reducedRowOrder,
839  bool reducedColOrder) const;
840 #endif
841 #ifdef ESYS_HAVE_TRILINOS
842  esys_trilinos::GraphType* createTrilinosGraph(bool reducedOrder) const;
843 #endif
844  void createColoring(const IndexVector& dofMap);
845  void distributeByRankOfDOF(const IndexVector& distribution);
846  void markNodes(std::vector<short>& mask, index_t offset, bool useLinear) const;
847  void optimizeDOFDistribution(IndexVector& distribution);
848  void optimizeDOFLabeling(const IndexVector& distribution);
850  void findMatchingFaces(double safetyFactor, double tolerance, int* numPairs,
851  int* elem0, int* elem1, int* matchingNodes) const;
852  void updateTagList();
853  void printElementInfo(const ElementFile* e, const std::string& title,
854  const std::string& defaultType, bool full) const;
855 
856  void writeElementInfo(std::ostream& stream, const ElementFile* e,
857  const std::string& defaultType) const;
858 
862  std::string m_name;
878  TagMap m_tagMap;
879 #ifdef ESYS_HAVE_PASO
880  // pointer to the sparse matrix patterns
881  mutable paso::SystemMatrixPattern_ptr FullFullPattern;
882  mutable paso::SystemMatrixPattern_ptr FullReducedPattern;
883  mutable paso::SystemMatrixPattern_ptr ReducedFullPattern;
884  mutable paso::SystemMatrixPattern_ptr ReducedReducedPattern;
885 #endif
886 #ifdef ESYS_HAVE_TRILINOS
887  mutable esys_trilinos::TrilinosGraph_ptr m_fullGraph;
888  mutable esys_trilinos::TrilinosGraph_ptr m_reducedGraph;
889 #endif
890 
891  static FunctionSpaceNamesMapType m_functionSpaceTypeNames;
892 };
893 
894 } // end of namespace
895 
896 #endif // __FINLEY_DOMAIN_H__
897 
virtual void setToSize(escript::Data &out) const
copies the size of samples into out. The actual function space to be considered is defined by out...
Definition: finley/src/FinleyDomain.cpp:1296
ElementFile * getPoints() const
returns a pointer to this domain&#39;s point (nodal) element file
Definition: finley/src/FinleyDomain.h:332
AbstractContinuousDomain, base class for continuous domains.
Definition: AbstractContinuousDomain.h:45
Definition: FunctionSpace.h:34
void writeElementInfo(std::ostream &stream, const ElementFile *e, const std::string &defaultType) const
Definition: finley/src/Mesh_write.cpp:33
void relabelElementNodes(const IndexVector &newNode, index_t offset)
Definition: finley/src/FinleyDomain.cpp:185
Definition: finley/src/FinleyDomain.h:85
escript::ATP_ptr newTransportProblem(int blocksize, const escript::FunctionSpace &functionspace, int type) const
creates a TransportProblem
Definition: finley/src/FinleyDomain.cpp:1440
virtual int getSystemMatrixTypeId(const boost::python::object &options) const
return the identifier of the matrix type to be used for the global stiffness matrix when a particular...
Definition: finley/src/FinleyDomain.cpp:1780
double l2(dim_t n, const double *x, escript::JMPI mpiinfo)
returns the global L2 norm of x
Definition: PasoUtil.cpp:501
virtual void setNewX(const escript::Data &arg)
assigns new location to the domain
Definition: finley/src/FinleyDomain.cpp:1334
virtual escript::Data getX() const
returns locations in the FEM nodes
Definition: finley/src/FinleyDomain.cpp:1852
void Print_Mesh_Info(bool full=false) const
Definition: finley/src/Mesh_write.cpp:147
Definition: finley/src/FinleyDomain.h:82
std::map< std::string, int > TagMap
Definition: finley/src/FinleyDomain.h:79
int getTagFromSampleNo(int functionSpaceType, index_t sampleNo) const
Return the tag key for the given sample number.
Definition: finley/src/FinleyDomain.cpp:1908
ElementFile * getFaceElements() const
returns a pointer to this domain&#39;s face element file
Definition: finley/src/FinleyDomain.h:308
void optimizeDOFDistribution(IndexVector &distribution)
Definition: finley/src/Mesh_optimizeDOFDistribution.cpp:65
ElementFile * m_points
the table of points (treated as elements of dimension 0)
Definition: finley/src/FinleyDomain.h:876
void glueFaces(double safetyFactor, double tolerance, bool optimize)
Definition: Mesh_glueFaces.cpp:32
ElementFile * m_faceElements
the table of face elements
Definition: finley/src/FinleyDomain.h:872
virtual void setTags(int functionSpaceType, int newTag, const escript::Data &mask) const
assigns new tag newTag to all samples of functionspace with a positive value of mask for any its samp...
Definition: finley/src/FinleyDomain.cpp:1948
int numDim
number of spatial dimensions
Definition: finley/src/NodeFile.h:161
boost::shared_ptr< AbstractDomain > Domain_ptr
Definition: AbstractDomain.h:36
virtual void interpolateOnDomain(escript::Data &target, const escript::Data &source) const
interpolates data given on source onto target where source and target have to be given on the same do...
Definition: finley/src/FinleyDomain.cpp:878
FinleyDomain implements the AbstractContinuousDomain interface for the Finley library.
Definition: finley/src/FinleyDomain.h:93
NodeFile * getNodes() const
returns a pointer to this domain&#39;s node file
Definition: finley/src/FinleyDomain.h:284
void joinFaces(double safetyFactor, double tolerance, bool optimize)
Definition: Mesh_joinFaces.cpp:32
escript::JMPI m_mpiInfo
MPI information.
Definition: finley/src/FinleyDomain.h:860
void prepare(bool optimize)
prepares the mesh for further use
Definition: finley/src/FinleyDomain.cpp:2161
int integrationOrder
Definition: finley/src/FinleyDomain.h:865
virtual std::string functionSpaceTypeAsString(int functionSpaceType) const
Return a description for the given function space type code.
Definition: finley/src/FinleyDomain.cpp:559
ElementFile * getElements() const
returns a pointer to this domain&#39;s element file
Definition: finley/src/FinleyDomain.h:296
virtual bool isValidFunctionSpaceType(int functionSpaceType) const
Returns true if the given integer is a valid function space type for this domain. ...
Definition: finley/src/FinleyDomain.cpp:570
virtual int getDiracDeltaFunctionsCode() const
Return a DiracDeltaFunctions code.
Definition: finley/src/FinleyDomain.cpp:667
std::vector< index_t > IndexVector
Definition: DataTypes.h:62
virtual void addPDEToSystem(escript::AbstractSystemMatrix &mat, escript::Data &rhs, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y, const escript::Data &d, const escript::Data &y, const escript::Data &d_contact, const escript::Data &y_contact, const escript::Data &d_dirac, const escript::Data &y_dirac) const
adds a PDE onto the stiffness matrix mat and a rhs
Definition: finley/src/FinleyDomain.cpp:766
boost::shared_ptr< JMPI_ > JMPI
Definition: EsysMPI.h:71
virtual int getFunctionCode() const
Return a function FunctionSpace code.
Definition: finley/src/FinleyDomain.cpp:617
virtual void setToNormal(escript::Data &out) const
copies the surface normals at data points into out. The actual function space to be considered is def...
Definition: finley/src/FinleyDomain.cpp:1164
static escript::Domain_ptr createHex20(dim_t NE0, dim_t NE1, dim_t NE2, double l0, double l1, double l2, bool periodic0, bool periodic1, bool periodic2, int order, int reducedOrder, bool useElementsOnFace, bool useFullElementOrder, bool useMacroElements, bool optimize, escript::JMPI jmpi)
Creates a 3-dimensional rectangular domain with second order (Hex20 or Hex27) elements.
Definition: Mesh_hex20.cpp:38
virtual std::string showTagNames() const
Returns all tag names in a single string sperated by commas.
Definition: finley/src/FinleyDomain.cpp:2006
void resolveNodeIds()
Definition: finley/src/FinleyDomain.cpp:2349
void markNodes(std::vector< short > &mask, index_t offset, bool useLinear) const
Definition: finley/src/FinleyDomain.cpp:176
std::vector< int > ShapeType
The shape of a single datapoint.
Definition: DataTypes.h:42
virtual void interpolateAcross(escript::Data &target, const escript::Data &source) const
interpolates data given on source onto target where source and target are given on different domains...
Definition: finley/src/FinleyDomain.cpp:1188
static dim_t M
Definition: SparseMatrix_saveHB.cpp:37
virtual int getSolutionCode() const
Return a Solution code.
Definition: finley/src/FinleyDomain.cpp:657
static FunctionSpaceNamesMapType m_functionSpaceTypeNames
Definition: finley/src/FinleyDomain.h:891
virtual bool probeInterpolationAcross(int functionSpaceType_source, const escript::AbstractDomain &targetDomain, int functionSpaceType_target) const
determines whether interpolation from source to target is possible.
Definition: finley/src/FinleyDomain.cpp:1756
virtual std::string getDescription() const
Return a description for this domain.
Definition: finley/src/FinleyDomain.cpp:554
A suite of factory methods for creating various finley domains.
Definition: finley/src/Assemble.h:31
virtual int getMPISize() const
returns the number of processors used for this domain
Definition: finley/src/FinleyDomain.h:344
virtual void addPDEToTransportProblem(escript::AbstractTransportProblem &tp, escript::Data &source, const escript::Data &M, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y, const escript::Data &d, const escript::Data &y, const escript::Data &d_contact, const escript::Data &y_contact, const escript::Data &d_dirac, const escript::Data &y_dirac) const
adds a PDE onto a transport problem
Definition: finley/src/FinleyDomain.cpp:837
boost::shared_ptr< SystemMatrixPattern > SystemMatrixPattern_ptr
Definition: SystemMatrixPattern.h:39
boost::shared_ptr< AbstractTransportProblem > ATP_ptr
Definition: AbstractTransportProblem.h:160
virtual int getReducedFunctionCode() const
Return a function with reduced integration order FunctionSpace code.
Definition: finley/src/FinleyDomain.cpp:622
static escript::Domain_ptr readGmsh(escript::JMPI mpiInfo, const std::string &filename, int numDim, int integrationOrder=-1, int reducedIntegrationOrder=-1, bool optimize=false, bool useMacroElements=false)
reads a gmsh mesh file.
Definition: finley/src/Mesh_readGmsh.cpp:1149
void optimizeElementOrdering()
redistributes elements to minimize communication during assemblage
Definition: finley/src/FinleyDomain.cpp:2440
bool commonFunctionSpace(const std::vector< int > &fs, int &resultcode) const
given a vector of FunctionSpace typecodes, pass back a code which then can all be interpolated to...
Definition: finley/src/FinleyDomain.cpp:1497
virtual escript::Data getSize() const
returns the element size
Definition: finley/src/FinleyDomain.cpp:1862
virtual int getFunctionOnBoundaryCode() const
Return a function on boundary FunctionSpace code.
Definition: finley/src/FinleyDomain.cpp:627
virtual bool onMasterProcessor() const
returns true if on MPI processor 0, else false
Definition: finley/src/FinleyDomain.h:362
virtual const int * borrowListOfTagsInUse(int functionSpaceCode) const
Definition: finley/src/FinleyDomain.cpp:2050
int index_t
type for array/matrix indices used both globally and on each rank
Definition: DataTypes.h:59
void distributeByRankOfDOF(const IndexVector &distribution)
Definition: finley/src/FinleyDomain.cpp:2218
void findMatchingFaces(double safetyFactor, double tolerance, int *numPairs, int *elem0, int *elem1, int *matchingNodes) const
Definition: Mesh_findMatchingFaces.cpp:68
virtual void setTagMap(const std::string &name, int tag)
sets a map from a clear tag name to a tag key
Definition: finley/src/FinleyDomain.cpp:1985
std::string m_name
domain description
Definition: finley/src/FinleyDomain.h:862
static FinleyDomain * merge(const std::vector< const FinleyDomain *> &meshes)
Definition: Mesh_merge.cpp:24
virtual bool isValidTagName(const std::string &name) const
Returns true if name is a defined tag name.
Definition: finley/src/FinleyDomain.cpp:2001
virtual bool operator==(const escript::AbstractDomain &other) const
comparison operators
Definition: finley/src/FinleyDomain.cpp:1762
virtual bool canTag(int functionSpaceCode) const
Checks if this domain allows tags for the specified functionSpace code.
Definition: finley/src/FinleyDomain.cpp:2096
virtual bool probeInterpolationOnDomain(int functionSpaceType_source, int functionSpaceType_target) const
True if interpolation is possible from source to target.
Definition: finley/src/FinleyDomain.cpp:1618
void setFaceElements(ElementFile *elements)
replaces the face element file by elements
Definition: finley/src/FinleyDomain.cpp:110
virtual int getReducedFunctionOnBoundaryCode() const
Return a function on boundary with reduced integration order FunctionSpace code.
Definition: finley/src/FinleyDomain.cpp:632
virtual int getDim() const
returns the dimensionality of this domain
Definition: finley/src/FinleyDomain.h:513
void addDiracPoints(const std::vector< double > &points, const std::vector< int > &tags)
adds Dirac delta points. Do NOT call this at any time other than construction! Using them later creat...
Definition: Mesh_addPoints.cpp:48
std::map< int, std::string > FunctionSpaceNamesMapType
Definition: finley/src/FinleyDomain.h:508
virtual escript::JMPI getMPI() const
returns a reference to the MPI information wrapper for this domain
Definition: finley/src/FinleyDomain.h:338
int StatusType
Definition: AbstractDomain.h:48
virtual int getApproximationOrder(int functionSpaceCode) const
returns the approximation order used for a function space functionSpaceCode
Definition: finley/src/FinleyDomain.cpp:2120
ElementFile * getContactElements() const
returns a pointer to this domain&#39;s contact element file
Definition: finley/src/FinleyDomain.h:320
virtual bool isCellOriented(int functionSpaceCode) const
returns true if data on this domain and a function space of type functionSpaceCode has to considered ...
Definition: finley/src/FinleyDomain.cpp:1472
int reducedApproximationOrder
Definition: finley/src/FinleyDomain.h:864
Data represents a collection of datapoints.
Definition: Data.h:63
void createMappings(const IndexVector &dofDistribution, const IndexVector &nodeDistribution)
Definition: finley/src/FinleyDomain.cpp:167
virtual std::pair< int, dim_t > getDataShape(int functionSpaceCode) const
Return the number of data points per sample, and the number of samples as a pair. ...
Definition: finley/src/FinleyDomain.cpp:684
virtual int getContinuousFunctionCode() const
Return a continuous FunctionSpace code.
Definition: finley/src/FinleyDomain.cpp:607
Definition: finley/src/NodeFile.h:40
virtual int getTag(const std::string &name) const
Return the tag key for tag name.
Definition: finley/src/FinleyDomain.cpp:1990
int reducedIntegrationOrder
Definition: finley/src/FinleyDomain.h:866
void printElementInfo(const ElementFile *e, const std::string &title, const std::string &defaultType, bool full) const
Definition: finley/src/Mesh_write.cpp:52
virtual int getReducedSolutionCode() const
Return a ReducedSolution code.
Definition: finley/src/FinleyDomain.cpp:662
const TagMap & getTagMap() const
returns a reference to the tag name->value map
Definition: finley/src/FinleyDomain.h:795
void setElements(ElementFile *elements)
replaces the element file by elements
Definition: finley/src/FinleyDomain.cpp:104
virtual int getMPIRank() const
returns the number MPI rank of this processor
Definition: finley/src/FinleyDomain.h:350
static escript::Domain_ptr createRec4(dim_t NE0, dim_t NE1, double L0, double L1, bool periodic0, bool periodic1, int order, int reducedOrder, bool useElementsOnFace, bool optimize, escript::JMPI jmpi)
Creates a 2-dimensional rectangular domain with first order (Rec4) elements in the rectangle [0...
Definition: Mesh_rec4.cpp:25
void optimizeDOFLabeling(const IndexVector &distribution)
optimizes the labeling of the DOFs on each processor
Definition: finley/src/FinleyDomain.cpp:2268
Definition: finley/src/FinleyDomain.h:84
void setFunctionSpaceTypeNames()
Build the table of function space type names.
Definition: finley/src/FinleyDomain.cpp:577
void setOrders()
Definition: finley/src/FinleyDomain.cpp:128
virtual void MPIBarrier() const
If compiled for MPI then execute an MPI_Barrier, else do nothing.
Definition: finley/src/FinleyDomain.cpp:97
virtual void addPDEToLumpedSystem(escript::Data &mat, const escript::Data &D, const escript::Data &d, const escript::Data &d_dirac, bool useHRZ) const
adds a PDE onto the lumped stiffness matrix matrix
Definition: finley/src/FinleyDomain.cpp:799
int MPI_Comm
Definition: EsysMPI.h:41
static escript::Domain_ptr load(const std::string &filename)
recovers domain from a dump file
Definition: finley/src/DomainFactory.cpp:62
MPI_Comm getMPIComm() const
get the communicator for this domain. Returns an integer on non-MPI builds Routine must be implemente...
Definition: finley/src/FinleyDomain.h:364
static escript::Domain_ptr read(escript::JMPI mpiInfo, const std::string &fileName, int integrationOrder=-1, int reducedIntegrationOrder=-1, bool optimize=false)
reads a mesh from a fly file. For MPI parallel runs fans out the mesh to multiple processes...
Definition: finley/src/Mesh_read.cpp:152
virtual bool ownSample(int fsCode, index_t id) const
Definition: finley/src/FinleyDomain.cpp:1347
virtual int getReducedFunctionOnContactZeroCode() const
Return a FunctionOnContactZero code with reduced integration order.
Definition: finley/src/FinleyDomain.cpp:642
void setContactElements(ElementFile *elements)
replaces the contact element file by elements
Definition: finley/src/FinleyDomain.cpp:116
Give a short description of what AbstractTransportProblem does.
Definition: AbstractTransportProblem.h:43
void setPoints(ElementFile *elements)
replaces the point element file by elements
Definition: finley/src/FinleyDomain.cpp:122
virtual int getTransportTypeId(int solver, int preconditioner, int package, bool symmetry) const
return the identifier of the transport problem type to be used when a particular solver, perconditioner, package and symmetric matrix is used.
Definition: finley/src/FinleyDomain.cpp:1840
virtual void setToX(escript::Data &arg) const
copies the location of data points into arg. The domain of arg has to match this. has to be implement...
Definition: finley/src/FinleyDomain.cpp:1145
virtual int getFunctionOnContactZeroCode() const
Return a FunctionOnContactZero code.
Definition: finley/src/FinleyDomain.cpp:637
virtual StatusType getStatus() const
Returns a status indicator of the domain. The status identifier should be unique over the live time i...
Definition: finley/src/FinleyDomain.cpp:2115
virtual void setToGradient(escript::Data &grad, const escript::Data &arg) const
copies the gradient of arg into grad. The actual function space to be considered for the gradient is ...
Definition: finley/src/FinleyDomain.cpp:1241
escript::ASM_ptr newSystemMatrix(int row_blocksize, const escript::FunctionSpace &row_functionspace, int column_blocksize, const escript::FunctionSpace &column_functionspace, int type) const
creates a stiffness matrix and initializes it with zeros
Definition: finley/src/FinleyDomain.cpp:1377
ElementFile * m_contactElements
the table of contact elements
Definition: finley/src/FinleyDomain.h:874
Base class for escript system matrices.
Definition: AbstractSystemMatrix.h:42
ElementFile * m_elements
the table of the elements
Definition: finley/src/FinleyDomain.h:870
static escript::Domain_ptr createRec8(dim_t NE0, dim_t NE1, double l0, double l1, bool periodic0, bool periodic1, int order, int reducedOrder, bool useElementsOnFace, bool useFullElementOrder, bool useMacroElements, bool optimize, escript::JMPI jmpi)
Creates a 2-dimensional rectangular domain with second order (Rec8 or Rec9) elements in the rectangle...
Definition: Mesh_rec8.cpp:25
TagMap m_tagMap
the tag map mapping names to tag keys
Definition: finley/src/FinleyDomain.h:878
void createColoring(const IndexVector &dofMap)
tries to reduce the number of colours for all element files
Definition: finley/src/FinleyDomain.cpp:2431
void updateTagList()
regenerates list of tags in use for node file and element files
Definition: finley/src/FinleyDomain.cpp:2449
virtual escript::Data getNormal() const
returns boundary normals at the quadrature point on the face elements
Definition: finley/src/FinleyDomain.cpp:1857
void dump(const std::string &fileName) const
dumps the mesh to a file with the given name.
Definition: finley/src/FinleyDomain.cpp:193
boost::shared_ptr< AbstractSystemMatrix > ASM_ptr
Definition: AbstractSystemMatrix.h:32
static escript::Domain_ptr createHex8(dim_t NE0, dim_t NE1, dim_t NE2, double l0, double l1, double l2, bool periodic0, bool periodic1, bool periodic2, int order, int reducedOrder, bool useElementsOnFace, bool optimize, escript::JMPI jmpi)
Creates a 3-dimensional rectangular domain with first order (Hex8) elements.
Definition: Mesh_hex8.cpp:39
FinleyDomain(const std::string &name, int numDim, escript::JMPI jmpi)
Constructor for FinleyDomain.
Definition: finley/src/FinleyDomain.cpp:55
Definition: finley/src/FinleyDomain.h:83
const index_t * borrowSampleReferenceIDs(int functionSpaceType) const
Return the reference number of the given sample number.
Definition: finley/src/FinleyDomain.cpp:1867
SystemMatrixType
Definition: finley/src/FinleyDomain.h:81
Base class for all escript domains.
Definition: AbstractDomain.h:45
virtual dim_t getNumDataPointsGlobal() const
Return the number of data points summed across all MPI processes.
Definition: finley/src/FinleyDomain.cpp:675
virtual int getReducedFunctionOnContactOneCode() const
Return a FunctionOnContactOne code with reduced integration order.
Definition: finley/src/FinleyDomain.cpp:652
void write(const std::string &fileName) const
writes the current mesh to a file with the given name in the fly file format.
Definition: finley/src/Mesh_write.cpp:87
virtual bool supportsContactElements() const
Definition: finley/src/FinleyDomain.h:785
virtual void setToIntegrals(std::vector< double > &integrals, const escript::Data &arg) const
copies the integrals of the function defined by arg into integrals. arg has to be defined on this...
Definition: finley/src/FinleyDomain.cpp:1198
virtual int getNumberOfTagsInUse(int functionSpaceCode) const
returns the number of tags in use and a pointer to an array with the number of tags in use ...
Definition: finley/src/FinleyDomain.cpp:2019
virtual int getFunctionOnContactOneCode() const
Return a FunctionOnContactOne code.
Definition: finley/src/FinleyDomain.cpp:647
virtual void addPDEToRHS(escript::Data &rhs, const escript::Data &X, const escript::Data &Y, const escript::Data &y, const escript::Data &y_contact, const escript::Data &y_dirac) const
adds a PDE onto the stiffness matrix mat and a rhs
Definition: finley/src/FinleyDomain.cpp:813
int approximationOrder
Definition: finley/src/FinleyDomain.h:863
virtual int getReducedContinuousFunctionCode() const
Return a continuous on reduced order nodes FunctionSpace code.
Definition: finley/src/FinleyDomain.cpp:612
virtual signed char preferredInterpolationOnDomain(int functionSpaceType_source, int functionSpaceType_target) const
Preferred direction of interpolation. If you really need to test for a particular direction...
Definition: finley/src/FinleyDomain.cpp:1745
NodeFile * m_nodes
the table of the nodes
Definition: finley/src/FinleyDomain.h:868
~FinleyDomain()
Destructor for FinleyDomain.
Definition: finley/src/FinleyDomain.cpp:88
index_t dim_t
Definition: DataTypes.h:64
Definition: finley/src/ElementFile.h:61
virtual bool operator!=(const escript::AbstractDomain &other) const
Return true if given domains are not equal.
Definition: finley/src/FinleyDomain.cpp:1775
virtual escript::Data randomFill(const escript::DataTypes::ShapeType &shape, const escript::FunctionSpace &what, long seed, const boost::python::tuple &filter) const
Fills the data object with filtered random values.
Definition: finley/src/FinleyDomain.cpp:2147