escript  Revision_
RipleyDomain.h
Go to the documentation of this file.
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