escript  Revision_4925
ripley/src/RipleyDomain.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2014 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16 
17 #ifndef __RIPLEY_DOMAIN_H__
18 #define __RIPLEY_DOMAIN_H__
19 
20 #ifdef BADPYTHONMACROS
21 // This hack is required for BSD/OSX builds with python 2.7
22 // (and possibly others). It must be the first include.
23 // From bug reports online it seems that python redefines
24 // some c macros that are functions in c++.
25 // c++ doesn't like that!
26 #include <Python.h>
27 #undef BADPYTHONMACROS
28 #endif
29 
30 #include <boost/python/tuple.hpp>
31 #include <boost/python/list.hpp>
32 
33 #include <ripley/Ripley.h>
34 #include <ripley/RipleyException.h>
35 #include <ripley/AbstractAssembler.h>
36 
37 #include <escript/AbstractContinuousDomain.h>
38 #include <escript/Data.h>
39 #include <escript/FunctionSpace.h>
40 
41 #include <paso/SystemMatrix.h>
42 
43 namespace ripley {
44 
49 };
50 
51 /* There is no particular significance to this type,
52 It is here as a typedef because a bug in clang++ prevents
53 that compiler from recognising it as a valid part of
54 a constant expression.
55 */
56 typedef std::map<std::string, int> simap_t;
57 
58 
64 {
66  std::vector<int> first;
68  std::vector<int> numValues;
71  std::vector<int> multiplier;
73  std::vector<int> reverse;
75  int byteOrder;
77  int dataType;
78 };
79 
84 struct DiracPoint
85 {
86  int node;
87  int tag;
88 };
89 
97 {
98 public:
103  RipleyDomain(dim_t dim);
104 
109  ~RipleyDomain();
110 
115  virtual int getMPISize() const { return m_mpiInfo->size; }
116 
121  virtual int getMPIRank() const { return m_mpiInfo->rank; }
122 
127  virtual void MPIBarrier() const {
128 #ifdef ESYS_MPI
129  MPI_Barrier(m_mpiInfo->comm);
130 #endif
131  }
132 
137  virtual bool onMasterProcessor() const { return getMPIRank()==0; }
138 
143 #ifdef ESYS_MPI
144  MPI_Comm
145 #else
146  unsigned int
147 #endif
148  getMPIComm() const {
149 #ifdef ESYS_MPI
150  return m_mpiInfo->comm;
151 #else
152  return 0;
153 #endif
154  }
155 
161  virtual bool isValidFunctionSpaceType(int fsType) const;
162 
167  virtual std::string functionSpaceTypeAsString(int fsType) const;
168 
173  virtual int getDim() const { return m_numDim; }
174 
178  virtual bool operator==(const escript::AbstractDomain& other) const;
179 
183  virtual bool operator!=(const escript::AbstractDomain& other) const {
184  return !(operator==(other));
185  }
186 
193  virtual std::pair<int,int> getDataShape(int fsType) const;
194 
201  int getTagFromSampleNo(int fsType, int sampleNo) const;
202 
209  virtual void setTagMap(const std::string& name, int tag) {
210  m_tagMap[name] = tag;
211  }
212 
218  virtual int getTag(const std::string& name) const {
219  if (m_tagMap.find(name) != m_tagMap.end()) {
220  return m_tagMap.find(name)->second;
221  } else {
222  throw RipleyException("getTag: invalid tag name");
223  }
224  }
225 
231  virtual bool isValidTagName(const std::string& name) const {
232  return (m_tagMap.find(name)!=m_tagMap.end());
233  }
234 
239  virtual std::string showTagNames() const;
240 
246  virtual void setNewX(const escript::Data& arg);
247 
253  virtual void interpolateOnDomain(escript::Data& target, const escript::Data& source) const;
254 
260  virtual bool probeInterpolationOnDomain(int fsType_source, int fsType_target) const;
261 
269  virtual signed char preferredInterpolationOnDomain(int fsType_source, int fsType_target) const;
270 
277  bool
278  commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const;
279 
285  virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
286 
291  virtual bool probeInterpolationACross(int, const escript::AbstractDomain&, int) const;
292 
297  virtual escript::Data getX() const;
298 
303  virtual escript::Data getNormal() const;
304 
308  virtual escript::Data getSize() const;
309 
315  virtual void setToX(escript::Data& arg) const;
316 
323  virtual void setToGradient(escript::Data& out, const escript::Data& in) const;
324 
330  virtual void setTags(const int fsType, const int newTag, const escript::Data& mask) const;
331 
337  virtual bool isCellOriented(int fsType) const;
338 
345  virtual StatusType getStatus() const { return m_status; }
346 
351  virtual int getNumberOfTagsInUse(int fsType) const;
352 
357  virtual const int* borrowListOfTagsInUse(int fsType) const;
358 
363  virtual bool canTag(int fsType) const;
364 
369  virtual int getApproximationOrder(const int fsType) const { return 1; }
370 
375  virtual bool supportsContactElements() const { return false; }
376 
381  virtual int getContinuousFunctionCode() const { return Nodes; }
382 
387  virtual int getReducedContinuousFunctionCode() const { return ReducedNodes; }
388 
393  virtual int getFunctionCode() const { return Elements; }
394 
399  virtual int getReducedFunctionCode() const { return ReducedElements; }
400 
405  virtual int getFunctionOnBoundaryCode() const { return FaceElements; }
406 
413 
418  virtual int getFunctionOnContactZeroCode() const {
419  throw RipleyException("Ripley does not support contact elements");
420  }
421 
427  throw RipleyException("Ripley does not support contact elements");
428  }
429 
434  virtual int getFunctionOnContactOneCode() const {
435  throw RipleyException("Ripley does not support contact elements");
436  }
437 
443  throw RipleyException("Ripley does not support contact elements");
444  }
445 
450  virtual int getSolutionCode() const { return DegreesOfFreedom; }
451 
456  virtual int getReducedSolutionCode() const { return ReducedDegreesOfFreedom; }
457 
462  virtual int getDiracDeltaFunctionsCode() const { return Points; }
463 
474  virtual int getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
475 
485  virtual int getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
486 
492  virtual void setToIntegrals(DoubleVector& integrals, const escript::Data& arg) const;
493 
498  virtual void addPDEToSystem(escript::AbstractSystemMatrix& mat,
499  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
500  const escript::Data& C, const escript::Data& D,
501  const escript::Data& X, const escript::Data& Y,
502  const escript::Data& d, const escript::Data& y,
503  const escript::Data& d_contact, const escript::Data& y_contact,
504  const escript::Data& d_dirac, const escript::Data& y_dirac) const;
505 
511  virtual void addToSystem(escript::AbstractSystemMatrix& mat,
512  escript::Data& rhs,std::map<std::string, escript::Data> data) const;
513 
518  virtual void addToSystemFromPython(escript::AbstractSystemMatrix& mat,
519  escript::Data& rhs,boost::python::list data) const;
520 
525  virtual void addPDEToRHS(escript::Data& rhs, const escript::Data& X,
526  const escript::Data& Y, const escript::Data& y,
527  const escript::Data& y_contact, const escript::Data& y_dirac) const;
528 
534  virtual void addToRHS(escript::Data& rhs,
535  std::map<std::string, escript::Data> data) const;
536 
541  virtual void addToRHSFromPython(escript::Data& rhs,
542  boost::python::list data) const;
543 
548  virtual void addPDEToTransportProblem(escript::AbstractTransportProblem& tp,
549  escript::Data& source, const escript::Data& M,
550  const escript::Data& A, const escript::Data& B,
551  const escript::Data& C, const escript::Data& D,
552  const escript::Data& X, const escript::Data& Y,
553  const escript::Data& d, const escript::Data& y,
554  const escript::Data& d_contact, const escript::Data& y_contact,
555  const escript::Data& d_dirac, const escript::Data& y_dirac) const;
556 
557 
562  virtual escript::ASM_ptr newSystemMatrix(const int row_blocksize,
563  const escript::FunctionSpace& row_functionspace,
564  const int column_blocksize,
565  const escript::FunctionSpace& column_functionspace, const int type) const;
566 
571  virtual escript::ATP_ptr newTransportProblem(
572  const int blocksize, const escript::FunctionSpace& functionspace,
573  const int type) const;
574 
580  virtual void Print_Mesh_Info(const bool full=false) const;
581 
582 
583  /************************************************************************/
584 
590  //void write(const std::string& filename) const = 0;
591 
596  virtual std::string getDescription() const = 0;
597 
603  void dump(const std::string& filename) const = 0;
604 
610  const int* borrowSampleReferenceIDs(int fsType) const = 0;
611 
618  virtual void setToNormal(escript::Data& out) const = 0;
619 
625  virtual void setToSize(escript::Data& out) const = 0;
626 
629  virtual void readNcGrid(escript::Data& out, std::string filename,
630  std::string varname, const ReaderParameters& params) const = 0;
631 
634  virtual void readBinaryGrid(escript::Data& out, std::string filename,
635  const ReaderParameters& params) const = 0;
636 #ifdef USE_BOOSTIO
637  virtual void readBinaryGridFromZipped(escript::Data& out, std::string filename,
638  const ReaderParameters& params) const = 0;
639 #endif
640 
643  virtual void writeBinaryGrid(const escript::Data& in, std::string filename,
644  int byteOrder, int dataType) const = 0;
645 
650  virtual bool ownSample(int fsType, index_t id) const = 0;
651 
656  virtual int getNumDataPointsGlobal() const = 0;
657 
662  virtual const int* getNumNodesPerDim() const = 0;
663 
668  virtual const int* getNumElementsPerDim() const = 0;
669 
675  virtual const int* getNumFacesPerBoundary() const = 0;
676 
681  virtual IndexVector getNodeDistribution() const = 0;
682 
687  virtual const int* getNumSubdivisionsPerDim() const = 0;
688 
693  virtual double getLocalCoordinate(int index, int dim) const = 0;
694 
699  virtual boost::python::tuple getGridParameters() const = 0;
700 
701 
706  virtual bool supportsFilter(const boost::python::tuple& t) const;
707 
708  virtual void setAssembler(std::string type,
709  std::map<std::string, escript::Data> options) {
710  throw RipleyException("Domain does not support custom assemblers");
711  }
712 
713  void setAssemblerFromPython(std::string type,
714  boost::python::list options);
715 protected:
720  mutable IndexVector m_nodeTags, m_nodeTagsInUse;
721  mutable IndexVector m_elementTags, m_elementTagsInUse;
722  mutable IndexVector m_faceTags, m_faceTagsInUse;
724  std::vector<struct DiracPoint> m_diracPoints;
725  IndexVector m_diracPointNodeIDs; //for borrowSampleID
727 
729  void copyData(escript::Data& out, const escript::Data& in) const;
730 
732  void averageData(escript::Data& out, const escript::Data& in) const;
733 
735  void multiplyData(escript::Data& out, const escript::Data& in) const;
736 
737  // this is const because setTags is const
738  void updateTagsInUse(int fsType) const;
739 
741  paso::Pattern_ptr createPasoPattern(const IndexVector& ptr,
742  const IndexVector& index, const dim_t M, const dim_t N) const;
743 
745  paso::Pattern_ptr createMainPattern() const;
746 
750  void createCouplePatterns(const std::vector<IndexVector>& colIndices,
751  const std::vector<IndexVector>& rowIndices,
752  const dim_t N, paso::Pattern_ptr& colPattern,
753  paso::Pattern_ptr& rowPattern) const;
754 
755  void addToSystemMatrix(paso::SystemMatrix_ptr in,
756  const IndexVector& nodes_Eq, dim_t num_Eq,
757  const IndexVector& nodes_Sol, dim_t num_Sol,
758  const DoubleVector& array) const;
759 
760  void addPoints(int numPoints, const double* points_ptr,
761  const int* tags_ptr);
762 
763  /***********************************************************************/
764 
766  virtual dim_t getNumNodes() const = 0;
767 
769  virtual dim_t getNumElements() const = 0;
770 
772  virtual dim_t getNumDOF() const = 0;
773 
775  virtual dim_t getNumFaceElements() const = 0;
776 
779  virtual dim_t insertNeighbourNodes(IndexVector& index, index_t node) const = 0;
780 
782  virtual void assembleCoordinates(escript::Data& arg) const = 0;
783 
785  virtual void assembleGradient(escript::Data& out, const escript::Data& in) const = 0;
786 
788  virtual void assembleIntegrate(DoubleVector& integrals, const escript::Data& arg) const = 0;
789 
791  virtual paso::SystemMatrixPattern_ptr getPattern(bool reducedRowOrder,
792  bool reducedColOrder) const = 0;
793 
795  virtual void interpolateNodesOnElements(escript::Data& out,
796  const escript::Data& in,
797  bool reduced) const = 0;
798 
800  virtual void interpolateNodesOnFaces(escript::Data& out,
801  const escript::Data& in,
802  bool reduced) const = 0;
803 
805  virtual void nodesToDOF(escript::Data& out, const escript::Data& in) const = 0;
806 
808  virtual void dofToNodes(escript::Data& out, const escript::Data& in) const = 0;
809 
810  virtual int getDofOfNode(int node) const = 0;
811 
812 private:
814  void assemblePDE(paso::SystemMatrix_ptr mat, escript::Data& rhs,
815  std::map<std::string, escript::Data> coefs) const;
816 
819  void assemblePDEBoundary(paso::SystemMatrix_ptr mat, escript::Data& rhs,
820  std::map<std::string, escript::Data> coefs) const;
821 
822  void assemblePDEDirac(paso::SystemMatrix_ptr mat, escript::Data& rhs,
823  std::map<std::string, escript::Data> coefs) const;
824 
825  // finds the node that the given point belongs to
826  virtual int findNode(const double *coords) const = 0;
827 };
828 
829 } // end of namespace ripley
830 
831 #endif // __RIPLEY_DOMAIN_H__
832 
AbstractContinuousDomain, base class for continuous domains.
Definition: AbstractContinuousDomain.h:46
Give a short description of what FunctionSpace does.
Definition: FunctionSpace.h:46
boost::shared_ptr< Pattern > Pattern_ptr
Definition: Pattern.h:36
int dataType
data type in the file (used by binary reader only)
Definition: ripley/src/RipleyDomain.h:77
virtual void setAssembler(std::string type, std::map< std::string, escript::Data > options)
Definition: ripley/src/RipleyDomain.h:708
std::vector< struct DiracPoint > m_diracPoints
Definition: ripley/src/RipleyDomain.h:724
IndexVector m_faceTagsInUse
Definition: ripley/src/RipleyDomain.h:722
virtual bool onMasterProcessor() const
returns true if on MPI processor 0, else false
Definition: ripley/src/RipleyDomain.h:137
virtual int getFunctionOnBoundaryCode() const
returns a function on boundary FunctionSpace code
Definition: ripley/src/RipleyDomain.h:405
Definition: Ripley.h:49
escript::Data readNcGrid(std::string filename, std::string varname, escript::FunctionSpace fs, const object &pyShape, double fill, const object &pyFirst, const object &pyNum, const object &pyMultiplier, const object &pyReverse)
Definition: ripleycpp.cpp:109
virtual int getFunctionOnContactOneCode() const
returns a FunctionOnContactOne code
Definition: ripley/src/RipleyDomain.h:434
unsigned int getMPIComm() const
returns the MPI communicator
Definition: ripley/src/RipleyDomain.h:148
virtual bool operator!=(const escript::AbstractDomain &other) const
inequality operator
Definition: ripley/src/RipleyDomain.h:183
Definition: ripley/src/RipleyDomain.h:47
Definition: AbstractAssembler.h:35
int node
Definition: ripley/src/RipleyDomain.h:86
Struct that holds MPI communicator, rank, size and a tag counter.
Definition: Esys_MPI.h:48
int byteOrder
byte order in the file (used by binary reader only)
Definition: ripley/src/RipleyDomain.h:75
A struct to contain a dirac point's information.
Definition: ripley/src/RipleyDomain.h:84
Definition: ripley/src/RipleyDomain.h:46
virtual int getReducedContinuousFunctionCode() const
returns a continuous on reduced order nodes FunctionSpace code
Definition: ripley/src/RipleyDomain.h:387
Definition: Ripley.h:52
assembler_t assembler_type
Definition: ripley/src/RipleyDomain.h:726
RipleyException exception class.
Definition: RipleyException.h:29
Esys_MPIInfo * m_mpiInfo
Definition: ripley/src/RipleyDomain.h:718
Definition: Ripley.h:50
static dim_t M
Definition: SparseMatrix_saveHB.cpp:36
virtual bool supportsContactElements() const
returns true if this domain supports contact elements, false otherwise
Definition: ripley/src/RipleyDomain.h:375
std::vector< int > numValues
the number of values to read from file
Definition: ripley/src/RipleyDomain.h:68
std::vector< int > reverse
if non-zero, values are written from last index to first index
Definition: ripley/src/RipleyDomain.h:73
escript::Data readBinaryGrid(std::string filename, escript::FunctionSpace fs, const object &pyShape, double fill, int byteOrder, int dataType, const object &pyFirst, const object &pyNum, const object &pyMultiplier, const object &pyReverse)
Definition: ripleycpp.cpp:57
assembler_t
Definition: ripley/src/RipleyDomain.h:45
virtual int getFunctionOnContactZeroCode() const
return a FunctionOnContactZero code
Definition: ripley/src/RipleyDomain.h:418
virtual int getApproximationOrder(const int fsType) const
returns the approximation order used for a function space
Definition: ripley/src/RipleyDomain.h:369
boost::shared_ptr< SystemMatrixPattern > SystemMatrixPattern_ptr
Definition: SystemMatrixPattern.h:38
boost::shared_ptr< AbstractTransportProblem > ATP_ptr
Definition: AbstractTransportProblem.h:162
int MPI_Comm
Definition: Esys_MPI.h:29
virtual int getDim() const
returns the number of spatial dimensions of the domain
Definition: ripley/src/RipleyDomain.h:173
boost::shared_ptr< SystemMatrix > SystemMatrix_ptr
Definition: SystemMatrix.h:38
IndexVector m_diracPointNodeIDs
Definition: ripley/src/RipleyDomain.h:725
virtual int getMPIRank() const
returns the MPI rank of this processor
Definition: ripley/src/RipleyDomain.h:121
std::vector< int > first
the (global) offset into the data object to start writing into
Definition: ripley/src/RipleyDomain.h:66
Definition: Ripley.h:45
std::map< std::string, int > simap_t
Definition: ripley/src/RipleyDomain.h:56
std::vector< index_t > IndexVector
Definition: Ripley.h:38
StatusType m_status
Definition: ripley/src/RipleyDomain.h:717
IndexVector m_nodeTagsInUse
Definition: ripley/src/RipleyDomain.h:720
virtual int getReducedFunctionCode() const
returns a function with reduced integration order FunctionSpace code
Definition: ripley/src/RipleyDomain.h:399
Structure that wraps parameters for the grid reading routines.
Definition: ripley/src/RipleyDomain.h:63
std::vector< int > multiplier
Definition: ripley/src/RipleyDomain.h:71
virtual int getReducedFunctionOnContactOneCode() const
returns a FunctionOnContactOne code with reduced integration order
Definition: ripley/src/RipleyDomain.h:442
int StatusType
Definition: AbstractDomain.h:78
virtual int getDiracDeltaFunctionsCode() const
returns a DiracDeltaFunctions FunctionSpace code
Definition: ripley/src/RipleyDomain.h:462
Data represents a collection of datapoints.
Definition: Data.h:71
virtual int getReducedFunctionOnContactZeroCode() const
returns a FunctionOnContactZero code with reduced integration order
Definition: ripley/src/RipleyDomain.h:426
virtual void setTagMap(const std::string &name, int tag)
sets a map from a clear tag name to a tag key
Definition: ripley/src/RipleyDomain.h:209
static dim_t N
Definition: SparseMatrix_saveHB.cpp:36
RipleyDomain extends the AbstractContinuousDomain interface for the Ripley library and is the base cl...
Definition: ripley/src/RipleyDomain.h:96
virtual bool isValidTagName(const std::string &name) const
returns true if name is a defined tag name
Definition: ripley/src/RipleyDomain.h:231
Definition: Ripley.h:48
virtual int getReducedFunctionOnBoundaryCode() const
returns a function on boundary with reduced integration order FunctionSpace code
Definition: ripley/src/RipleyDomain.h:412
virtual int getMPISize() const
returns the number of processors used for this domain
Definition: ripley/src/RipleyDomain.h:115
virtual int getTag(const std::string &name) const
returns the tag key for tag name
Definition: ripley/src/RipleyDomain.h:218
int index_t
Definition: types.h:25
AbstractAssembler * assembler
Definition: ripley/src/RipleyDomain.h:723
virtual int getReducedSolutionCode() const
returns a ReducedSolution FunctionSpace code
Definition: ripley/src/RipleyDomain.h:456
virtual int getFunctionCode() const
returns a function FunctionSpace code
Definition: ripley/src/RipleyDomain.h:393
int tag
Definition: ripley/src/RipleyDomain.h:87
Give a short description of what AbstractTransportProblem does.
Definition: AbstractTransportProblem.h:45
Definition: ripley/src/RipleyDomain.h:48
Definition: Ripley.h:46
Definition: Ripley.h:51
Give a short description of what AbstractSystemMatrix does.
Definition: AbstractSystemMatrix.h:44
IndexVector m_elementTagsInUse
Definition: ripley/src/RipleyDomain.h:721
virtual StatusType getStatus() const
returns a status indicator of the domain. The status identifier should be unique over the lifetime of...
Definition: ripley/src/RipleyDomain.h:345
boost::shared_ptr< AbstractSystemMatrix > ASM_ptr
Definition: AbstractSystemMatrix.h:201
Definition: AbstractDomain.h:62
TagMap m_tagMap
Definition: ripley/src/RipleyDomain.h:719
virtual int getSolutionCode() const
returns a Solution FunctionSpace code
Definition: ripley/src/RipleyDomain.h:450
std::map< std::string, index_t > TagMap
Definition: Ripley.h:41
Definition: Ripley.h:44
dim_t m_numDim
Definition: ripley/src/RipleyDomain.h:716
#define RIPLEY_DLL_API
Definition: ripley/src/system_dep.h:35
int dim_t
Definition: types.h:24
Definition: Ripley.h:47
std::vector< double > DoubleVector
Definition: Ripley.h:39
virtual int getContinuousFunctionCode() const
returns a continuous FunctionSpace code
Definition: ripley/src/RipleyDomain.h:381
virtual void MPIBarrier() const
if compiled for MPI then executes an MPI_Barrier, else does nothing
Definition: ripley/src/RipleyDomain.h:127