escript
Revision_
|
00001 00002 /***************************************************************************** 00003 * 00004 * Copyright (c) 2003-2014 by University of Queensland 00005 * http://www.uq.edu.au 00006 * 00007 * Primary Business: Queensland, Australia 00008 * Licensed under the Open Software License version 3.0 00009 * http://www.opensource.org/licenses/osl-3.0.php 00010 * 00011 * Development until 2012 by Earth Systems Science Computational Center (ESSCC) 00012 * Development 2012-2013 by School of Earth Sciences 00013 * Development from 2014 by Centre for Geoscience Computing (GeoComp) 00014 * 00015 *****************************************************************************/ 00016 00017 #ifndef __RIPLEY_DOMAIN_H__ 00018 #define __RIPLEY_DOMAIN_H__ 00019 00020 #ifdef BADPYTHONMACROS 00021 // This hack is required for BSD/OSX builds with python 2.7 00022 // (and possibly others). It must be the first include. 00023 // From bug reports online it seems that python redefines 00024 // some c macros that are functions in c++. 00025 // c++ doesn't like that! 00026 #include <Python.h> 00027 #undef BADPYTHONMACROS 00028 #endif 00029 00030 #include <boost/python/tuple.hpp> 00031 #include <boost/python/list.hpp> 00032 00033 #include <ripley/Ripley.h> 00034 #include <ripley/RipleyException.h> 00035 #include <ripley/AbstractAssembler.h> 00036 00037 #include <escript/AbstractContinuousDomain.h> 00038 #include <escript/Data.h> 00039 #include <escript/FunctionSpace.h> 00040 00041 #include <paso/SystemMatrix.h> 00042 00043 namespace ripley { 00044 00045 enum assembler_t { 00046 DEFAULT_ASSEMBLER, 00047 WAVE_ASSEMBLER, 00048 LAME_ASSEMBLER 00049 }; 00050 00051 /* There is no particular significance to this type, 00052 It is here as a typedef because a bug in clang++ prevents 00053 that compiler from recognising it as a valid part of 00054 a constant expression. 00055 */ 00056 typedef std::map<std::string, int> simap_t; 00057 00058 00063 struct ReaderParameters 00064 { 00066 std::vector<int> first; 00068 std::vector<int> numValues; 00071 std::vector<int> multiplier; 00073 std::vector<int> reverse; 00075 int byteOrder; 00077 int dataType; 00078 }; 00079 00084 struct DiracPoint 00085 { 00086 int node; 00087 int tag; 00088 }; 00089 00096 class RIPLEY_DLL_API RipleyDomain : public escript::AbstractContinuousDomain 00097 { 00098 public: 00103 RipleyDomain(dim_t dim); 00104 00109 ~RipleyDomain(); 00110 00115 virtual int getMPISize() const { return m_mpiInfo->size; } 00116 00121 virtual int getMPIRank() const { return m_mpiInfo->rank; } 00122 00127 virtual void MPIBarrier() const { 00128 #ifdef ESYS_MPI 00129 MPI_Barrier(m_mpiInfo->comm); 00130 #endif 00131 } 00132 00137 virtual bool onMasterProcessor() const { return getMPIRank()==0; } 00138 00143 #ifdef ESYS_MPI 00144 MPI_Comm 00145 #else 00146 unsigned int 00147 #endif 00148 getMPIComm() const { 00149 #ifdef ESYS_MPI 00150 return m_mpiInfo->comm; 00151 #else 00152 return 0; 00153 #endif 00154 } 00155 00161 virtual bool isValidFunctionSpaceType(int fsType) const; 00162 00167 virtual std::string functionSpaceTypeAsString(int fsType) const; 00168 00173 virtual int getDim() const { return m_numDim; } 00174 00178 virtual bool operator==(const escript::AbstractDomain& other) const; 00179 00183 virtual bool operator!=(const escript::AbstractDomain& other) const { 00184 return !(operator==(other)); 00185 } 00186 00193 virtual std::pair<int,int> getDataShape(int fsType) const; 00194 00201 int getTagFromSampleNo(int fsType, int sampleNo) const; 00202 00209 virtual void setTagMap(const std::string& name, int tag) { 00210 m_tagMap[name] = tag; 00211 } 00212 00218 virtual int getTag(const std::string& name) const { 00219 if (m_tagMap.find(name) != m_tagMap.end()) { 00220 return m_tagMap.find(name)->second; 00221 } else { 00222 throw RipleyException("getTag: invalid tag name"); 00223 } 00224 } 00225 00231 virtual bool isValidTagName(const std::string& name) const { 00232 return (m_tagMap.find(name)!=m_tagMap.end()); 00233 } 00234 00239 virtual std::string showTagNames() const; 00240 00246 virtual void setNewX(const escript::Data& arg); 00247 00253 virtual void interpolateOnDomain(escript::Data& target, const escript::Data& source) const; 00254 00260 virtual bool probeInterpolationOnDomain(int fsType_source, int fsType_target) const; 00261 00269 virtual signed char preferredInterpolationOnDomain(int fsType_source, int fsType_target) const; 00270 00277 bool 00278 commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const; 00279 00285 virtual void interpolateACross(escript::Data& target, const escript::Data& source) const; 00286 00291 virtual bool probeInterpolationACross(int, const escript::AbstractDomain&, int) const; 00292 00297 virtual escript::Data getX() const; 00298 00303 virtual escript::Data getNormal() const; 00304 00308 virtual escript::Data getSize() const; 00309 00315 virtual void setToX(escript::Data& arg) const; 00316 00323 virtual void setToGradient(escript::Data& out, const escript::Data& in) const; 00324 00330 virtual void setTags(const int fsType, const int newTag, const escript::Data& mask) const; 00331 00337 virtual bool isCellOriented(int fsType) const; 00338 00345 virtual StatusType getStatus() const { return m_status; } 00346 00351 virtual int getNumberOfTagsInUse(int fsType) const; 00352 00357 virtual const int* borrowListOfTagsInUse(int fsType) const; 00358 00363 virtual bool canTag(int fsType) const; 00364 00369 virtual int getApproximationOrder(const int fsType) const { return 1; } 00370 00375 virtual bool supportsContactElements() const { return false; } 00376 00381 virtual int getContinuousFunctionCode() const { return Nodes; } 00382 00387 virtual int getReducedContinuousFunctionCode() const { return ReducedNodes; } 00388 00393 virtual int getFunctionCode() const { return Elements; } 00394 00399 virtual int getReducedFunctionCode() const { return ReducedElements; } 00400 00405 virtual int getFunctionOnBoundaryCode() const { return FaceElements; } 00406 00412 virtual int getReducedFunctionOnBoundaryCode() const { return ReducedFaceElements; } 00413 00418 virtual int getFunctionOnContactZeroCode() const { 00419 throw RipleyException("Ripley does not support contact elements"); 00420 } 00421 00426 virtual int getReducedFunctionOnContactZeroCode() const { 00427 throw RipleyException("Ripley does not support contact elements"); 00428 } 00429 00434 virtual int getFunctionOnContactOneCode() const { 00435 throw RipleyException("Ripley does not support contact elements"); 00436 } 00437 00442 virtual int getReducedFunctionOnContactOneCode() const { 00443 throw RipleyException("Ripley does not support contact elements"); 00444 } 00445 00450 virtual int getSolutionCode() const { return DegreesOfFreedom; } 00451 00456 virtual int getReducedSolutionCode() const { return ReducedDegreesOfFreedom; } 00457 00462 virtual int getDiracDeltaFunctionsCode() const { return Points; } 00463 00474 virtual int getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const; 00475 00485 virtual int getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const; 00486 00492 virtual void setToIntegrals(DoubleVector& integrals, const escript::Data& arg) const; 00493 00498 virtual void addPDEToSystem(escript::AbstractSystemMatrix& mat, 00499 escript::Data& rhs, const escript::Data& A, const escript::Data& B, 00500 const escript::Data& C, const escript::Data& D, 00501 const escript::Data& X, const escript::Data& Y, 00502 const escript::Data& d, const escript::Data& y, 00503 const escript::Data& d_contact, const escript::Data& y_contact, 00504 const escript::Data& d_dirac, const escript::Data& y_dirac) const; 00505 00511 virtual void addToSystem(escript::AbstractSystemMatrix& mat, 00512 escript::Data& rhs,std::map<std::string, escript::Data> data) const; 00513 00518 virtual void addToSystemFromPython(escript::AbstractSystemMatrix& mat, 00519 escript::Data& rhs,boost::python::list data) const; 00520 00525 virtual void addPDEToRHS(escript::Data& rhs, const escript::Data& X, 00526 const escript::Data& Y, const escript::Data& y, 00527 const escript::Data& y_contact, const escript::Data& y_dirac) const; 00528 00534 virtual void addToRHS(escript::Data& rhs, 00535 std::map<std::string, escript::Data> data) const; 00536 00541 virtual void addToRHSFromPython(escript::Data& rhs, 00542 boost::python::list data) const; 00543 00548 virtual void addPDEToTransportProblem(escript::AbstractTransportProblem& tp, 00549 escript::Data& source, const escript::Data& M, 00550 const escript::Data& A, const escript::Data& B, 00551 const escript::Data& C, const escript::Data& D, 00552 const escript::Data& X, const escript::Data& Y, 00553 const escript::Data& d, const escript::Data& y, 00554 const escript::Data& d_contact, const escript::Data& y_contact, 00555 const escript::Data& d_dirac, const escript::Data& y_dirac) const; 00556 00557 00562 virtual escript::ASM_ptr newSystemMatrix(const int row_blocksize, 00563 const escript::FunctionSpace& row_functionspace, 00564 const int column_blocksize, 00565 const escript::FunctionSpace& column_functionspace, const int type) const; 00566 00571 virtual escript::ATP_ptr newTransportProblem( 00572 const int blocksize, const escript::FunctionSpace& functionspace, 00573 const int type) const; 00574 00580 virtual void Print_Mesh_Info(const bool full=false) const; 00581 00582 00583 /************************************************************************/ 00584 00590 //void write(const std::string& filename) const = 0; 00591 00596 virtual std::string getDescription() const = 0; 00597 00603 void dump(const std::string& filename) const = 0; 00604 00610 const int* borrowSampleReferenceIDs(int fsType) const = 0; 00611 00618 virtual void setToNormal(escript::Data& out) const = 0; 00619 00625 virtual void setToSize(escript::Data& out) const = 0; 00626 00629 virtual void readNcGrid(escript::Data& out, std::string filename, 00630 std::string varname, const ReaderParameters& params) const = 0; 00631 00634 virtual void readBinaryGrid(escript::Data& out, std::string filename, 00635 const ReaderParameters& params) const = 0; 00636 #ifdef USE_BOOSTIO 00637 virtual void readBinaryGridFromZipped(escript::Data& out, std::string filename, 00638 const ReaderParameters& params) const = 0; 00639 #endif 00640 00643 virtual void writeBinaryGrid(const escript::Data& in, std::string filename, 00644 int byteOrder, int dataType) const = 0; 00645 00650 virtual bool ownSample(int fsType, index_t id) const = 0; 00651 00656 virtual int getNumDataPointsGlobal() const = 0; 00657 00662 virtual const int* getNumNodesPerDim() const = 0; 00663 00668 virtual const int* getNumElementsPerDim() const = 0; 00669 00675 virtual const int* getNumFacesPerBoundary() const = 0; 00676 00681 virtual IndexVector getNodeDistribution() const = 0; 00682 00687 virtual const int* getNumSubdivisionsPerDim() const = 0; 00688 00693 virtual double getLocalCoordinate(int index, int dim) const = 0; 00694 00699 virtual boost::python::tuple getGridParameters() const = 0; 00700 00701 00706 virtual bool supportsFilter(const boost::python::tuple& t) const; 00707 00708 virtual void setAssembler(std::string type, 00709 std::map<std::string, escript::Data> options) { 00710 throw RipleyException("Domain does not support custom assemblers"); 00711 } 00712 00713 void setAssemblerFromPython(std::string type, 00714 boost::python::list options); 00715 protected: 00716 dim_t m_numDim; 00717 StatusType m_status; 00718 Esys_MPIInfo *m_mpiInfo; 00719 TagMap m_tagMap; 00720 mutable IndexVector m_nodeTags, m_nodeTagsInUse; 00721 mutable IndexVector m_elementTags, m_elementTagsInUse; 00722 mutable IndexVector m_faceTags, m_faceTagsInUse; 00723 AbstractAssembler *assembler; 00724 std::vector<struct DiracPoint> m_diracPoints; 00725 IndexVector m_diracPointNodeIDs; //for borrowSampleID 00726 assembler_t assembler_type; 00727 00729 void copyData(escript::Data& out, const escript::Data& in) const; 00730 00732 void averageData(escript::Data& out, const escript::Data& in) const; 00733 00735 void multiplyData(escript::Data& out, const escript::Data& in) const; 00736 00737 // this is const because setTags is const 00738 void updateTagsInUse(int fsType) const; 00739 00741 paso::Pattern_ptr createPasoPattern(const IndexVector& ptr, 00742 const IndexVector& index, const dim_t M, const dim_t N) const; 00743 00745 paso::Pattern_ptr createMainPattern() const; 00746 00750 void createCouplePatterns(const std::vector<IndexVector>& colIndices, 00751 const std::vector<IndexVector>& rowIndices, 00752 const dim_t N, paso::Pattern_ptr& colPattern, 00753 paso::Pattern_ptr& rowPattern) const; 00754 00755 void addToSystemMatrix(paso::SystemMatrix_ptr in, 00756 const IndexVector& nodes_Eq, dim_t num_Eq, 00757 const IndexVector& nodes_Sol, dim_t num_Sol, 00758 const DoubleVector& array) const; 00759 00760 void addPoints(int numPoints, const double* points_ptr, 00761 const int* tags_ptr); 00762 00763 /***********************************************************************/ 00764 00766 virtual dim_t getNumNodes() const = 0; 00767 00769 virtual dim_t getNumElements() const = 0; 00770 00772 virtual dim_t getNumDOF() const = 0; 00773 00775 virtual dim_t getNumFaceElements() const = 0; 00776 00779 virtual dim_t insertNeighbourNodes(IndexVector& index, index_t node) const = 0; 00780 00782 virtual void assembleCoordinates(escript::Data& arg) const = 0; 00783 00785 virtual void assembleGradient(escript::Data& out, const escript::Data& in) const = 0; 00786 00788 virtual void assembleIntegrate(DoubleVector& integrals, const escript::Data& arg) const = 0; 00789 00791 virtual paso::SystemMatrixPattern_ptr getPattern(bool reducedRowOrder, 00792 bool reducedColOrder) const = 0; 00793 00795 virtual void interpolateNodesOnElements(escript::Data& out, 00796 const escript::Data& in, 00797 bool reduced) const = 0; 00798 00800 virtual void interpolateNodesOnFaces(escript::Data& out, 00801 const escript::Data& in, 00802 bool reduced) const = 0; 00803 00805 virtual void nodesToDOF(escript::Data& out, const escript::Data& in) const = 0; 00806 00808 virtual void dofToNodes(escript::Data& out, const escript::Data& in) const = 0; 00809 00810 virtual int getDofOfNode(int node) const = 0; 00811 00812 private: 00814 void assemblePDE(paso::SystemMatrix_ptr mat, escript::Data& rhs, 00815 std::map<std::string, escript::Data> coefs) const; 00816 00819 void assemblePDEBoundary(paso::SystemMatrix_ptr mat, escript::Data& rhs, 00820 std::map<std::string, escript::Data> coefs) const; 00821 00822 void assemblePDEDirac(paso::SystemMatrix_ptr mat, escript::Data& rhs, 00823 std::map<std::string, escript::Data> coefs) const; 00824 00825 // finds the node that the given point belongs to 00826 virtual int findNode(const double *coords) const = 0; 00827 }; 00828 00829 } // end of namespace ripley 00830 00831 #endif // __RIPLEY_DOMAIN_H__ 00832