16 #ifndef __RIPLEY_DOMAIN_H__
17 #define __RIPLEY_DOMAIN_H__
19 #include <boost/python/tuple.hpp>
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>
58 virtual int getMPISize()
const {
return m_mpiInfo->size; }
64 virtual int getMPIRank()
const {
return m_mpiInfo->rank; }
72 MPI_Barrier(m_mpiInfo->comm);
93 return m_mpiInfo->comm;
104 virtual bool isValidFunctionSpaceType(
int fsType)
const;
110 virtual std::string functionSpaceTypeAsString(
int fsType)
const;
116 virtual int getDim()
const {
return m_numDim; }
121 virtual bool operator==(
const escript::AbstractDomain& other)
const;
126 virtual bool operator!=(
const escript::AbstractDomain& other)
const {
127 return !(operator==(other));
136 virtual std::pair<int,int> getDataShape(
int fsType)
const;
144 int getTagFromSampleNo(
int fsType,
int sampleNo)
const;
152 virtual void setTagMap(
const std::string& name,
int tag) {
153 m_tagMap[name] = tag;
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;
175 return (m_tagMap.find(name)!=m_tagMap.end());
182 virtual std::string showTagNames()
const;
189 virtual void setNewX(
const escript::Data& arg);
196 virtual void interpolateOnDomain(escript::Data& target,
const escript::Data& source)
const;
203 virtual bool probeInterpolationOnDomain(
int fsType_source,
int fsType_target)
const;
212 virtual signed char preferredInterpolationOnDomain(
int fsType_source,
int fsType_target)
const;
221 commonFunctionSpace(
const std::vector<int>& fs,
int& resultcode)
const;
228 virtual void interpolateACross(escript::Data& target,
const escript::Data& source)
const;
234 virtual bool probeInterpolationACross(
int,
const escript::AbstractDomain&,
int)
const;
240 virtual escript::Data getX()
const;
246 virtual escript::Data getNormal()
const;
251 virtual escript::Data getSize()
const;
258 virtual void setToX(escript::Data& arg)
const;
266 virtual void setToGradient(escript::Data& out,
const escript::Data& in)
const;
273 virtual void setTags(
const int fsType,
const int newTag,
const escript::Data& mask)
const;
280 virtual bool isCellOriented(
int fsType)
const;
288 virtual StatusType
getStatus()
const {
return m_status; }
294 virtual int getNumberOfTagsInUse(
int fsType)
const;
300 virtual const int* borrowListOfTagsInUse(
int fsType)
const;
306 virtual bool canTag(
int fsType)
const;
417 virtual int getSystemMatrixTypeId(
const int solver,
const int preconditioner,
const int package,
const bool symmetry)
const;
428 virtual int getTransportTypeId(
const int solver,
const int preconditioner,
const int package,
const bool symmetry)
const;
435 virtual void setToIntegrals(
DoubleVector& integrals,
const escript::Data& arg)
const;
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;
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;
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;
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;
484 virtual escript::ATP_ptr newTransportProblem(
485 const int blocksize,
const escript::FunctionSpace& functionspace,
486 const int type)
const;
493 virtual void Print_Mesh_Info(
const bool full=
false)
const;
509 virtual std::string getDescription()
const = 0;
516 void dump(
const std::string& filename)
const = 0;
523 const int* borrowSampleReferenceIDs(
int fsType)
const = 0;
531 virtual void setToNormal(escript::Data& out)
const = 0;
538 virtual void setToSize(escript::Data& out)
const = 0;
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;
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;
557 virtual void writeBinaryGrid(
const escript::Data& in, std::string filename,
558 int byteOrder,
int dataType)
const = 0;
564 virtual bool ownSample(
int fsType,
index_t id)
const = 0;
570 virtual int getNumDataPointsGlobal()
const = 0;
576 virtual const int* getNumNodesPerDim()
const = 0;
582 virtual const int* getNumElementsPerDim()
const = 0;
589 virtual const int* getNumFacesPerBoundary()
const = 0;
595 virtual IndexVector getNodeDistribution()
const = 0;
601 virtual const int* getNumSubdivisionsPerDim()
const = 0;
607 virtual double getLocalCoordinate(
int index,
int dim)
const = 0;
613 virtual boost::python::tuple getGridParameters()
const = 0;
619 virtual bool supportsFilter(
const boost::python::tuple& t)
const;
624 virtual escript::Data randomFill(
long seed,
const boost::python::tuple& filter)
const;
637 void copyData(escript::Data& out, escript::Data& in)
const;
640 void averageData(escript::Data& out, escript::Data& in)
const;
643 void multiplyData(escript::Data& out, escript::Data& in)
const;
646 void updateTagsInUse(
int fsType)
const;
658 void createCouplePatterns(
const std::vector<IndexVector>& colIndices,
669 virtual dim_t getNumNodes()
const = 0;
672 virtual dim_t getNumElements()
const = 0;
675 virtual dim_t getNumDOF()
const = 0;
678 virtual dim_t getNumFaceElements()
const = 0;
685 virtual void assembleCoordinates(escript::Data& arg)
const = 0;
688 virtual void assembleGradient(escript::Data& out, escript::Data& in)
const = 0;
691 virtual void assembleIntegrate(
DoubleVector& integrals, escript::Data& arg)
const = 0;
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;
703 escript::Data& rhs,
const escript::Data& d,
704 const escript::Data& y)
const = 0;
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;
716 escript::Data& rhs,
const escript::Data& d,
717 const escript::Data& y)
const = 0;
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;
729 escript::Data& rhs,
const escript::Data& d,
730 const escript::Data& y)
const = 0;
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;
742 escript::Data& rhs,
const escript::Data& d,
743 const escript::Data& y)
const = 0;
747 bool reducedColOrder)
const = 0;
750 virtual void interpolateNodesOnElements(escript::Data& out,
751 escript::Data& in,
bool reduced)
const = 0;
754 virtual void interpolateNodesOnFaces(escript::Data& out, escript::Data& in,
755 bool reduced)
const = 0;
758 virtual void nodesToDOF(escript::Data& out, escript::Data& in)
const = 0;
761 virtual void dofToNodes(escript::Data& out, escript::Data& in)
const = 0;
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;
773 const escript::Data& d,
const escript::Data& y)
const;
778 #endif // __RIPLEY_DOMAIN_H__