escript  Revision_
DudleyDomain.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2020 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-2017 by Centre for Geoscience Computing (GeoComp)
14 * Development from 2019 by School of Earth and Environmental Sciences
15 **
16 *****************************************************************************/
17 
18 #ifndef __DUDLEY_DOMAIN_H__
19 #define __DUDLEY_DOMAIN_H__
20 
21 /****************************************************************************
22 
23  Dudley: Domain
24 
25  A mesh is built from nodes and elements which describe the domain, surface,
26  and point sources (the latter are needed to establish links with other
27  codes, in particular particle codes). The nodes are stored in a NodeFile
28  and elements in ElementFiles. Dudley domains have three ElementFiles
29  containing the elements, surface and point sources, respectively.
30  Notice that the surface elements do not necessarily cover the entire
31  surface of the domain.
32 
33  The element type is either Tri3 or Tet4 depending on dimensionality
34  and also determines the type of surface elements to be used.
35 
36  The numbering of the nodes starts with 0.
37 
38  Important: it is assumed that every node appears in at least one element or
39  surface element and that any node used in an element, surface element or as
40  a point is specified in the NodeFile, see also resolveNodeIds.
41 
42  All nodes and elements are tagged. The tag allows to group nodes and
43  elements. A typical application is to mark surface elements on a
44  certain portion of the domain with the same tag. All these surface
45  elements can then be assigned the same value e.g. for the pressure.
46 
47 *****************************************************************************/
48 
49 #include "system_dep.h"
50 
51 #include <dudley/Dudley.h>
52 #include <dudley/ElementFile.h>
53 #include <dudley/NodeFile.h>
54 #include <dudley/Util.h>
55 
56 #include <escript/AbstractContinuousDomain.h>
57 #include <escript/FunctionSpace.h>
58 #include <escript/FunctionSpaceFactory.h>
59 
60 #ifdef ESYS_HAVE_PASO
61 #include <paso/SystemMatrixPattern.h>
62 #endif
63 #ifdef ESYS_HAVE_TRILINOS
64 #include <trilinoswrap/types.h>
65 #endif
66 
67 #include <map>
68 #include <string>
69 #include <vector>
70 
71 namespace dudley {
72 
73 typedef std::map<std::string, int> TagMap;
74 
76  SMT_PASO = 1<<8,
77  SMT_TRILINOS = 1<<10,
78  SMT_COMPLEX = 1<<16,
79  SMT_UNROLL = 1<<17
80 };
81 
88 {
89 public:
95  static escript::Domain_ptr load(const std::string& filename);
96 
105  static escript::Domain_ptr read(escript::JMPI mpiInfo,
106  const std::string& filename, bool optimize);
107 
117  const std::string& filename, int numDim,
118  bool optimize);
119 
131  static escript::Domain_ptr create2D(dim_t NE0, dim_t NE1, double l0,
132  double l1, bool optimize,
133  escript::JMPI jmpi);
134 
148  static escript::Domain_ptr create3D(dim_t NE0, dim_t NE1, dim_t NE2,
149  double l0, double l1, double l2,
150  bool optimize, escript::JMPI jmpi);
151 
160  DudleyDomain(const std::string& name, int numDim, escript::JMPI jmpi);
161 
166  DudleyDomain(const DudleyDomain& in);
167 
172  ~DudleyDomain();
173 
178  NodeFile* getNodes() const { return m_nodes; }
179 
184  void setElements(ElementFile* elements);
185 
190  ElementFile* getElements() const { return m_elements; }
191 
196  void setFaceElements(ElementFile* elements);
197 
202  ElementFile* getFaceElements() const { return m_faceElements; }
203 
208  void setPoints(ElementFile* elements);
209 
214  ElementFile* getPoints() const { return m_points; }
215 
220  virtual escript::JMPI getMPI() const { return m_mpiInfo; }
221 
226  virtual int getMPISize() const { return m_mpiInfo->size; }
227 
232  virtual int getMPIRank() const { return m_mpiInfo->rank; }
233 
238  virtual void MPIBarrier() const;
239 
244  virtual bool onMasterProcessor() const { return getMPIRank() == 0; }
245 
246  MPI_Comm getMPIComm() const { return m_mpiInfo->comm; }
247 
254  void write(const std::string& fileName) const;
255 
260  void Print_Mesh_Info(bool full=false) const;
261 
267  void dump(const std::string& fileName) const;
268 
275  int getTagFromSampleNo(int functionSpaceType, index_t sampleNo) const;
276 
282  const index_t* borrowSampleReferenceIDs(int functionSpaceType) const;
283 
289  virtual bool isValidFunctionSpaceType(int functionSpaceType) const;
290 
295  virtual std::string getDescription() const;
296 
301  virtual std::string functionSpaceTypeAsString(int functionSpaceType) const;
302 
307  void setFunctionSpaceTypeNames();
308 
313  virtual int getContinuousFunctionCode() const;
314 
319  virtual int getReducedContinuousFunctionCode() const;
320 
325  virtual int getFunctionCode() const;
326 
331  virtual int getReducedFunctionCode() const;
332 
337  virtual int getFunctionOnBoundaryCode() const;
338 
343  virtual int getReducedFunctionOnBoundaryCode() const;
344 
349  virtual int getFunctionOnContactZeroCode() const;
350 
355  virtual int getReducedFunctionOnContactZeroCode() const;
356 
361  virtual int getFunctionOnContactOneCode() const;
362 
367  virtual int getReducedFunctionOnContactOneCode() const;
368 
373  virtual int getSolutionCode() const;
374 
379  virtual int getReducedSolutionCode() const;
380 
385  virtual int getDiracDeltaFunctionsCode() const;
386 
390  typedef std::map<int, std::string> FunctionSpaceNamesMapType;
391 
395  virtual int getDim() const { return m_nodes->numDim; }
396 
403  virtual StatusType getStatus() const;
404 
409  virtual dim_t getNumDataPointsGlobal() const;
410 
416  virtual std::pair<int,dim_t> getDataShape(int functionSpaceCode) const;
417 
423  virtual void setToX(escript::Data& arg) const;
424 
431  virtual void setTagMap(const std::string& name, int tag);
432 
438  virtual int getTag(const std::string& name) const;
439 
445  virtual bool isValidTagName(const std::string& name) const;
446 
451  virtual std::string showTagNames() const;
452 
457  virtual void setNewX(const escript::Data& arg);
458 
463  virtual void interpolateOnDomain(escript::Data& target,
464  const escript::Data& source) const;
465 
466  virtual bool probeInterpolationOnDomain(int functionSpaceType_source,
467  int functionSpaceType_target) const;
468 
469  virtual signed char preferredInterpolationOnDomain(int functionSpaceType_source, int functionSpaceType_target) const;
470 
475  bool commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const;
476 
481  virtual void interpolateAcross(escript::Data& target, const escript::Data& source) const;
482 
486  virtual bool probeInterpolationAcross(int functionSpaceType_source,
487  const escript::AbstractDomain& targetDomain,
488  int functionSpaceType_target) const;
489 
495  virtual void setToNormal(escript::Data& out) const;
496 
502  virtual void setToSize(escript::Data& out) const;
503 
509  virtual void setToGradient(escript::Data& grad, const escript::Data& arg) const;
510 
516  virtual void setToIntegrals(std::vector<escript::DataTypes::real_t>& integrals,
517  const escript::Data& arg) const;
518  virtual void setToIntegrals(std::vector<escript::DataTypes::cplx_t>& integrals,
519  const escript::Data& arg) const;
520 
529  virtual int getSystemMatrixTypeId(const boost::python::object& options) const;
530 
540  virtual int getTransportTypeId(int solver, int preconditioner, int package,
541  bool symmetry) const;
542 
548  virtual bool isCellOriented(int functionSpaceCode) const;
549 
550  virtual bool ownSample(int fsCode, index_t id) const;
551 
556  virtual void addPDEToSystem(
558  const escript::Data& A, const escript::Data& B,
559  const escript::Data& C, const escript::Data& D,
560  const escript::Data& X, const escript::Data& Y,
561  const escript::Data& d, const escript::Data& y,
562  const escript::Data& d_contact,
563  const escript::Data& y_contact,
564  const escript::Data& d_dirac,
565  const escript::Data& y_dirac) const;
566 
571  virtual void addPDEToLumpedSystem(escript::Data& mat,
572  const escript::Data& D,
573  const escript::Data& d,
574  const escript::Data& d_dirac,
575  bool useHRZ) const;
576 
581  virtual void addPDEToRHS(escript::Data& rhs, const escript::Data& X,
582  const escript::Data& Y, const escript::Data& y,
583  const escript::Data& y_contact,
584  const escript::Data& y_dirac) const;
585 
590  virtual void addPDEToTransportProblem(
592  escript::Data& source, const escript::Data& M,
593  const escript::Data& A, const escript::Data& B,
594  const escript::Data& C, const escript::Data& D,
595  const escript::Data& X, const escript::Data& Y,
596  const escript::Data& d, const escript::Data& y,
597  const escript::Data& d_contact,
598  const escript::Data& y_contact,
599  const escript::Data& d_dirac,
600  const escript::Data& y_dirac) const;
601 
606  escript::ASM_ptr newSystemMatrix(
607  int row_blocksize,
608  const escript::FunctionSpace& row_functionspace,
609  int column_blocksize,
610  const escript::FunctionSpace& column_functionspace,
611  int type) const;
612 
617  escript::ATP_ptr newTransportProblem(int blocksize,
618  const escript::FunctionSpace& functionspace,
619  int type) const;
620 
624  virtual escript::Data getX() const;
625 
629 #ifdef ESYS_HAVE_BOOST_NUMPY
630  boost::python::numpy::ndarray getNumpyX() const;
631 
632  boost::python::numpy::ndarray getConnectivityInfo() const;
633 #endif
634 
639  virtual escript::Data getNormal() const;
640 
644  virtual escript::Data getSize() const;
645 
649  virtual bool operator==(const escript::AbstractDomain& other) const;
650  virtual bool operator!=(const escript::AbstractDomain& other) const;
651 
656  virtual void setTags(int functionSpaceType, int newTag,
657  const escript::Data& mask) const;
658 
664  virtual int getNumberOfTagsInUse(int functionSpaceCode) const;
665 
666  virtual const int* borrowListOfTagsInUse(int functionSpaceCode) const;
667 
672  virtual bool canTag(int functionSpaceCode) const;
673 
677  virtual int getApproximationOrder(int functionSpaceCode) const;
678 
679  virtual bool supportsContactElements() const { return false; }
680 
681  virtual escript::Data randomFill(const escript::DataTypes::ShapeType& shape,
682  const escript::FunctionSpace& what, long seed,
683  const boost::python::tuple& filter) const;
684 
685  void createMappings(const std::vector<index_t>& dofDistribution,
686  const std::vector<index_t>& nodeDistribution);
687 
690  void relabelElementNodes(const index_t* newNode, index_t offset);
691 
692  template<typename Scalar>
693  void setToIntegralsWorker(std::vector<Scalar>& integrals,
694  const escript::Data& arg) const;
695 
696 #ifdef ESYS_HAVE_PASO
698  paso::SystemMatrixPattern_ptr getPasoPattern() const;
699 #endif
700 
701 #ifdef ESYS_HAVE_TRILINOS
703  esys_trilinos::const_TrilinosGraph_ptr getTrilinosGraph() const {
704  return m_nodes->getTrilinosGraph();
705  }
706 #endif
707 
708 private:
709  void prepare(bool optimize);
710 
718  void resolveNodeIds();
719 
720 #ifdef ESYS_HAVE_PASO
721  paso::SystemMatrixPattern_ptr makePasoPattern() const;
722 #endif
723 
724  void createColoring(const index_t* dofMap);
725  void distributeByRankOfDOF(const IndexVector& distribution);
726  void markNodes(std::vector<short>& mask, index_t offset) const;
727  void optimizeDOFDistribution(IndexVector& distribution);
728  void optimizeDOFLabeling(const IndexVector& distribution);
729  void optimizeElementOrdering();
730  void updateTagList();
731  void printElementInfo(const ElementFile* e, const std::string& title,
732  const std::string& defaultType, bool full) const;
733 
734  void writeElementInfo(std::ostream& stream, const ElementFile* e,
735  const std::string& defaultType) const;
736 
740  std::string m_name;
751 #ifdef ESYS_HAVE_PASO
752  // pointer to the sparse matrix pattern
753  mutable paso::SystemMatrixPattern_ptr pasoPattern;
754 #endif
755 
757 };
758 
759 } // end of namespace
760 
761 #endif // __DUDLEY_DOMAIN_H__
762 
int MPI_Comm
Definition: EsysMPI.h:44
int getTag(unsigned char sourcex, unsigned char sourcey, unsigned char sourcez, unsigned char targetx, unsigned char targety, unsigned char targetz)
Definition: blocktools.cpp:413
DudleyDomain implements the AbstractContinuousDomain interface for the Dudley library.
Definition: DudleyDomain.h:88
virtual int getDim() const
Returns the spatial dimension of the domain.
Definition: DudleyDomain.h:395
virtual escript::JMPI getMPI() const
returns a reference to the MPI information wrapper for this domain
Definition: DudleyDomain.h:220
ElementFile * getFaceElements() const
returns a pointer to this domain's face element file
Definition: DudleyDomain.h:202
virtual int getMPISize() const
returns the number of processors used for this domain
Definition: DudleyDomain.h:226
ElementFile * m_points
the table of points (treated as elements of dimension 0)
Definition: DudleyDomain.h:748
ElementFile * getPoints() const
returns a pointer to this domain's point (nodal) element file
Definition: DudleyDomain.h:214
virtual bool onMasterProcessor() const
returns true if on MPI processor 0, else false
Definition: DudleyDomain.h:244
ElementFile * m_elements
the table of the elements
Definition: DudleyDomain.h:744
virtual void setToIntegrals(std::vector< escript::DataTypes::real_t > &integrals, const escript::Data &arg) const
copies the integrals of the function defined by arg into integrals. arg has to be defined on this.
ElementFile * m_faceElements
the table of face elements
Definition: DudleyDomain.h:746
TagMap m_tagMap
the tag map mapping names to tag keys
Definition: DudleyDomain.h:750
static FunctionSpaceNamesMapType m_functionSpaceTypeNames
Definition: DudleyDomain.h:756
virtual bool supportsContactElements() const
Definition: DudleyDomain.h:679
std::string m_name
domain description
Definition: DudleyDomain.h:740
virtual int getMPIRank() const
returns the number MPI rank of this processor
Definition: DudleyDomain.h:232
std::map< int, std::string > FunctionSpaceNamesMapType
Definition: DudleyDomain.h:390
NodeFile * m_nodes
the table of the nodes
Definition: DudleyDomain.h:742
MPI_Comm getMPIComm() const
get the communicator for this domain. Returns an integer on non-MPI builds Routine must be implemente...
Definition: DudleyDomain.h:246
escript::JMPI m_mpiInfo
MPI information.
Definition: DudleyDomain.h:738
NodeFile * getNodes() const
returns a pointer to this domain's node file
Definition: DudleyDomain.h:178
virtual void setToIntegrals(std::vector< escript::DataTypes::cplx_t > &integrals, const escript::Data &arg) const
ElementFile * getElements() const
returns a pointer to this domain's element file
Definition: DudleyDomain.h:190
Definition: dudley/src/ElementFile.h:53
Definition: dudley/src/NodeFile.h:40
AbstractContinuousDomain, base class for continuous domains.
Definition: AbstractContinuousDomain.h:47
Base class for all escript domains.
Definition: AbstractDomain.h:51
Base class for escript system matrices.
Definition: AbstractSystemMatrix.h:44
Give a short description of what AbstractTransportProblem does.
Definition: AbstractTransportProblem.h:45
Data represents a collection of datapoints.
Definition: Data.h:64
Definition: FunctionSpace.h:36
#define DUDLEY_DLL_API
Definition: dudley/src/system_dep.h:29
A suite of factory methods for creating 2D and 3D dudley domains.
Definition: dudley/src/Assemble.h:32
Domain_ptr readGmsh(const string &fileName, int numDim, int, int, bool optimize)
reads a gmsh mesh file
Definition: dudley/src/DomainFactory.cpp:681
SystemMatrixType
Definition: DudleyDomain.h:75
@ SMT_PASO
Definition: DudleyDomain.h:76
@ SMT_TRILINOS
Definition: DudleyDomain.h:77
@ SMT_COMPLEX
Definition: DudleyDomain.h:78
@ SMT_UNROLL
Definition: DudleyDomain.h:79
std::map< std::string, int > TagMap
Definition: DudleyDomain.h:73
std::vector< index_t > IndexVector
Definition: DataTypes.h:64
std::vector< int > ShapeType
The shape of a single datapoint.
Definition: DataTypes.h:44
index_t dim_t
Definition: DataTypes.h:66
int index_t
type for array/matrix indices used both globally and on each rank
Definition: DataTypes.h:61
Data load(const std::string fileName, const AbstractDomain &domain)
reads Data on domain from file in netCDF format
Definition: DataFactory.cpp:708
boost::shared_ptr< AbstractTransportProblem > ATP_ptr
Definition: AbstractTransportProblem.h:161
boost::shared_ptr< AbstractSystemMatrix > ASM_ptr
Definition: AbstractSystemMatrix.h:33
boost::shared_ptr< AbstractDomain > Domain_ptr
Definition: AbstractDomain.h:41
boost::shared_ptr< JMPI_ > JMPI
Definition: EsysMPI.h:74
double l2(dim_t n, const double *x, escript::JMPI mpiinfo)
returns the global L2 norm of x
Definition: PasoUtil.cpp:501
boost::shared_ptr< SystemMatrixPattern > SystemMatrixPattern_ptr
Definition: SystemMatrixPattern.h:40
static dim_t M
Definition: SparseMatrix_saveHB.cpp:37
bool probeInterpolationAcross(int fsType_source, const escript::AbstractDomain &domain, int fsType_target, int dim)
Definition: CrossDomainCoupler.cpp:32