Escript  Revision_4320
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 <ripley/Ripley.h>
20 #include <ripley/RipleyException.h>
21 #include <escript/AbstractContinuousDomain.h>
22 #include <escript/Data.h>
23 #include <escript/FunctionSpace.h>
24 
25 struct Paso_Pattern;
27 struct Paso_SystemMatrix;
28 
29 namespace ripley {
30 
38 {
39 public:
44  RipleyDomain(dim_t dim);
45 
50  ~RipleyDomain();
51 
56  virtual int getMPISize() const { return m_mpiInfo->size; }
57 
62  virtual int getMPIRank() const { return m_mpiInfo->rank; }
63 
68  virtual void MPIBarrier() const {
69 #ifdef ESYS_MPI
70  MPI_Barrier(m_mpiInfo->comm);
71 #endif
72  }
73 
78  virtual bool onMasterProcessor() const { return getMPIRank()==0; }
79 
84 #ifdef ESYS_MPI
85  MPI_Comm
86 #else
87  unsigned int
88 #endif
89  getMPIComm() const {
90 #ifdef ESYS_MPI
91  return m_mpiInfo->comm;
92 #else
93  return 0;
94 #endif
95  }
96 
102  virtual bool isValidFunctionSpaceType(int fsType) const;
103 
108  virtual std::string functionSpaceTypeAsString(int fsType) const;
109 
114  virtual int getDim() const { return m_numDim; }
115 
119  virtual bool operator==(const escript::AbstractDomain& other) const;
120 
124  virtual bool operator!=(const escript::AbstractDomain& other) const {
125  return !(operator==(other));
126  }
127 
134  virtual std::pair<int,int> getDataShape(int fsType) const;
135 
142  int getTagFromSampleNo(int fsType, int sampleNo) const;
143 
150  virtual void setTagMap(const std::string& name, int tag) {
151  m_tagMap[name] = tag;
152  }
153 
159  virtual int getTag(const std::string& name) const {
160  if (m_tagMap.find(name) != m_tagMap.end()) {
161  return m_tagMap.find(name)->second;
162  } else {
163  throw RipleyException("getTag: invalid tag name");
164  }
165  }
166 
172  virtual bool isValidTagName(const std::string& name) const {
173  return (m_tagMap.find(name)!=m_tagMap.end());
174  }
175 
180  virtual std::string showTagNames() const;
181 
187  virtual void setNewX(const escript::Data& arg);
188 
194  virtual void interpolateOnDomain(escript::Data& target, const escript::Data& source) const;
195 
201  virtual bool probeInterpolationOnDomain(int fsType_source, int fsType_target) const;
202 
210  virtual signed char preferredInterpolationOnDomain(int fsType_source, int fsType_target) const;
211 
218  bool
219  commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const;
220 
226  virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
227 
232  virtual bool probeInterpolationACross(int, const escript::AbstractDomain&, int) const;
233 
238  virtual escript::Data getX() const;
239 
244  virtual escript::Data getNormal() const;
245 
249  virtual escript::Data getSize() const;
250 
256  virtual void setToX(escript::Data& arg) const;
257 
264  virtual void setToGradient(escript::Data& out, const escript::Data& in) const;
265 
271  virtual void setTags(const int fsType, const int newTag, const escript::Data& mask) const;
272 
278  virtual bool isCellOriented(int fsType) const;
279 
286  virtual StatusType getStatus() const { return m_status; }
287 
292  virtual int getNumberOfTagsInUse(int fsType) const;
293 
298  virtual const int* borrowListOfTagsInUse(int fsType) const;
299 
304  virtual bool canTag(int fsType) const;
305 
310  virtual int getApproximationOrder(const int fsType) const { return 1; }
311 
316  virtual bool supportsContactElements() const { return false; }
317 
322  virtual int getContinuousFunctionCode() const { return Nodes; }
323 
328  virtual int getReducedContinuousFunctionCode() const { return ReducedNodes; }
329 
334  virtual int getFunctionCode() const { return Elements; }
335 
340  virtual int getReducedFunctionCode() const { return ReducedElements; }
341 
346  virtual int getFunctionOnBoundaryCode() const { return FaceElements; }
347 
354 
359  virtual int getFunctionOnContactZeroCode() const {
360  throw RipleyException("Ripley does not support contact elements");
361  }
362 
368  throw RipleyException("Ripley does not support contact elements");
369  }
370 
375  virtual int getFunctionOnContactOneCode() const {
376  throw RipleyException("Ripley does not support contact elements");
377  }
378 
384  throw RipleyException("Ripley does not support contact elements");
385  }
386 
391  virtual int getSolutionCode() const { return DegreesOfFreedom; }
392 
397  virtual int getReducedSolutionCode() const { return ReducedDegreesOfFreedom; }
398 
403  virtual int getDiracDeltaFunctionsCode() const { return Points; }
404 
415  virtual int getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
416 
426  virtual int getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
427 
433  virtual void setToIntegrals(std::vector<double>& integrals, const escript::Data& arg) const;
434 
439  virtual void addPDEToSystem(escript::AbstractSystemMatrix& mat,
440  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
441  const escript::Data& C, const escript::Data& D,
442  const escript::Data& X, const escript::Data& Y,
443  const escript::Data& d, const escript::Data& y,
444  const escript::Data& d_contact, const escript::Data& y_contact,
445  const escript::Data& d_dirac, const escript::Data& y_dirac) const;
446 
451  virtual void addPDEToRHS(escript::Data& rhs, const escript::Data& X,
452  const escript::Data& Y, const escript::Data& y,
453  const escript::Data& y_contact, const escript::Data& y_dirac) const;
454 
459  virtual void addPDEToTransportProblem(escript::AbstractTransportProblem& tp,
460  escript::Data& source, const escript::Data& M,
461  const escript::Data& A, const escript::Data& B,
462  const escript::Data& C, const escript::Data& D,
463  const escript::Data& X, const escript::Data& Y,
464  const escript::Data& d, const escript::Data& y,
465  const escript::Data& d_contact, const escript::Data& y_contact,
466  const escript::Data& d_dirac, const escript::Data& y_dirac) const;
467 
468 
473  virtual escript::ASM_ptr newSystemMatrix(const int row_blocksize,
474  const escript::FunctionSpace& row_functionspace,
475  const int column_blocksize,
476  const escript::FunctionSpace& column_functionspace, const int type) const;
477 
482  virtual escript::ATP_ptr newTransportProblem(
483  const int blocksize, const escript::FunctionSpace& functionspace,
484  const int type) const;
485 
491  virtual void Print_Mesh_Info(const bool full=false) const;
492 
493 
494  /************************************************************************/
495 
501  //void write(const std::string& filename) const = 0;
502 
507  virtual std::string getDescription() const = 0;
508 
514  void dump(const std::string& filename) const = 0;
515 
521  const int* borrowSampleReferenceIDs(int fsType) const = 0;
522 
529  virtual void setToNormal(escript::Data& out) const = 0;
530 
536  virtual void setToSize(escript::Data& out) const = 0;
537 
540  virtual void readBinaryGrid(escript::Data& out, std::string filename,
541  const std::vector<int>& first,
542  const std::vector<int>& numValues,
543  const std::vector<int>& multiplier) const = 0;
544 
547  virtual void readNcGrid(escript::Data& out, std::string filename,
548  std::string varname, const std::vector<int>& first,
549  const std::vector<int>& numValues,
550  const std::vector<int>& multiplier) const = 0;
551 
556  virtual bool ownSample(int fsType, index_t id) const = 0;
557 
562  virtual int getNumDataPointsGlobal() const = 0;
563 
568  virtual IndexVector getNumNodesPerDim() const = 0;
569 
574  virtual IndexVector getNumElementsPerDim() const = 0;
575 
581  virtual IndexVector getNumFacesPerBoundary() const = 0;
582 
587  virtual IndexVector getNodeDistribution() const = 0;
588 
593  virtual IndexVector getNumSubdivisionsPerDim() const = 0;
594 
600  virtual std::pair<double,double> getFirstCoordAndSpacing(dim_t dim) const = 0;
601 
602 protected:
607  mutable IndexVector m_nodeTags, m_nodeTagsInUse;
608  mutable IndexVector m_elementTags, m_elementTagsInUse;
609  mutable IndexVector m_faceTags, m_faceTagsInUse;
610 
612  void copyData(escript::Data& out, escript::Data& in) const;
613 
615  void averageData(escript::Data& out, escript::Data& in) const;
616 
618  void multiplyData(escript::Data& out, escript::Data& in) const;
619 
620  // this is const because setTags is const
621  void updateTagsInUse(int fsType) const;
622 
624  Paso_Pattern* createPasoPattern(const IndexVector& ptr,
625  const IndexVector& index, const dim_t M, const dim_t N) const;
626 
628  Paso_Pattern* createMainPattern() const;
629 
633  void createCouplePatterns(const std::vector<IndexVector>& colIndices,
634  const dim_t N, Paso_Pattern** colPattern,
635  Paso_Pattern** rowPattern) const;
636 
637  void addToSystemMatrix(Paso_SystemMatrix* in, const IndexVector& nodes_Eq,
638  dim_t num_Eq, const IndexVector& nodes_Sol, dim_t num_Sol,
639  const std::vector<double>& array) const;
640 
641  /***********************************************************************/
642 
644  virtual dim_t getNumNodes() const = 0;
645 
647  virtual dim_t getNumElements() const = 0;
648 
650  virtual dim_t getNumDOF() const = 0;
651 
653  virtual dim_t getNumFaceElements() const = 0;
654 
657  virtual dim_t insertNeighbourNodes(IndexVector& index, index_t node) const = 0;
658 
660  virtual void assembleCoordinates(escript::Data& arg) const = 0;
661 
663  virtual void assembleGradient(escript::Data& out, escript::Data& in) const = 0;
664 
666  virtual void assembleIntegrate(std::vector<double>& integrals, escript::Data& arg) const = 0;
667 
670  virtual void assemblePDESingle(Paso_SystemMatrix* mat, escript::Data& rhs,
671  const escript::Data& A, const escript::Data& B,
672  const escript::Data& C, const escript::Data& D,
673  const escript::Data& X, const escript::Data& Y) const = 0;
674 
677  virtual void assemblePDEBoundarySingle(Paso_SystemMatrix* mat,
678  escript::Data& rhs, const escript::Data& d,
679  const escript::Data& y) const = 0;
680 
683  virtual void assemblePDESingleReduced(Paso_SystemMatrix* mat,
684  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
685  const escript::Data& C, const escript::Data& D,
686  const escript::Data& X, const escript::Data& Y) const = 0;
687 
690  virtual void assemblePDEBoundarySingleReduced(Paso_SystemMatrix* mat,
691  escript::Data& rhs, const escript::Data& d,
692  const escript::Data& y) const = 0;
693 
696  virtual void assemblePDESystem(Paso_SystemMatrix* mat, escript::Data& rhs,
697  const escript::Data& A, const escript::Data& B,
698  const escript::Data& C, const escript::Data& D,
699  const escript::Data& X, const escript::Data& Y) const = 0;
700 
703  virtual void assemblePDEBoundarySystem(Paso_SystemMatrix* mat,
704  escript::Data& rhs, const escript::Data& d,
705  const escript::Data& y) const = 0;
706 
709  virtual void assemblePDESystemReduced(Paso_SystemMatrix* mat,
710  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
711  const escript::Data& C, const escript::Data& D,
712  const escript::Data& X, const escript::Data& Y) const = 0;
713 
716  virtual void assemblePDEBoundarySystemReduced(Paso_SystemMatrix* mat,
717  escript::Data& rhs, const escript::Data& d,
718  const escript::Data& y) const = 0;
719 
721  virtual Paso_SystemMatrixPattern* getPattern(bool reducedRowOrder,
722  bool reducedColOrder) const = 0;
723 
725  virtual void interpolateNodesOnElements(escript::Data& out,
726  escript::Data& in, bool reduced) const = 0;
727 
729  virtual void interpolateNodesOnFaces(escript::Data& out, escript::Data& in,
730  bool reduced) const = 0;
731 
733  virtual void nodesToDOF(escript::Data& out, escript::Data& in) const = 0;
734 
736  virtual void dofToNodes(escript::Data& out, escript::Data& in) const = 0;
737 
738 private:
740  void assemblePDE(Paso_SystemMatrix* mat, escript::Data& rhs,
741  const escript::Data& A, const escript::Data& B,
742  const escript::Data& C, const escript::Data& D,
743  const escript::Data& X, const escript::Data& Y) const;
744 
747  void assemblePDEBoundary(Paso_SystemMatrix* mat, escript::Data& rhs,
748  const escript::Data& d, const escript::Data& y) const;
749 };
750 
751 } // end of namespace ripley
752 
753 #endif // __RIPLEY_DOMAIN_H__
754