escript  Revision_
Public Types | Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
dudley::DudleyDomain Class Reference

DudleyDomain implements the AbstractContinuousDomain interface for the Dudley library. More...

#include <DudleyDomain.h>

Inheritance diagram for dudley::DudleyDomain:
escript::AbstractContinuousDomain escript::AbstractDomain

Public Types

typedef std::map< int, std::string > FunctionSpaceNamesMapType
 
- Public Types inherited from escript::AbstractDomain
typedef int StatusType
 

Public Member Functions

 DudleyDomain (const std::string &name, int numDim, escript::JMPI jmpi)
 Constructor for DudleyDomain. More...
 
 DudleyDomain (const DudleyDomain &in)
 Copy constructor. More...
 
 ~DudleyDomain ()
 Destructor for DudleyDomain. More...
 
NodeFilegetNodes () const
 returns a pointer to this domain's node file More...
 
void setElements (ElementFile *elements)
 replaces the element file by elements More...
 
ElementFilegetElements () const
 returns a pointer to this domain's element file More...
 
void setFaceElements (ElementFile *elements)
 replaces the face element file by elements More...
 
ElementFilegetFaceElements () const
 returns a pointer to this domain's face element file More...
 
void setPoints (ElementFile *elements)
 replaces the point element file by elements More...
 
ElementFilegetPoints () const
 returns a pointer to this domain's point (nodal) element file More...
 
virtual escript::JMPI getMPI () const
 returns a reference to the MPI information wrapper for this domain More...
 
virtual int getMPISize () const
 returns the number of processors used for this domain More...
 
virtual int getMPIRank () const
 returns the number MPI rank of this processor More...
 
virtual void MPIBarrier () const
 If compiled for MPI then execute an MPI_Barrier, else do nothing. More...
 
virtual bool onMasterProcessor () const
 returns true if on MPI processor 0, else false More...
 
MPI_Comm getMPIComm () const
 get the communicator for this domain. Returns an integer on non-MPI builds Routine must be implemented by the DomainAdapter. More...
 
void write (const std::string &fileName) const
 writes the current mesh to a file with the given name in the fly file format. More...
 
void Print_Mesh_Info (bool full=false) const
 
void dump (const std::string &fileName) const
 dumps the mesh to a file with the given name. More...
 
int getTagFromSampleNo (int functionSpaceType, index_t sampleNo) const
 Return the tag key for the given sample number. More...
 
const index_t * borrowSampleReferenceIDs (int functionSpaceType) const
 Return the reference number of the given sample number. More...
 
virtual bool isValidFunctionSpaceType (int functionSpaceType) const
 Returns true if the given integer is a valid function space type for this domain. More...
 
virtual std::string getDescription () const
 Return a description for this domain. More...
 
virtual std::string functionSpaceTypeAsString (int functionSpaceType) const
 Return a description for the given function space type code. More...
 
void setFunctionSpaceTypeNames ()
 Build the table of function space type names. More...
 
virtual int getContinuousFunctionCode () const
 Return a continuous FunctionSpace code. More...
 
virtual int getReducedContinuousFunctionCode () const
 Return a continuous on reduced order nodes FunctionSpace code. More...
 
virtual int getFunctionCode () const
 Return a function FunctionSpace code. More...
 
virtual int getReducedFunctionCode () const
 Return a function with reduced integration order FunctionSpace code. More...
 
virtual int getFunctionOnBoundaryCode () const
 Return a function on boundary FunctionSpace code. More...
 
virtual int getReducedFunctionOnBoundaryCode () const
 Return a function on boundary with reduced integration order FunctionSpace code. More...
 
virtual int getFunctionOnContactZeroCode () const
 Return a FunctionOnContactZero code. More...
 
virtual int getReducedFunctionOnContactZeroCode () const
 Return a FunctionOnContactZero code with reduced integration order. More...
 
virtual int getFunctionOnContactOneCode () const
 Return a FunctionOnContactOne code. More...
 
virtual int getReducedFunctionOnContactOneCode () const
 Return a FunctionOnContactOne code with reduced integration order. More...
 
virtual int getSolutionCode () const
 Return a Solution code. More...
 
virtual int getReducedSolutionCode () const
 Return a ReducedSolution code. More...
 
virtual int getDiracDeltaFunctionsCode () const
 Return a DiracDeltaFunctions code. More...
 
virtual int getDim () const
 Returns the spatial dimension of the domain. More...
 
virtual StatusType getStatus () const
 Returns a status indicator of the domain. The status identifier should be unique over the live time if the object but may be updated if changes to the domain happen, e.g. modifications to its geometry. More...
 
virtual dim_t getNumDataPointsGlobal () const
 Return the number of data points summed across all MPI processes. More...
 
virtual std::pair< int, dim_t > getDataShape (int functionSpaceCode) const
 Return the number of data points per sample, and the number of samples as a pair. More...
 
virtual void setToX (escript::Data &arg) const
 copies the location of data points into arg. The domain of arg has to match this. has to be implemented by the actual Domain adapter. More...
 
virtual void setTagMap (const std::string &name, int tag)
 sets a map from a clear tag name to a tag key More...
 
virtual int getTag (const std::string &name) const
 Return the tag key for tag name. More...
 
virtual bool isValidTagName (const std::string &name) const
 Returns true if name is a defined tag name. More...
 
virtual std::string showTagNames () const
 Returns all tag names in a single string sperated by commas. More...
 
virtual void setNewX (const escript::Data &arg)
 assigns new location to the domain More...
 
virtual void interpolateOnDomain (escript::Data &target, const escript::Data &source) const
 interpolates data given on source onto target where source and target have to be given on the same domain. More...
 
virtual bool probeInterpolationOnDomain (int functionSpaceType_source, int functionSpaceType_target) const
 True if interpolation is possible from source to target. More...
 
virtual signed char preferredInterpolationOnDomain (int functionSpaceType_source, int functionSpaceType_target) const
 Preferred direction of interpolation. If you really need to test for a particular direction, then use probeInterpolation. More...
 
bool commonFunctionSpace (const std::vector< int > &fs, int &resultcode) const
 given a vector of FunctionSpace typecodes, pass back a code which then can all be interpolated to. More...
 
virtual void interpolateAcross (escript::Data &target, const escript::Data &source) const
 interpolates data given on source onto target where source and target are given on different domains. More...
 
virtual bool probeInterpolationAcross (int functionSpaceType_source, const escript::AbstractDomain &targetDomain, int functionSpaceType_target) const
 determines whether interpolation from source to target is possible. More...
 
virtual void setToNormal (escript::Data &out) const
 copies the surface normals at data points into out. The actual function space to be considered is defined by out. out has to be defined on this. More...
 
virtual void setToSize (escript::Data &out) const
 copies the size of samples into out. The actual function space to be considered is defined by out. out has to be defined on this. More...
 
virtual void setToGradient (escript::Data &grad, const escript::Data &arg) const
 copies the gradient of arg into grad. The actual function space to be considered for the gradient is defined by grad. arg and grad have to be defined on this. More...
 
virtual void setToIntegrals (std::vector< double > &integrals, const escript::Data &arg) const
 copies the integrals of the function defined by arg into integrals. arg has to be defined on this. More...
 
virtual int getSystemMatrixTypeId (const boost::python::object &options) const
 return the identifier of the matrix type to be used for the global stiffness matrix when a particular solver, package, preconditioner, and symmetric matrix is used. More...
 
virtual int getTransportTypeId (int solver, int preconditioner, int package, bool symmetry) const
 return the identifier of the transport problem type to be used when a particular solver, perconditioner, package and symmetric matrix is used. More...
 
virtual bool isCellOriented (int functionSpaceCode) const
 returns true if data on this domain and a function space of type functionSpaceCode has to considered as cell centered data. More...
 
virtual bool ownSample (int fsCode, index_t id) const
 
virtual void addPDEToSystem (escript::AbstractSystemMatrix &mat, escript::Data &rhs, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y, const escript::Data &d, const escript::Data &y, const escript::Data &d_contact, const escript::Data &y_contact, const escript::Data &d_dirac, const escript::Data &y_dirac) const
 adds a PDE onto the stiffness matrix mat and a rhs More...
 
virtual void addPDEToLumpedSystem (escript::Data &mat, const escript::Data &D, const escript::Data &d, const escript::Data &d_dirac, bool useHRZ) const
 adds a PDE onto the lumped stiffness matrix matrix More...
 
virtual void addPDEToRHS (escript::Data &rhs, const escript::Data &X, const escript::Data &Y, const escript::Data &y, const escript::Data &y_contact, const escript::Data &y_dirac) const
 adds a PDE onto the stiffness matrix mat and a rhs More...
 
virtual void addPDEToTransportProblem (escript::AbstractTransportProblem &tp, escript::Data &source, const escript::Data &M, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y, const escript::Data &d, const escript::Data &y, const escript::Data &d_contact, const escript::Data &y_contact, const escript::Data &d_dirac, const escript::Data &y_dirac) const
 adds a PDE onto a transport problem More...
 
escript::ASM_ptr newSystemMatrix (int row_blocksize, const escript::FunctionSpace &row_functionspace, int column_blocksize, const escript::FunctionSpace &column_functionspace, int type) const
 creates a stiffness matrix and initializes it with zeros More...
 
escript::ATP_ptr newTransportProblem (int blocksize, const escript::FunctionSpace &functionspace, int type) const
 creates a TransportProblem More...
 
virtual escript::Data getX () const
 returns locations in the FEM nodes More...
 
virtual escript::Data getNormal () const
 returns boundary normals at the quadrature point on the face elements More...
 
virtual escript::Data getSize () const
 returns the element size More...
 
virtual bool operator== (const escript::AbstractDomain &other) const
 comparison operators More...
 
virtual bool operator!= (const escript::AbstractDomain &other) const
 Return true if given domains are not equal. More...
 
virtual void setTags (int functionSpaceType, int newTag, const escript::Data &mask) const
 assigns new tag newTag to all samples of functionspace with a positive value of mask for any its sample point. More...
 
virtual int getNumberOfTagsInUse (int functionSpaceCode) const
 returns the number of tags in use and a pointer to an array with the number of tags in use More...
 
virtual const int * borrowListOfTagsInUse (int functionSpaceCode) const
 
virtual bool canTag (int functionSpaceCode) const
 Checks if this domain allows tags for the specified functionSpace code. More...
 
virtual int getApproximationOrder (int functionSpaceCode) const
 returns the approximation order used for a function space functionSpaceCode More...
 
virtual bool supportsContactElements () const
 
virtual escript::Data randomFill (const escript::DataTypes::ShapeType &shape, const escript::FunctionSpace &what, long seed, const boost::python::tuple &filter) const
 Fills the data object with filtered random values. More...
 
void createMappings (const std::vector< index_t > &dofDistribution, const std::vector< index_t > &nodeDistribution)
 
void relabelElementNodes (const index_t *newNode, index_t offset)
 
- Public Member Functions inherited from escript::AbstractContinuousDomain
 AbstractContinuousDomain ()
 Default constructor for AbstractContinuousDomain. More...
 
virtual ~AbstractContinuousDomain ()
 Destructor for AbstractContinuousDomain. More...
 
- Public Member Functions inherited from escript::AbstractDomain
Domain_ptr getPtr ()
 Returns smart pointer which is managing this object. If one does not exist yet it creates one. More...
 
const_Domain_ptr getPtr () const
 
virtual ~AbstractDomain ()
 Destructor for AbstractDomain. More...
 
virtual int getTagFromSampleNo (int functionSpaceType, DataTypes::index_t sampleNo) const =0
 Return the tag key for the given sample number. More...
 
virtual bool ownSample (int fs_code, DataTypes::index_t id) const =0
 True if this rank owns the sample(id) Must be implemented by the Domain adapter. More...
 
void throwStandardException (const std::string &functionName) const
 Throw a standard exception. This function is called if any attempt is made to use a base class function. More...
 
virtual bool supportsFilter (const boost::python::tuple &t) const
 true if this domain can handle to specified tuple of filter options. More...
 

Static Public Member Functions

static escript::Domain_ptr load (const std::string &filename)
 recovers domain from a dump file More...
 
static escript::Domain_ptr read (escript::JMPI mpiInfo, const std::string &filename, bool optimize)
 reads a mesh from a fly file. For MPI parallel runs fans out the mesh to multiple processes. More...
 
static escript::Domain_ptr readGmsh (escript::JMPI mpiInfo, const std::string &filename, int numDim, bool optimize)
 reads a gmsh mesh file. More...
 
static escript::Domain_ptr create2D (dim_t NE0, dim_t NE1, double l0, double l1, bool optimize, escript::JMPI jmpi)
 Creates a 2-dimensional rectangular domain. More...
 
static escript::Domain_ptr create3D (dim_t NE0, dim_t NE1, dim_t NE2, double l0, double l1, double l2, bool optimize, escript::JMPI jmpi)
 Creates a 3-dimensional rectangular domain. More...
 

Private Member Functions

void prepare (bool optimize)
 prepares the mesh for further use More...
 
void resolveNodeIds ()
 
void createColoring (const index_t *dofMap)
 tries to reduce the number of colours for all element files More...
 
void distributeByRankOfDOF (const IndexVector &distribution)
 
void markNodes (std::vector< short > &mask, index_t offset) const
 
void optimizeDOFDistribution (IndexVector &distribution)
 
void optimizeDOFLabeling (const IndexVector &distribution)
 optimizes the labeling of the DOFs on each processor More...
 
void optimizeElementOrdering ()
 redistributes elements to minimize communication during assemblage More...
 
void updateTagList ()
 regenerates list of tags in use for node file and element files More...
 
void printElementInfo (const ElementFile *e, const std::string &title, const std::string &defaultType, bool full) const
 
void writeElementInfo (std::ostream &stream, const ElementFile *e, const std::string &defaultType) const
 

Private Attributes

escript::JMPI m_mpiInfo
 MPI information. More...
 
std::string m_name
 domain description More...
 
NodeFilem_nodes
 the table of the nodes More...
 
ElementFilem_elements
 the table of the elements More...
 
ElementFilem_faceElements
 the table of face elements More...
 
ElementFilem_points
 the table of points (treated as elements of dimension 0) More...
 
TagMap m_tagMap
 the tag map mapping names to tag keys More...
 

Static Private Attributes

static FunctionSpaceNamesMapType m_functionSpaceTypeNames
 

Detailed Description

DudleyDomain implements the AbstractContinuousDomain interface for the Dudley library.

Member Typedef Documentation

◆ FunctionSpaceNamesMapType

typedef std::map<int, std::string> dudley::DudleyDomain::FunctionSpaceNamesMapType

Constructor & Destructor Documentation

◆ DudleyDomain() [1/2]

dudley::DudleyDomain::DudleyDomain ( const std::string &  name,
int  numDim,
escript::JMPI  jmpi 
)

Constructor for DudleyDomain.

Parameters
namea descriptive name for the domain
numDimdimensionality of the domain (2 or 3)
jmpishared pointer to MPI Information to be used

References m_mpiInfo, m_nodes, and setFunctionSpaceTypeNames().

Referenced by create2D(), create3D(), read(), and readGmsh().

◆ DudleyDomain() [2/2]

dudley::DudleyDomain::DudleyDomain ( const DudleyDomain in)

Copy constructor.

References setFunctionSpaceTypeNames().

◆ ~DudleyDomain()

dudley::DudleyDomain::~DudleyDomain ( )

Destructor for DudleyDomain.

References m_elements, m_faceElements, m_nodes, and m_points.

Member Function Documentation

◆ addPDEToLumpedSystem()

void dudley::DudleyDomain::addPDEToLumpedSystem ( escript::Data mat,
const escript::Data D,
const escript::Data d,
const escript::Data d_dirac,
bool  useHRZ 
) const
virtual

adds a PDE onto the lumped stiffness matrix matrix

References dudley::Assemble_LumpedSystem(), m_elements, m_faceElements, m_nodes, and m_points.

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ addPDEToRHS()

void dudley::DudleyDomain::addPDEToRHS ( escript::Data rhs,
const escript::Data X,
const escript::Data Y,
const escript::Data y,
const escript::Data y_contact,
const escript::Data y_dirac 
) const
virtual

adds a PDE onto the stiffness matrix mat and a rhs

Reimplemented from escript::AbstractContinuousDomain.

References dudley::Assemble_PDE(), escript::Data::isEmpty(), m_elements, m_faceElements, m_nodes, and m_points.

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ addPDEToSystem()

void dudley::DudleyDomain::addPDEToSystem ( escript::AbstractSystemMatrix mat,
escript::Data rhs,
const escript::Data A,
const escript::Data B,
const escript::Data C,
const escript::Data D,
const escript::Data X,
const escript::Data Y,
const escript::Data d,
const escript::Data y,
const escript::Data d_contact,
const escript::Data y_contact,
const escript::Data d_dirac,
const escript::Data y_dirac 
) const
virtual

adds a PDE onto the stiffness matrix mat and a rhs

Reimplemented from escript::AbstractContinuousDomain.

References dudley::Assemble_PDE(), escript::AbstractSystemMatrix::getPtr(), escript::Data::isEmpty(), m_elements, m_faceElements, m_nodes, and m_points.

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ addPDEToTransportProblem()

void dudley::DudleyDomain::addPDEToTransportProblem ( escript::AbstractTransportProblem tp,
escript::Data source,
const escript::Data M,
const escript::Data A,
const escript::Data B,
const escript::Data C,
const escript::Data D,
const escript::Data X,
const escript::Data Y,
const escript::Data d,
const escript::Data y,
const escript::Data d_contact,
const escript::Data y_contact,
const escript::Data d_dirac,
const escript::Data y_dirac 
) const
virtual

◆ borrowListOfTagsInUse()

const int * dudley::DudleyDomain::borrowListOfTagsInUse ( int  functionSpaceCode) const
virtual

◆ borrowSampleReferenceIDs()

const index_t * dudley::DudleyDomain::borrowSampleReferenceIDs ( int  functionSpaceType) const
virtual

◆ canTag()

bool dudley::DudleyDomain::canTag ( int  functionSpaceCode) const
virtual

Checks if this domain allows tags for the specified functionSpace code.

Implements escript::AbstractDomain.

References dudley::Elements, dudley::FaceElements, dudley::Nodes, dudley::Points, dudley::ReducedElements, and dudley::ReducedFaceElements.

Referenced by getDim().

◆ commonFunctionSpace()

bool dudley::DudleyDomain::commonFunctionSpace ( const std::vector< int > &  fs,
int &  resultcode 
) const
virtual

given a vector of FunctionSpace typecodes, pass back a code which then can all be interpolated to.

Returns
true is result is valid, false if not

Implements escript::AbstractDomain.

References dudley::DegreesOfFreedom, dudley::Elements, dudley::FaceElements, dudley::Nodes, dudley::Points, dudley::ReducedElements, and dudley::ReducedFaceElements.

Referenced by getDim().

◆ create2D()

escript::Domain_ptr dudley::DudleyDomain::create2D ( dim_t  NE0,
dim_t  NE1,
double  l0,
double  l1,
bool  optimize,
escript::JMPI  jmpi 
)
static

Creates a 2-dimensional rectangular domain.

Parameters
NE0Input - number of elements in first dimension
NE1Input - number of elements in second dimension
l0Input - length of domain in first dimension (width)
l1Input - length of domain in second dimension (height)
optimizeInput - whether to optimize node/DOF labelling
jmpiInput - Shared pointer to MPI Information to be used

References dudley::NodeFile::allocTable(), dudley::ElementFile::allocTable(), dudley::NodeFile::Coordinates, DIM, dudley::Dudley_Line2, dudley::Dudley_Point1, dudley::Dudley_Tri3, DudleyDomain(), getNodes(), escript::AbstractDomain::getPtr(), dudley::NodeFile::globalDegreesOfFreedom, dudley::ElementFile::Id, dudley::NodeFile::Id, INDEX2, dudley::ElementFile::Nodes, dudley::ElementFile::numNodes, dudley::ElementFile::Owner, prepare(), resolveNodeIds(), setElements(), setFaceElements(), setPoints(), setTagMap(), dudley::ElementFile::Tag, and dudley::NodeFile::Tag.

◆ create3D()

escript::Domain_ptr dudley::DudleyDomain::create3D ( dim_t  NE0,
dim_t  NE1,
dim_t  NE2,
double  l0,
double  l1,
double  l2,
bool  optimize,
escript::JMPI  jmpi 
)
static

Creates a 3-dimensional rectangular domain.

Parameters
NE0Input - number of elements in first dimension
NE1Input - number of elements in second dimension
NE2Input - number of elements in third dimension
l0Input - length of domain in first dimension (width)
l1Input - length of domain in second dimension (height)
l2Input - length of domain in third dimension (depth)
optimizeInput - whether to optimize node/DOF labelling
jmpiInput - Shared pointer to MPI Information to be used

References dudley::ElementFile::allocTable(), DIM, dudley::Dudley_Point1, dudley::Dudley_Tet4, dudley::Dudley_Tri3, DudleyDomain(), getNodes(), escript::AbstractDomain::getPtr(), dudley::ElementFile::Id, INDEX2, MAX3, dudley::ElementFile::Nodes, dudley::ElementFile::numNodes, dudley::ElementFile::Owner, prepare(), resolveNodeIds(), setElements(), setFaceElements(), setPoints(), setTagMap(), and dudley::ElementFile::Tag.

◆ createColoring()

void dudley::DudleyDomain::createColoring ( const index_t *  dofMap)
private

tries to reduce the number of colours for all element files

References dudley::ElementFile::createColoring(), dudley::NodeFile::getNumNodes(), m_elements, m_faceElements, m_nodes, and m_points.

Referenced by distributeByRankOfDOF(), and supportsContactElements().

◆ createMappings()

void dudley::DudleyDomain::createMappings ( const std::vector< index_t > &  dofDistribution,
const std::vector< index_t > &  nodeDistribution 
)

◆ distributeByRankOfDOF()

void dudley::DudleyDomain::distributeByRankOfDOF ( const IndexVector &  distribution)
private

redistributes the Nodes and Elements including overlap according to the DOF distribution. It will create an element colouring but will not create any mappings.

References dudley::NodeFile::assignMPIRankToDOFs(), createColoring(), dudley::ElementFile::distributeByRankOfDOF(), ESYS_ASSERT, dudley::NodeFile::getDOFRange(), dudley::NodeFile::getNumNodes(), dudley::NodeFile::globalDegreesOfFreedom, dudley::NodeFile::Id, m_elements, m_faceElements, m_nodes, m_points, and resolveNodeIds().

Referenced by prepare(), and supportsContactElements().

◆ dump()

void dudley::DudleyDomain::dump ( const std::string &  fileName) const
virtual

◆ functionSpaceTypeAsString()

string dudley::DudleyDomain::functionSpaceTypeAsString ( int  functionSpaceType) const
virtual

Return a description for the given function space type code.

Implements escript::AbstractDomain.

References m_functionSpaceTypeNames.

Referenced by getMPIComm().

◆ getApproximationOrder()

int dudley::DudleyDomain::getApproximationOrder ( int  functionSpaceCode) const
virtual

returns the approximation order used for a function space functionSpaceCode

Implements escript::AbstractDomain.

References dudley::DegreesOfFreedom, dudley::Elements, dudley::FaceElements, dudley::Nodes, dudley::Points, dudley::ReducedElements, and dudley::ReducedFaceElements.

Referenced by getDim().

◆ getContinuousFunctionCode()

int dudley::DudleyDomain::getContinuousFunctionCode ( ) const
virtual

Return a continuous FunctionSpace code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::Nodes.

Referenced by getMPIComm().

◆ getDataShape()

pair< int, dim_t > dudley::DudleyDomain::getDataShape ( int  functionSpaceCode) const
virtual

◆ getDescription()

string dudley::DudleyDomain::getDescription ( ) const
virtual

Return a description for this domain.

Reimplemented from escript::AbstractContinuousDomain.

Referenced by BOOST_PYTHON_MODULE(), borrowSampleReferenceIDs(), getDataShape(), getMPIComm(), and getTagFromSampleNo().

◆ getDim()

virtual int dudley::DudleyDomain::getDim ( ) const
inlinevirtual

◆ getDiracDeltaFunctionsCode()

int dudley::DudleyDomain::getDiracDeltaFunctionsCode ( ) const
virtual

Return a DiracDeltaFunctions code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::Points.

Referenced by getMPIComm().

◆ getElements()

ElementFile* dudley::DudleyDomain::getElements ( ) const
inline

returns a pointer to this domain's element file

References m_elements, and setFaceElements().

Referenced by weipa::FinleyDomain::initFromEscript().

◆ getFaceElements()

ElementFile* dudley::DudleyDomain::getFaceElements ( ) const
inline

returns a pointer to this domain's face element file

References m_faceElements, and setPoints().

Referenced by weipa::FinleyDomain::initFromEscript().

◆ getFunctionCode()

int dudley::DudleyDomain::getFunctionCode ( ) const
virtual

Return a function FunctionSpace code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::Elements.

Referenced by getMPIComm().

◆ getFunctionOnBoundaryCode()

int dudley::DudleyDomain::getFunctionOnBoundaryCode ( ) const
virtual

Return a function on boundary FunctionSpace code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::FaceElements.

Referenced by getMPIComm().

◆ getFunctionOnContactOneCode()

int dudley::DudleyDomain::getFunctionOnContactOneCode ( ) const
virtual

Return a FunctionOnContactOne code.

Reimplemented from escript::AbstractContinuousDomain.

Referenced by getMPIComm().

◆ getFunctionOnContactZeroCode()

int dudley::DudleyDomain::getFunctionOnContactZeroCode ( ) const
virtual

Return a FunctionOnContactZero code.

Reimplemented from escript::AbstractContinuousDomain.

Referenced by getMPIComm().

◆ getMPI()

virtual escript::JMPI dudley::DudleyDomain::getMPI ( ) const
inlinevirtual

returns a reference to the MPI information wrapper for this domain

Implements escript::AbstractDomain.

References m_mpiInfo.

Referenced by getTransportTypeId().

◆ getMPIComm()

MPI_Comm dudley::DudleyDomain::getMPIComm ( ) const
inlinevirtual

◆ getMPIRank()

virtual int dudley::DudleyDomain::getMPIRank ( ) const
inlinevirtual

returns the number MPI rank of this processor

Implements escript::AbstractDomain.

References m_mpiInfo, and MPIBarrier().

Referenced by BOOST_PYTHON_MODULE(), dump(), and onMasterProcessor().

◆ getMPISize()

virtual int dudley::DudleyDomain::getMPISize ( ) const
inlinevirtual

returns the number of processors used for this domain

Implements escript::AbstractDomain.

References m_mpiInfo.

Referenced by BOOST_PYTHON_MODULE(), dump(), getSystemMatrixTypeId(), interpolateOnDomain(), ownSample(), and setToGradient().

◆ getNodes()

NodeFile* dudley::DudleyDomain::getNodes ( ) const
inline

returns a pointer to this domain's node file

References m_nodes, and setElements().

Referenced by create2D(), create3D(), weipa::FinleyDomain::initFromEscript(), load(), read(), and readGmsh().

◆ getNormal()

escript::Data dudley::DudleyDomain::getNormal ( ) const
virtual

returns boundary normals at the quadrature point on the face elements

Implements escript::AbstractDomain.

References escript::functionOnBoundary(), and escript::FunctionSpace::getNormal().

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ getNumberOfTagsInUse()

int dudley::DudleyDomain::getNumberOfTagsInUse ( int  functionSpaceCode) const
virtual

◆ getNumDataPointsGlobal()

dim_t dudley::DudleyDomain::getNumDataPointsGlobal ( ) const
virtual

Return the number of data points summed across all MPI processes.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::NodeFile::getGlobalNumNodes(), and m_nodes.

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ getPoints()

ElementFile* dudley::DudleyDomain::getPoints ( ) const
inline

returns a pointer to this domain's point (nodal) element file

References m_points.

◆ getReducedContinuousFunctionCode()

int dudley::DudleyDomain::getReducedContinuousFunctionCode ( ) const
virtual

Return a continuous on reduced order nodes FunctionSpace code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::Nodes.

Referenced by getMPIComm().

◆ getReducedFunctionCode()

int dudley::DudleyDomain::getReducedFunctionCode ( ) const
virtual

Return a function with reduced integration order FunctionSpace code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::ReducedElements.

Referenced by getMPIComm().

◆ getReducedFunctionOnBoundaryCode()

int dudley::DudleyDomain::getReducedFunctionOnBoundaryCode ( ) const
virtual

Return a function on boundary with reduced integration order FunctionSpace code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::ReducedFaceElements.

Referenced by getMPIComm().

◆ getReducedFunctionOnContactOneCode()

int dudley::DudleyDomain::getReducedFunctionOnContactOneCode ( ) const
virtual

Return a FunctionOnContactOne code with reduced integration order.

Reimplemented from escript::AbstractContinuousDomain.

Referenced by getMPIComm().

◆ getReducedFunctionOnContactZeroCode()

int dudley::DudleyDomain::getReducedFunctionOnContactZeroCode ( ) const
virtual

Return a FunctionOnContactZero code with reduced integration order.

Reimplemented from escript::AbstractContinuousDomain.

Referenced by getMPIComm().

◆ getReducedSolutionCode()

int dudley::DudleyDomain::getReducedSolutionCode ( ) const
virtual

Return a ReducedSolution code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::DegreesOfFreedom.

Referenced by getMPIComm().

◆ getSize()

escript::Data dudley::DudleyDomain::getSize ( ) const
virtual

returns the element size

Implements escript::AbstractDomain.

References escript::function(), and escript::FunctionSpace::getSize().

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ getSolutionCode()

int dudley::DudleyDomain::getSolutionCode ( ) const
virtual

Return a Solution code.

Reimplemented from escript::AbstractContinuousDomain.

References dudley::DegreesOfFreedom.

Referenced by getMPIComm().

◆ getStatus()

DudleyDomain::StatusType dudley::DudleyDomain::getStatus ( ) const
virtual

Returns a status indicator of the domain. The status identifier should be unique over the live time if the object but may be updated if changes to the domain happen, e.g. modifications to its geometry.

Reimplemented from escript::AbstractDomain.

References m_nodes, and dudley::NodeFile::status.

Referenced by getDim().

◆ getSystemMatrixTypeId()

int dudley::DudleyDomain::getSystemMatrixTypeId ( const boost::python::object &  options) const
virtual

◆ getTag()

int dudley::DudleyDomain::getTag ( const std::string &  name) const
virtual

Return the tag key for tag name.

Parameters
nameInput - tag name

Implements escript::AbstractDomain.

References m_tagMap.

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ getTagFromSampleNo()

int dudley::DudleyDomain::getTagFromSampleNo ( int  functionSpaceType,
index_t  sampleNo 
) const

Return the tag key for the given sample number.

Parameters
functionSpaceTypeInput - The function space type.
sampleNoInput - The sample number.

References dudley::DegreesOfFreedom, dudley::Elements, dudley::FaceElements, getDescription(), m_elements, m_faceElements, m_nodes, m_points, dudley::Nodes, dudley::Points, dudley::ReducedElements, dudley::ReducedFaceElements, dudley::ElementFile::Tag, and dudley::NodeFile::Tag.

Referenced by getMPIComm().

◆ getTransportTypeId()

int dudley::DudleyDomain::getTransportTypeId ( int  solver,
int  preconditioner,
int  package,
bool  symmetry 
) const
virtual

return the identifier of the transport problem type to be used when a particular solver, perconditioner, package and symmetric matrix is used.

Parameters
solver
preconditioner
package
symmetry

Reimplemented from escript::AbstractContinuousDomain.

References getMPI(), and paso::TransportProblem::getTypeId().

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ getX()

escript::Data dudley::DudleyDomain::getX ( ) const
virtual

returns locations in the FEM nodes

Implements escript::AbstractDomain.

References escript::continuousFunction(), and escript::FunctionSpace::getX().

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ interpolateAcross()

void dudley::DudleyDomain::interpolateAcross ( escript::Data target,
const escript::Data source 
) const
virtual

interpolates data given on source onto target where source and target are given on different domains.

Implements escript::AbstractDomain.

Referenced by getDim().

◆ interpolateOnDomain()

void dudley::DudleyDomain::interpolateOnDomain ( escript::Data target,
const escript::Data source 
) const
virtual

◆ isCellOriented()

bool dudley::DudleyDomain::isCellOriented ( int  functionSpaceCode) const
virtual

returns true if data on this domain and a function space of type functionSpaceCode has to considered as cell centered data.

Implements escript::AbstractDomain.

References dudley::DegreesOfFreedom, dudley::Elements, dudley::FaceElements, dudley::Nodes, dudley::Points, dudley::ReducedElements, and dudley::ReducedFaceElements.

Referenced by getDim().

◆ isValidFunctionSpaceType()

bool dudley::DudleyDomain::isValidFunctionSpaceType ( int  functionSpaceType) const
virtual

Returns true if the given integer is a valid function space type for this domain.

Reimplemented from escript::AbstractContinuousDomain.

References m_functionSpaceTypeNames.

Referenced by getMPIComm().

◆ isValidTagName()

bool dudley::DudleyDomain::isValidTagName ( const std::string &  name) const
virtual

Returns true if name is a defined tag name.

Parameters
nameInput - tag name to be checked.

Reimplemented from escript::AbstractDomain.

References m_tagMap.

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ load()

Domain_ptr dudley::DudleyDomain::load ( const std::string &  filename)
static

◆ markNodes()

void dudley::DudleyDomain::markNodes ( std::vector< short > &  mask,
index_t  offset 
) const
private

◆ MPIBarrier()

void dudley::DudleyDomain::MPIBarrier ( ) const
virtual

If compiled for MPI then execute an MPI_Barrier, else do nothing.

Implements escript::AbstractDomain.

References getMPIComm().

Referenced by BOOST_PYTHON_MODULE(), and getMPIRank().

◆ newSystemMatrix()

escript::ASM_ptr dudley::DudleyDomain::newSystemMatrix ( int  row_blocksize,
const escript::FunctionSpace row_functionspace,
int  column_blocksize,
const escript::FunctionSpace column_functionspace,
int  type 
) const
virtual

◆ newTransportProblem()

escript::ATP_ptr dudley::DudleyDomain::newTransportProblem ( int  blocksize,
const escript::FunctionSpace functionspace,
int  type 
) const
virtual

◆ onMasterProcessor()

virtual bool dudley::DudleyDomain::onMasterProcessor ( ) const
inlinevirtual

returns true if on MPI processor 0, else false

Implements escript::AbstractDomain.

References getMPIRank().

Referenced by BOOST_PYTHON_MODULE().

◆ operator!=()

bool dudley::DudleyDomain::operator!= ( const escript::AbstractDomain other) const
virtual

Return true if given domains are not equal.

Implements escript::AbstractDomain.

References operator==().

Referenced by getDim().

◆ operator==()

bool dudley::DudleyDomain::operator== ( const escript::AbstractDomain other) const
virtual

comparison operators

Implements escript::AbstractDomain.

References m_elements, m_faceElements, m_nodes, and m_points.

Referenced by getDim(), and operator!=().

◆ optimizeDOFDistribution()

void dudley::DudleyDomain::optimizeDOFDistribution ( IndexVector &  distribution)
private

optimizes the distribution of DOFs across processors using ParMETIS. On return a new distribution is given and the globalDOF are relabeled accordingly but the mesh has not been redistributed yet

References dudley::NodeFile::Coordinates, dudley::NodeFile::getNumNodes(), dudley::NodeFile::globalDegreesOfFreedom, INDEX2, dudley::IndexList_insertElementsWithRowRangeNoMainDiagonal(), m_elements, m_faceElements, m_mpiInfo, m_nodes, m_points, and dudley::NodeFile::numDim.

Referenced by prepare(), and supportsContactElements().

◆ optimizeDOFLabeling()

void dudley::DudleyDomain::optimizeDOFLabeling ( const IndexVector &  distribution)
private

◆ optimizeElementOrdering()

void dudley::DudleyDomain::optimizeElementOrdering ( )
private

redistributes elements to minimize communication during assemblage

References m_elements, m_faceElements, m_points, and dudley::ElementFile::optimizeOrdering().

Referenced by prepare(), and supportsContactElements().

◆ ownSample()

bool dudley::DudleyDomain::ownSample ( int  fsCode,
index_t  id 
) const
virtual

◆ preferredInterpolationOnDomain()

signed char dudley::DudleyDomain::preferredInterpolationOnDomain ( int  functionSpaceType_source,
int  functionSpaceType_target 
) const
virtual

Preferred direction of interpolation. If you really need to test for a particular direction, then use probeInterpolation.

Returns
0 for not possible, 1 for possible and preferred, -1 other direction preferred (does not mean this direction is possible)

Implements escript::AbstractDomain.

References probeInterpolationOnDomain().

Referenced by getDim().

◆ prepare()

void dudley::DudleyDomain::prepare ( bool  optimize)
private

◆ Print_Mesh_Info()

void dudley::DudleyDomain::Print_Mesh_Info ( bool  full = false) const
virtual

◆ printElementInfo()

void dudley::DudleyDomain::printElementInfo ( const ElementFile e,
const std::string &  title,
const std::string &  defaultType,
bool  full 
) const
private

◆ probeInterpolationAcross()

bool dudley::DudleyDomain::probeInterpolationAcross ( int  functionSpaceType_source,
const escript::AbstractDomain targetDomain,
int  functionSpaceType_target 
) const
virtual

determines whether interpolation from source to target is possible.

Implements escript::AbstractDomain.

Referenced by getDim().

◆ probeInterpolationOnDomain()

bool dudley::DudleyDomain::probeInterpolationOnDomain ( int  functionSpaceType_source,
int  functionSpaceType_target 
) const
virtual

◆ randomFill()

escript::Data dudley::DudleyDomain::randomFill ( const escript::DataTypes::ShapeType shape,
const escript::FunctionSpace what,
long  seed,
const boost::python::tuple &  filter 
) const
virtual

Fills the data object with filtered random values.

Implements escript::AbstractDomain.

References escript::Data::getExpandedVectorReference(), and escript::randomFillArray().

Referenced by supportsContactElements().

◆ read()

escript::Domain_ptr dudley::DudleyDomain::read ( escript::JMPI  mpiInfo,
const std::string &  filename,
bool  optimize 
)
static

reads a mesh from a fly file. For MPI parallel runs fans out the mesh to multiple processes.

Parameters
mpiInfothe MPI information structure
fileNamethe name of the file
optimizewhether to optimize the node labels

References dudley::NodeFile::allocTable(), dudley::NodeFile::Coordinates, DudleyDomain(), getNodes(), escript::AbstractDomain::getPtr(), dudley::NodeFile::globalDegreesOfFreedom, dudley::NodeFile::Id, INDEX2, MPI_DOUBLE, MPI_INT, prepare(), resolveNodeIds(), setElements(), setFaceElements(), setPoints(), setTagMap(), and dudley::NodeFile::Tag.

◆ readGmsh()

escript::Domain_ptr dudley::DudleyDomain::readGmsh ( escript::JMPI  mpiInfo,
const std::string &  filename,
int  numDim,
bool  optimize 
)
static

◆ relabelElementNodes()

void dudley::DudleyDomain::relabelElementNodes ( const index_t *  newNode,
index_t  offset 
)

assigns new node reference numbers to all element files. If k is the old node, the new node is newNode[k-offset].

References m_elements, m_faceElements, m_points, and dudley::ElementFile::relabelNodes().

Referenced by resolveNodeIds(), and supportsContactElements().

◆ resolveNodeIds()

void dudley::DudleyDomain::resolveNodeIds ( )
private

Initially the element nodes refer to the numbering defined by the global id assigned to the nodes in the NodeFile. It is also not ensured that all nodes referred by an element are actually available on the process. At the output, a local node labeling is used and all nodes are available. In particular the numbering of the element nodes is between 0 and Nodes->numNodes. The function does not create a distribution of the degrees of freedom.

References dudley::NodeFile::allocTable(), ESYS_ASSERT, dudley::NodeFile::gather_global(), getDim(), dudley::ElementFile::getNodeRange(), escript::DataTypes::index_t_max(), m_elements, m_faceElements, m_mpiInfo, m_nodes, m_points, markNodes(), MPI_MAX, dudley::util::packMask(), and relabelElementNodes().

Referenced by create2D(), create3D(), distributeByRankOfDOF(), read(), readGmsh(), and supportsContactElements().

◆ setElements()

void dudley::DudleyDomain::setElements ( ElementFile elements)

replaces the element file by elements

References m_elements.

Referenced by create2D(), create3D(), getNodes(), load(), read(), and readGmsh().

◆ setFaceElements()

void dudley::DudleyDomain::setFaceElements ( ElementFile elements)

replaces the face element file by elements

References m_faceElements.

Referenced by create2D(), create3D(), getElements(), load(), read(), and readGmsh().

◆ setFunctionSpaceTypeNames()

void dudley::DudleyDomain::setFunctionSpaceTypeNames ( )

◆ setNewX()

void dudley::DudleyDomain::setNewX ( const escript::Data arg)
virtual

◆ setPoints()

void dudley::DudleyDomain::setPoints ( ElementFile elements)

replaces the point element file by elements

References m_points.

Referenced by create2D(), create3D(), getFaceElements(), load(), read(), and readGmsh().

◆ setTagMap()

void dudley::DudleyDomain::setTagMap ( const std::string &  name,
int  tag 
)
virtual

sets a map from a clear tag name to a tag key

Parameters
nameInput - tag name.
tagInput - tag key.

Implements escript::AbstractDomain.

References m_tagMap.

Referenced by BOOST_PYTHON_MODULE(), create2D(), create3D(), getDim(), load(), read(), and readGmsh().

◆ setTags()

void dudley::DudleyDomain::setTags ( int  functionSpaceType,
int  newTag,
const escript::Data mask 
) const
virtual

assigns new tag newTag to all samples of functionspace with a positive value of mask for any its sample point.

Implements escript::AbstractDomain.

References dudley::DegreesOfFreedom, dudley::Elements, dudley::FaceElements, m_elements, m_faceElements, m_nodes, m_points, dudley::Nodes, dudley::Points, dudley::ReducedElements, dudley::ReducedFaceElements, dudley::ElementFile::setTags(), and dudley::NodeFile::setTags().

Referenced by getDim().

◆ setToGradient()

void dudley::DudleyDomain::setToGradient ( escript::Data grad,
const escript::Data arg 
) const
virtual

◆ setToIntegrals()

void dudley::DudleyDomain::setToIntegrals ( std::vector< double > &  integrals,
const escript::Data arg 
) const
virtual

◆ setToNormal()

void dudley::DudleyDomain::setToNormal ( escript::Data out) const
virtual

copies the surface normals at data points into out. The actual function space to be considered is defined by out. out has to be defined on this.

Implements escript::AbstractDomain.

References dudley::Assemble_getNormal(), dudley::FaceElements, escript::FunctionSpace::getDomain(), escript::Data::getFunctionSpace(), escript::FunctionSpace::getTypeCode(), m_faceElements, m_nodes, and dudley::ReducedFaceElements.

Referenced by getDim().

◆ setToSize()

void dudley::DudleyDomain::setToSize ( escript::Data out) const
virtual

copies the size of samples into out. The actual function space to be considered is defined by out. out has to be defined on this.

Implements escript::AbstractDomain.

References dudley::Assemble_getSize(), dudley::DegreesOfFreedom, dudley::Elements, dudley::FaceElements, escript::Data::getFunctionSpace(), escript::FunctionSpace::getTypeCode(), m_elements, m_faceElements, m_nodes, dudley::Nodes, dudley::Points, dudley::ReducedElements, and dudley::ReducedFaceElements.

Referenced by getDim().

◆ setToX()

void dudley::DudleyDomain::setToX ( escript::Data arg) const
virtual

copies the location of data points into arg. The domain of arg has to match this. has to be implemented by the actual Domain adapter.

Implements escript::AbstractDomain.

References dudley::Assemble_NodeCoordinates(), escript::continuousFunction(), escript::FunctionSpace::getDomain(), escript::Data::getFunctionSpace(), escript::FunctionSpace::getTypeCode(), interpolateOnDomain(), m_nodes, dudley::Nodes, and escript::Vector().

Referenced by getDim().

◆ showTagNames()

string dudley::DudleyDomain::showTagNames ( ) const
virtual

Returns all tag names in a single string sperated by commas.

Implements escript::AbstractDomain.

References m_tagMap.

Referenced by BOOST_PYTHON_MODULE(), and getDim().

◆ supportsContactElements()

virtual bool dudley::DudleyDomain::supportsContactElements ( ) const
inlinevirtual

◆ updateTagList()

void dudley::DudleyDomain::updateTagList ( )
private

regenerates list of tags in use for node file and element files

References m_elements, m_faceElements, m_nodes, m_points, dudley::NodeFile::updateTagList(), and dudley::ElementFile::updateTagList().

Referenced by prepare(), and supportsContactElements().

◆ write()

void dudley::DudleyDomain::write ( const std::string &  fileName) const
virtual

writes the current mesh to a file with the given name in the fly file format.

Parameters
fileNameInput - The name of the file to write to.

Implements escript::AbstractDomain.

References dudley::NodeFile::Coordinates, getDim(), dudley::NodeFile::getNumNodes(), dudley::NodeFile::globalDegreesOfFreedom, dudley::NodeFile::Id, INDEX2, m_elements, m_faceElements, m_mpiInfo, m_name, m_nodes, m_points, m_tagMap, dudley::NodeFile::Tag, and writeElementInfo().

Referenced by BOOST_PYTHON_MODULE(), and getMPIComm().

◆ writeElementInfo()

void dudley::DudleyDomain::writeElementInfo ( std::ostream &  stream,
const ElementFile e,
const std::string &  defaultType 
) const
private

Member Data Documentation

◆ m_elements

ElementFile* dudley::DudleyDomain::m_elements
private

◆ m_faceElements

ElementFile* dudley::DudleyDomain::m_faceElements
private

◆ m_functionSpaceTypeNames

DudleyDomain::FunctionSpaceNamesMapType dudley::DudleyDomain::m_functionSpaceTypeNames
staticprivate

◆ m_mpiInfo

escript::JMPI dudley::DudleyDomain::m_mpiInfo
private

◆ m_name

std::string dudley::DudleyDomain::m_name
private

domain description

Referenced by dump(), Print_Mesh_Info(), and write().

◆ m_nodes

NodeFile* dudley::DudleyDomain::m_nodes
private

◆ m_points

ElementFile* dudley::DudleyDomain::m_points
private

◆ m_tagMap

TagMap dudley::DudleyDomain::m_tagMap
private

the tag map mapping names to tag keys

Referenced by dump(), getTag(), isValidTagName(), Print_Mesh_Info(), setTagMap(), showTagNames(), and write().


The documentation for this class was generated from the following files: