ESScript  Revision_
ripley/src/RipleyDomain.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2013 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 since 2012 by School of Earth Sciences
13 *
14 *****************************************************************************/
15 
16 #ifndef __RIPLEY_DOMAIN_H__
17 #define __RIPLEY_DOMAIN_H__
18 
19 #include <boost/python/tuple.hpp>
20 
21 #include <ripley/Ripley.h>
22 #include <ripley/RipleyException.h>
23 #include <escript/AbstractContinuousDomain.h>
24 #include <escript/Data.h>
25 #include <escript/FunctionSpace.h>
26 
27 struct Paso_Pattern;
29 struct Paso_SystemMatrix;
30 
31 namespace ripley {
32 
39 class RIPLEY_DLL_API RipleyDomain : public escript::AbstractContinuousDomain
40 {
41 public:
46  RipleyDomain(dim_t dim);
47 
52  ~RipleyDomain();
53 
58  virtual int getMPISize() const { return m_mpiInfo->size; }
59 
64  virtual int getMPIRank() const { return m_mpiInfo->rank; }
65 
70  virtual void MPIBarrier() const {
71 #ifdef ESYS_MPI
72  MPI_Barrier(m_mpiInfo->comm);
73 #endif
74  }
75 
80  virtual bool onMasterProcessor() const { return getMPIRank()==0; }
81 
86 #ifdef ESYS_MPI
87  MPI_Comm
88 #else
89  unsigned int
90 #endif
91  getMPIComm() const {
92 #ifdef ESYS_MPI
93  return m_mpiInfo->comm;
94 #else
95  return 0;
96 #endif
97  }
98 
104  virtual bool isValidFunctionSpaceType(int fsType) const;
105 
110  virtual std::string functionSpaceTypeAsString(int fsType) const;
111 
116  virtual int getDim() const { return m_numDim; }
117 
121  virtual bool operator==(const escript::AbstractDomain& other) const;
122 
126  virtual bool operator!=(const escript::AbstractDomain& other) const {
127  return !(operator==(other));
128  }
129 
136  virtual std::pair<int,int> getDataShape(int fsType) const;
137 
144  int getTagFromSampleNo(int fsType, int sampleNo) const;
145 
152  virtual void setTagMap(const std::string& name, int tag) {
153  m_tagMap[name] = tag;
154  }
155 
161  virtual int getTag(const std::string& name) const {
162  if (m_tagMap.find(name) != m_tagMap.end()) {
163  return m_tagMap.find(name)->second;
164  } else {
165  throw RipleyException("getTag: invalid tag name");
166  }
167  }
168 
174  virtual bool isValidTagName(const std::string& name) const {
175  return (m_tagMap.find(name)!=m_tagMap.end());
176  }
177 
182  virtual std::string showTagNames() const;
183 
189  virtual void setNewX(const escript::Data& arg);
190 
196  virtual void interpolateOnDomain(escript::Data& target, const escript::Data& source) const;
197 
203  virtual bool probeInterpolationOnDomain(int fsType_source, int fsType_target) const;
204 
212  virtual signed char preferredInterpolationOnDomain(int fsType_source, int fsType_target) const;
213 
220  bool
221  commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const;
222 
228  virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
229 
234  virtual bool probeInterpolationACross(int, const escript::AbstractDomain&, int) const;
235 
240  virtual escript::Data getX() const;
241 
246  virtual escript::Data getNormal() const;
247 
251  virtual escript::Data getSize() const;
252 
258  virtual void setToX(escript::Data& arg) const;
259 
266  virtual void setToGradient(escript::Data& out, const escript::Data& in) const;
267 
273  virtual void setTags(const int fsType, const int newTag, const escript::Data& mask) const;
274 
280  virtual bool isCellOriented(int fsType) const;
281 
288  virtual StatusType getStatus() const { return m_status; }
289 
294  virtual int getNumberOfTagsInUse(int fsType) const;
295 
300  virtual const int* borrowListOfTagsInUse(int fsType) const;
301 
306  virtual bool canTag(int fsType) const;
307 
312  virtual int getApproximationOrder(const int fsType) const { return 1; }
313 
318  virtual bool supportsContactElements() const { return false; }
319 
324  virtual int getContinuousFunctionCode() const { return Nodes; }
325 
330  virtual int getReducedContinuousFunctionCode() const { return ReducedNodes; }
331 
336  virtual int getFunctionCode() const { return Elements; }
337 
342  virtual int getReducedFunctionCode() const { return ReducedElements; }
343 
348  virtual int getFunctionOnBoundaryCode() const { return FaceElements; }
349 
356 
361  virtual int getFunctionOnContactZeroCode() const {
362  throw RipleyException("Ripley does not support contact elements");
363  }
364 
370  throw RipleyException("Ripley does not support contact elements");
371  }
372 
377  virtual int getFunctionOnContactOneCode() const {
378  throw RipleyException("Ripley does not support contact elements");
379  }
380 
386  throw RipleyException("Ripley does not support contact elements");
387  }
388 
393  virtual int getSolutionCode() const { return DegreesOfFreedom; }
394 
399  virtual int getReducedSolutionCode() const { return ReducedDegreesOfFreedom; }
400 
405  virtual int getDiracDeltaFunctionsCode() const { return Points; }
406 
417  virtual int getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
418 
428  virtual int getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
429 
435  virtual void setToIntegrals(DoubleVector& integrals, const escript::Data& arg) const;
436 
441  virtual void addPDEToSystem(escript::AbstractSystemMatrix& mat,
442  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
443  const escript::Data& C, const escript::Data& D,
444  const escript::Data& X, const escript::Data& Y,
445  const escript::Data& d, const escript::Data& y,
446  const escript::Data& d_contact, const escript::Data& y_contact,
447  const escript::Data& d_dirac, const escript::Data& y_dirac) const;
448 
453  virtual void addPDEToRHS(escript::Data& rhs, const escript::Data& X,
454  const escript::Data& Y, const escript::Data& y,
455  const escript::Data& y_contact, const escript::Data& y_dirac) const;
456 
461  virtual void addPDEToTransportProblem(escript::AbstractTransportProblem& tp,
462  escript::Data& source, const escript::Data& M,
463  const escript::Data& A, const escript::Data& B,
464  const escript::Data& C, const escript::Data& D,
465  const escript::Data& X, const escript::Data& Y,
466  const escript::Data& d, const escript::Data& y,
467  const escript::Data& d_contact, const escript::Data& y_contact,
468  const escript::Data& d_dirac, const escript::Data& y_dirac) const;
469 
470 
475  virtual escript::ASM_ptr newSystemMatrix(const int row_blocksize,
476  const escript::FunctionSpace& row_functionspace,
477  const int column_blocksize,
478  const escript::FunctionSpace& column_functionspace, const int type) const;
479 
484  virtual escript::ATP_ptr newTransportProblem(
485  const int blocksize, const escript::FunctionSpace& functionspace,
486  const int type) const;
487 
493  virtual void Print_Mesh_Info(const bool full=false) const;
494 
495 
496  /************************************************************************/
497 
503  //void write(const std::string& filename) const = 0;
504 
509  virtual std::string getDescription() const = 0;
510 
516  void dump(const std::string& filename) const = 0;
517 
523  const int* borrowSampleReferenceIDs(int fsType) const = 0;
524 
531  virtual void setToNormal(escript::Data& out) const = 0;
532 
538  virtual void setToSize(escript::Data& out) const = 0;
539 
542  virtual void readNcGrid(escript::Data& out, std::string filename,
543  std::string varname, const std::vector<int>& first,
544  const std::vector<int>& numValues,
545  const std::vector<int>& multiplier) const = 0;
546 
549  virtual void readBinaryGrid(escript::Data& out, std::string filename,
550  const std::vector<int>& first,
551  const std::vector<int>& numValues,
552  const std::vector<int>& multiplier,
553  int byteOrder, int dataType) const = 0;
554 
557  virtual void writeBinaryGrid(const escript::Data& in, std::string filename,
558  int byteOrder, int dataType) const = 0;
559 
564  virtual bool ownSample(int fsType, index_t id) const = 0;
565 
570  virtual int getNumDataPointsGlobal() const = 0;
571 
576  virtual const int* getNumNodesPerDim() const = 0;
577 
582  virtual const int* getNumElementsPerDim() const = 0;
583 
589  virtual const int* getNumFacesPerBoundary() const = 0;
590 
595  virtual IndexVector getNodeDistribution() const = 0;
596 
601  virtual const int* getNumSubdivisionsPerDim() const = 0;
602 
607  virtual double getLocalCoordinate(int index, int dim) const = 0;
608 
613  virtual boost::python::tuple getGridParameters() const = 0;
614 
615 
619  virtual bool supportsFilter(const boost::python::tuple& t) const;
620 
624  virtual escript::Data randomFill(long seed, const boost::python::tuple& filter) const;
625 
626 
627 protected:
629  StatusType m_status;
632  mutable IndexVector m_nodeTags, m_nodeTagsInUse;
633  mutable IndexVector m_elementTags, m_elementTagsInUse;
634  mutable IndexVector m_faceTags, m_faceTagsInUse;
635 
637  void copyData(escript::Data& out, escript::Data& in) const;
638 
640  void averageData(escript::Data& out, escript::Data& in) const;
641 
643  void multiplyData(escript::Data& out, escript::Data& in) const;
644 
645  // this is const because setTags is const
646  void updateTagsInUse(int fsType) const;
647 
649  Paso_Pattern* createPasoPattern(const IndexVector& ptr,
650  const IndexVector& index, const dim_t M, const dim_t N) const;
651 
653  Paso_Pattern* createMainPattern() const;
654 
658  void createCouplePatterns(const std::vector<IndexVector>& colIndices,
659  const dim_t N, Paso_Pattern** colPattern,
660  Paso_Pattern** rowPattern) const;
661 
662  void addToSystemMatrix(Paso_SystemMatrix* in, const IndexVector& nodes_Eq,
663  dim_t num_Eq, const IndexVector& nodes_Sol, dim_t num_Sol,
664  const DoubleVector& array) const;
665 
666  /***********************************************************************/
667 
669  virtual dim_t getNumNodes() const = 0;
670 
672  virtual dim_t getNumElements() const = 0;
673 
675  virtual dim_t getNumDOF() const = 0;
676 
678  virtual dim_t getNumFaceElements() const = 0;
679 
682  virtual dim_t insertNeighbourNodes(IndexVector& index, index_t node) const = 0;
683 
685  virtual void assembleCoordinates(escript::Data& arg) const = 0;
686 
688  virtual void assembleGradient(escript::Data& out, escript::Data& in) const = 0;
689 
691  virtual void assembleIntegrate(DoubleVector& integrals, escript::Data& arg) const = 0;
692 
695  virtual void assemblePDESingle(Paso_SystemMatrix* mat, escript::Data& rhs,
696  const escript::Data& A, const escript::Data& B,
697  const escript::Data& C, const escript::Data& D,
698  const escript::Data& X, const escript::Data& Y) const = 0;
699 
702  virtual void assemblePDEBoundarySingle(Paso_SystemMatrix* mat,
703  escript::Data& rhs, const escript::Data& d,
704  const escript::Data& y) const = 0;
705 
708  virtual void assemblePDESingleReduced(Paso_SystemMatrix* mat,
709  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
710  const escript::Data& C, const escript::Data& D,
711  const escript::Data& X, const escript::Data& Y) const = 0;
712 
715  virtual void assemblePDEBoundarySingleReduced(Paso_SystemMatrix* mat,
716  escript::Data& rhs, const escript::Data& d,
717  const escript::Data& y) const = 0;
718 
721  virtual void assemblePDESystem(Paso_SystemMatrix* mat, escript::Data& rhs,
722  const escript::Data& A, const escript::Data& B,
723  const escript::Data& C, const escript::Data& D,
724  const escript::Data& X, const escript::Data& Y) const = 0;
725 
728  virtual void assemblePDEBoundarySystem(Paso_SystemMatrix* mat,
729  escript::Data& rhs, const escript::Data& d,
730  const escript::Data& y) const = 0;
731 
734  virtual void assemblePDESystemReduced(Paso_SystemMatrix* mat,
735  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
736  const escript::Data& C, const escript::Data& D,
737  const escript::Data& X, const escript::Data& Y) const = 0;
738 
741  virtual void assemblePDEBoundarySystemReduced(Paso_SystemMatrix* mat,
742  escript::Data& rhs, const escript::Data& d,
743  const escript::Data& y) const = 0;
744 
746  virtual Paso_SystemMatrixPattern* getPattern(bool reducedRowOrder,
747  bool reducedColOrder) const = 0;
748 
750  virtual void interpolateNodesOnElements(escript::Data& out,
751  escript::Data& in, bool reduced) const = 0;
752 
754  virtual void interpolateNodesOnFaces(escript::Data& out, escript::Data& in,
755  bool reduced) const = 0;
756 
758  virtual void nodesToDOF(escript::Data& out, escript::Data& in) const = 0;
759 
761  virtual void dofToNodes(escript::Data& out, escript::Data& in) const = 0;
762 
763 private:
765  void assemblePDE(Paso_SystemMatrix* mat, escript::Data& rhs,
766  const escript::Data& A, const escript::Data& B,
767  const escript::Data& C, const escript::Data& D,
768  const escript::Data& X, const escript::Data& Y) const;
769 
772  void assemblePDEBoundary(Paso_SystemMatrix* mat, escript::Data& rhs,
773  const escript::Data& d, const escript::Data& y) const;
774 };
775 
776 } // end of namespace ripley
777 
778 #endif // __RIPLEY_DOMAIN_H__
779