Geophysical Inversion and Modeling Library
jenkins-pyGIMLi-201-SUCCESS-1-g74330efc
|
Public Types | |
typedef RegionMarkerPLC | RegionMarker |
typedef std::vector< RegionMarker > | RegionMarkerList |
typedef RVector3 | HoleMarker |
typedef std::vector< RVector3 > | HoleMarkerList |
Public Member Functions | |
Mesh (Index dim=2, bool isGeometry=false) | |
Mesh (const std::string &filename, bool createNeighbourInfos=true) | |
Mesh (const Mesh &mesh) | |
Mesh & | operator= (const Mesh &mesh) |
~Mesh () | |
void | clear () |
void | setStaticGeometry (bool stat) |
bool | staticGeometry () const |
void | setGeometry (bool b) |
bool | isGeometry () const |
void | setDimension (uint dim) |
uint | dimension () const |
uint | dim () const |
Node * | createNode (double x, double y, double z, int marker=0) |
Node * | createNode (const Node &node) |
Node * | createNode (const RVector3 &pos, int marker=0) |
Node * | createSecondaryNode (const RVector3 &pos, double tol=-1) |
Node * | createNodeWithCheck (const RVector3 &pos, double tol=1e-6, bool warn=false, bool edgeCheck=false) |
Boundary * | createBoundary (std::vector< Node * > &nodes, int marker=0, bool check=true) |
Boundary * | createBoundary (const IndexArray &nodes, int marker=0, bool check=true) |
Boundary * | createBoundary (const Boundary &bound, bool check=true) |
Boundary * | createBoundary (const Cell &cell, bool check=true) |
Boundary * | createNodeBoundary (Node &n1, int marker=0, bool check=true) |
Boundary * | createEdge (Node &n1, Node &n2, int marker=0, bool check=true) |
Boundary * | createEdge3 (Node &n1, Node &n2, Node &n3, int marker=0, bool check=true) |
Boundary * | createTriangleFace (Node &n1, Node &n2, Node &n3, int marker=0, bool check=true) |
Boundary * | createQuadrangleFace (Node &n1, Node &n2, Node &n3, Node &n4, int marker=0, bool check=true) |
Cell * | createCell (int marker=0) |
Cell * | createCell (std::vector< Node * > &nodes, int marker=0) |
Cell * | createCell (const IndexArray &nodes, int marker=0) |
Cell * | createCell (const Cell &cell) |
Cell * | createTriangle (Node &n1, Node &n2, Node &n3, int marker=0) |
Cell * | createQuadrangle (Node &n1, Node &n2, Node &n3, Node &n4, int marker=0) |
Cell * | createTetrahedron (Node &n1, Node &n2, Node &n3, Node &n4, int marker=0) |
Cell * | copyCell (const Cell &cell, double tol=1e-6) |
Boundary * | copyBoundary (const Boundary &bound, double tol=1e-6, bool check=true) |
void | deleteCells (const std::vector< Cell * > &cells) |
void | create1DGrid (const RVector &x) |
void | create2DGrid (const RVector &x, const RVector &y, int markerType=0, bool worldBoundaryMarker=false) |
void | create3DGrid (const RVector &x, const RVector &y, const RVector &z, int markerType=0, bool worldBoundaryMarker=false) |
void | createGrid (const RVector &x) |
void | createGrid (const RVector &x, const RVector &y, int markerType=0, bool worldBoundaryMarker=false) |
void | createGrid (const RVector &x, const RVector &y, const RVector &z, int markerType=0, bool worldBoundaryMarker=false) |
void | createHull_ (const Mesh &mesh) |
Mesh | createHull () const |
void | createClosedGeometry (const std::vector< RVector3 > &vPos, int nSegments, double dxInner) |
void | createClosedGeometryParaMesh (const std::vector< RVector3 > &vPos, int nSegments, double dxInner) |
void | createClosedGeometryParaMesh (const std::vector< RVector3 > &vPos, int nSegments, double dxInner, const std::vector< RVector3 > &addit) |
Mesh | createH2 () const |
Mesh | createP2 () const |
void | createMeshByCells (const Mesh &mesh, const std::vector< Cell * > &cells) |
void | createMeshByBoundaries (const Mesh &mesh, const std::vector< Boundary * > &bounds) |
Mesh | createMeshByCellIdx (const IndexArray &idxList) |
void | createMeshByCellIdx (const Mesh &mesh, const IndexArray &idxList) |
void | createMeshByMarker (const Mesh &mesh, int from, int to=-1) |
Mesh | createSubMesh (const std::vector< Cell * > &cells) const |
Mesh | createSubMesh (const std::vector< Boundary * > &bounds) const |
Mesh | createSubMesh (const std::vector< Node * > &nodes) const |
void | showInfos () |
Show some infos. | |
const std::vector< Node *> & | nodes () const |
const std::vector< Cell *> & | cells () const |
const std::vector< Boundary *> & | boundaries () const |
std::vector< Node *> | nodes (const IndexArray &ids) const |
std::vector< Node *> | nodes (const BVector &b) const |
std::vector< Cell *> | cells (const IndexArray &ids) const |
std::vector< Cell *> | cells (const BVector &b) const |
std::vector< Boundary *> | boundaries (const IndexArray &ids) const |
std::vector< Boundary *> | boundaries (const BVector &b) const |
Index | nodeCount (bool withSecNodes=false) const |
Node & | node (Index i) const |
Node & | node (Index i) |
Index | secondaryNodeCount () const |
Node & | secondaryNode (Index id) const |
Node & | secondaryNode (Index id) |
Index | cellCount () const |
Cell & | cell (Index i) const |
Cell & | cell (Index i) |
Index | boundaryCount () const |
Boundary & | boundary (Index i) const |
Boundary & | boundary (Index i) |
R3Vector | positions (bool withSecNodes=false) const |
R3Vector | positions (const IndexArray &idx) const |
R3Vector | nodeCenters () const |
R3Vector | cellCenters () const |
R3Vector | cellCenter () const |
R3Vector | boundaryCenters () const |
RVector & | cellSizes () const |
RVector & | boundarySizes () const |
R3Vector & | boundarySizedNormals () const |
void | setBoundaryMarker (const IndexArray &ids, int marker) |
void | setBoundaryMarkers (const IndexArray &ids, int marker) |
void | setBoundaryMarkers (const IVector &marker) |
IVector | boundaryMarkers () const |
IVector | boundaryMarker () const |
IVector | nodeMarkers () const |
IVector | nodeMarker () const |
IndexArray | findNodesIdxByMarker (int marker) const |
std::vector< Boundary *> | findBoundaryByMarker (int marker) const |
std::vector< Boundary *> | findBoundaryByMarker (int from, int to) const |
Cell * | findCell (const RVector3 &pos, size_t &counter, bool extensive) const |
Cell * | findCell (const RVector3 &pos, bool extensive=true) const |
Index | findNearestNode (const RVector3 &pos) |
std::vector< Cell *> | findCellByMarker (int from, int to=0) const |
std::vector< Cell *> | findCellByAttribute (double from, double to=0.0) const |
std::vector< Cell *> | findCellsAlongRay (const RVector3 &start, const RVector3 &dir, R3Vector &pos) const |
void | fillEmptyCells (const std::vector< Cell * > &emptyList, double background=-1.0) |
void | prolongateEmptyCellsValues (RVector &vals, double background=-1.0) const |
void | recountNodes () |
void | sortNodes (const IndexArray &perm) |
bool | neighboursKnown () const |
void | cleanNeighbourInfos () |
void | createNeighbourInfos (bool force=false) |
void | createNeighbourInfosCell_ (Cell *c) |
void | relax () |
void | smooth (bool nodeMoving=true, bool edgeSwapping=true, uint smoothFunction=1, uint smoothIteration=10) |
Mesh & | scale (const RVector3 &s) |
Mesh & | scale (const double &s) |
Mesh & | translate (const RVector3 &t) |
Mesh & | rotate (const RVector3 &r) |
void | swapCoordinates (Index i, Index j) |
template<class Matrix > | |
Mesh & | transform (const Matrix &mat) |
int | save (const std::string &fileName, IOFormat format=Binary) const |
int | saveAscii (const std::string &fileName) const |
int | saveBinary (const std::string &fileName) const |
void | load (const std::string &fileName, bool createNeighbours=true, IOFormat format=Binary) |
void | loadAscii (const std::string &fileName) |
void | importMod (const std::string &fileName) |
void | importVTK (const std::string &fbody) |
void | importVTU (const std::string &fbody) |
void | importSTL (const std::string &fileName, bool isBinary=false) |
void | loadBinary (const std::string &fileName) |
void | saveBinaryV2 (const std::string &fbody) const |
void | loadBinaryV2 (const std::string &fbody) |
int | exportSimple (const std::string &fbody, const RVector &data) const |
int | exportMidCellValue (const std::string &fileName, const RVector &data1, const RVector &data2) const |
int | exportMidCellValue (const std::string &fileName, const RVector &data1) const |
void | exportVTK (const std::string &fbody, const std::map< std::string, RVector > &data, const std::vector< RVector3 > &vec, bool writeCells=true) const |
void | exportVTK (const std::string &fbody, const std::map< std::string, RVector > &data, bool writeCells=true) const |
void | exportVTK (const std::string &fbody, bool writeCells=true) const |
void | exportVTK (const std::string &fbody, const std::vector< RVector3 > &vec, bool writeCells=true) const |
void | exportVTK (const std::string &fbody, const RVector &arr) const |
void | readVTKPoints_ (std::fstream &file, const std::vector< std::string > &row) |
void | readVTKCells_ (std::fstream &file, const std::vector< std::string > &row) |
void | readVTKScalars_ (std::fstream &file, const std::vector< std::string > &row) |
void | readVTKPolygons_ (std::fstream &file, const std::vector< std::string > &row) |
void | exportVTU (const std::string &filename, bool binary=false) const |
void | exportBoundaryVTU (const std::string &fbody, bool binary=false) const |
void | addVTUPiece_ (std::fstream &file, const Mesh &mesh, const std::map< std::string, RVector > &data) const |
void | exportAsTetgenPolyFile (const std::string &filename) |
void | addExportData (const std::string &name, const RVector &data) |
RVector | exportData (const std::string &name) const |
const std::map< std::string, RVector > & | exportDataMap () const |
void | setExportDataMap (const std::map< std::string, RVector > &eMap) |
void | clearExportData () |
void | addData (const std::string &name, const CVector &data) |
void | addData (const std::string &name, const RVector &data) |
RVector | data (const std::string &name) const |
bool | haveData (const std::string &name) const |
const std::map< std::string, RVector > & | dataMap () const |
void | setDataMap (const std::map< std::string, RVector > m) |
void | dataInfo () const |
void | setCommentString (const std::string &commentString) |
const std::string & | commentString () const |
void | mapCellAttributes (const std::map< float, float > &aMap) |
void | mapAttributeToParameter (const IndexArray &cellMapIndex, const RVector &attributeMap, double defaultVal) |
void | mapBoundaryMarker (const std::map< int, int > &aMap) |
void | setCellAttributes (const RVector &attribute) |
void | setCellAttributes (double attribute) |
RVector | cellAttributes () const |
void | setCellMarkers (const IndexArray &ids, int marker) |
void | setCellMarkers (const IVector &marker) |
void | setCellMarkers (const RVector &attribute) |
IVector | cellMarkers () const |
double | xmin () const |
double | ymin () const |
double | zmin () const |
double | xmax () const |
double | ymax () const |
double | zmax () const |
const BoundingBox | boundingBox () const |
RSparseMapMatrix | interpolationMatrix (const R3Vector &q) |
void | interpolationMatrix (const R3Vector &q, RSparseMapMatrix &I) |
RSparseMapMatrix & | cellToBoundaryInterpolation () const |
RVector | divergence (const R3Vector &V) const |
R3Vector | boundaryDataToCellGradient (const RVector &boundaryData) const |
R3Vector | cellDataToBoundaryGradient (const RVector &cellData) const |
R3Vector | cellDataToBoundaryGradient (const RVector &cellData, const R3Vector &cellGradient) const |
void | addRegionMarker (const RVector3 &pos, int marker, double area=0) |
void | addRegionMarker (const RegionMarker ®) |
const RegionMarkerList & | regionMarker () const |
void | addHoleMarker (const RVector3 &pos) |
const HoleMarkerList & | holeMarker () const |
Protected Member Functions | |
void | copy_ (const Mesh &mesh) |
void | findRange_ () const |
Node * | createNodeGC_ (const RVector3 &pos, int marker) |
Node * | createNode_ (const RVector3 &pos, int marker) |
Node * | createSecondaryNode_ (const RVector3 &pos) |
template<class B > | |
Boundary * | createBoundary_ (std::vector< Node * > &nodes, int marker, int id) |
template<class B > | |
Boundary * | createBoundaryChecked_ (std::vector< Node * > &nodes, int marker, bool check=true) |
template<class C > | |
Cell * | createCell_ (std::vector< Node * > &nodes, int marker, int id) |
Node * | createRefinementNode_ (Node *n0, Node *n1, std::map< std::pair< Index, Index >, Node * > &nodeMatrix) |
void | createRefined_ (const Mesh &mesh, bool p2, bool r2) |
Cell * | findCellBySlopeSearch_ (const RVector3 &pos, Cell *start, size_t &count, bool tagging) const |
void | fillKDTree_ () const |
Protected Attributes | |
std::vector< Node *> | nodeVector_ |
std::vector< Node *> | secNodeVector_ |
std::vector< Boundary *> | boundaryVector_ |
std::vector< Cell *> | cellVector_ |
uint | dimension_ |
RVector3 | minRange_ |
RVector3 | maxRange_ |
bool | rangesKnown_ |
bool | neighboursKnown_ |
KDTreeWrapper * | tree_ |
bool | staticGeometry_ |
bool | isGeometry_ |
RVector | cellSizesCache_ |
RVector | boundarySizesCache_ |
R3Vector | boundarySizedNormCache_ |
RSparseMapMatrix * | cellToBoundaryInterpolationCache_ |
bool | oldTet10NumberingStyle_ |
std::map< std::string, RVector > | exportDataMap_ |
std::string | commentString_ |
RegionMarkerList | regionMarker_ |
HoleMarkerList | holeMarker_ |
GIMLI::Mesh::Mesh | ( | Index | dim = 2 , |
bool | isGeometry = false |
||
) |
Default constructor, create empty mesh with dimension dim If this mesh is supposed to be a geometry definition, all created nodes will be checked for duplicates.
GIMLI::Mesh::Mesh | ( | const std::string & | filename, |
bool | createNeighbourInfos = true |
||
) |
Constructor, read mesh from filename
References load().
GIMLI::Mesh::Mesh | ( | const Mesh & | mesh | ) |
Copy constructor.
GIMLI::Mesh::~Mesh | ( | ) |
Default destructor.
References clear().
|
inline |
Add data to the mesh that will be saved with by using the binary mesh format v.2. or exported with the appropriate name in VTK format, if the size of data equals the amount of nodes, cells or boundaries.
void GIMLI::Mesh::addExportData | ( | const std::string & | name, |
const RVector & | data | ||
) |
DEPRECATED will be removed. Add data to the mesh that will be saved with by using the binary mesh format v.2. or exported with the appropriate name in VTK format, if the size of data equals the amount of nodes, cells or boundaries.
References data().
Referenced by GIMLI::interpolate(), GIMLI::Inversion< ModelValType >::oneStep(), and GIMLI::DCMultiElectrodeModelling::response().
void GIMLI::Mesh::addHoleMarker | ( | const RVector3 & | pos | ) |
Add a hole marker for tetgen or triangle creation if the mesh is a PLC
Referenced by addRegionMarker(), and operator=().
void GIMLI::Mesh::addRegionMarker | ( | const RVector3 & | pos, |
int | marker, | ||
double | area = 0 |
||
) |
Add a region marker for tetgen or triangle creation if the mesh is a PLC, if area is < 0 a hole is added.
References addHoleMarker().
Referenced by operator=(), and smooth().
void GIMLI::Mesh::addVTUPiece_ | ( | std::fstream & | file, |
const Mesh & | mesh, | ||
const std::map< std::string, RVector > & | data | ||
) | const |
Internal function for exporting VTU
References GIMLI::Matrix< ValueType >::cols(), create2DGrid(), GIMLI::BaseEntity::id(), GIMLI::loadMatrixCol(), GIMLI::MeshEntity::rtti(), GIMLI::unique(), GIMLI::x(), and GIMLI::y().
Referenced by exportBoundaryVTU(), and exportVTU().
|
inline |
Return const reference to all boundaries
Referenced by GIMLI::addTriangleBoundary(), boundaries(), GIMLI::calcGBounds(), and exportBoundaryVTU().
std::vector< Boundary *> GIMLI::Mesh::boundaries | ( | const IndexArray & | ids | ) | const |
Return vector of boundaries from index list
Return a vector of boundary ptrs matching BVector b.
References boundaries(), and GIMLI::find().
R3Vector GIMLI::Mesh::boundaryCenters | ( | ) | const |
Return a vector of all center positions for all boundaries
References GIMLI::MeshEntity::center().
Interpolate boundary based values to cell based gradients.
References boundarySizedNormals(), GIMLI::BaseEntity::id(), and GIMLI::Boundary::leftCell().
Referenced by cellDataToBoundaryGradient().
IVector GIMLI::Mesh::boundaryMarker | ( | ) | const |
DEPRECATED
References boundaryMarkers().
IVector GIMLI::Mesh::boundaryMarkers | ( | ) | const |
Return a vector of all boundary marker
Referenced by boundaryMarker().
R3Vector & GIMLI::Mesh::boundarySizedNormals | ( | ) | const |
Return the reference to the vector of scaled normal directions for each boundary. Cached for static geometry and will be build on first call. Not thread safe, perhaps not python GC safe. Return . Where
is the size and
the normal direction for the i-th boundary. If you want to use this, i.e. for the calculation of inside or outside flow through the boundary, you need to recognize the orientation of this boundary to the cell the flow goes into or comes from. For the left cell neighbor the normal direction should be always the outer normal.
References createNeighbourInfos(), dim(), GIMLI::Vector< ValueType >::resize(), and staticGeometry_.
Referenced by boundaryDataToCellGradient(), and divergence().
RVector & GIMLI::Mesh::boundarySizes | ( | ) | const |
Return the reference to a RVector of all boundary sizes. Cached for static geometry.
References GIMLI::Vector< ValueType >::resize(), GIMLI::MeshEntity::size(), and staticGeometry_.
RVector GIMLI::Mesh::cellAttributes | ( | ) | const |
Return a RVector of all cell attributes
Referenced by GIMLI::DCMultiElectrodeModelling::createJacobian(), GIMLI::dcfemDomainAssembleStiffnessMatrix(), exportVTK(), exportVTU(), fillEmptyCells(), GIMLI::DCMultiElectrodeModelling::mapERTModel(), GIMLI::Inversion< ModelValType >::oneStep(), operator=(), and GIMLI::DCMultiElectrodeModelling::response().
R3Vector GIMLI::Mesh::cellCenters | ( | ) | const |
Return a vector of all cell center positions
References GIMLI::MeshEntity::center().
Interpolate cell based values to boundary based gradients.
References boundaryDataToCellGradient(), and cellToBoundaryInterpolation().
R3Vector GIMLI::Mesh::cellDataToBoundaryGradient | ( | const RVector & | cellData, |
const R3Vector & | cellGradient | ||
) | const |
Interpolate cell based values to boundary based gradients with a given cell Gradient.
References GIMLI::MeshEntity::center(), GIMLI::BaseEntity::id(), GIMLI::Boundary::leftCell(), and GIMLI::Boundary::norm().
IVector GIMLI::Mesh::cellMarkers | ( | ) | const |
Return a vector of all cell marker
Referenced by GIMLI::CreateSensitivityColMT< ValueType >::calc1(), createH2(), createRefined_(), and GIMLI::prepExportSensitivityData().
std::vector< Cell *> GIMLI::Mesh::cells | ( | const IndexArray & | ids | ) | const |
Return vector of cells from index list
Return a vector of cells ptrs matching BVector b.
References GIMLI::find().
RVector & GIMLI::Mesh::cellSizes | ( | ) | const |
Return the reference to a RVector of all cell sizes. Cached for static geometry.
References GIMLI::Vector< ValueType >::resize(), GIMLI::MeshEntity::size(), and staticGeometry_.
Referenced by divergence(), and GIMLI::prepExportSensitivityData().
RSparseMapMatrix & GIMLI::Mesh::cellToBoundaryInterpolation | ( | ) | const |
Return the reference to the matrix for cell value to boundary value interpolation matrix.
References GIMLI::MeshEntity::center(), GIMLI::BaseEntity::id(), GIMLI::Boundary::leftCell(), and staticGeometry_.
Referenced by cellDataToBoundaryGradient().
void GIMLI::Mesh::cleanNeighbourInfos | ( | ) |
Remove from each boundary the ptr to the corresponding left and right cell
Referenced by createNeighbourInfos().
void GIMLI::Mesh::clear | ( | ) |
Clear all data, inclusive all caches.
Referenced by create2DGrid(), create3DGrid(), createHull(), createMeshByBoundaries(), createMeshByCellIdx(), createMeshByCells(), createNeighbourInfosCell_(), createRefined_(), exportVTK(), loadBinary(), loadBinaryV2(), operator=(), GIMLI::TriangleWrapper::transformTriangleToMesh_(), and ~Mesh().
void GIMLI::Mesh::clearExportData | ( | ) |
Empty the data map.
Referenced by GIMLI::Inversion< ModelValType >::oneStep().
|
inline |
Return comment for VTK Ascii export headline.
Boundary * GIMLI::Mesh::copyBoundary | ( | const Boundary & | bound, |
double | tol = 1e-6 , |
||
bool | check = true |
||
) |
Create a boundary as a copy of a boundary from an alternative mesh. Each node of the new cell will be created with duplication check and reused if there is already a node withing the tolerance distance tol. tol=-1 disables this duplication check.
References createNodeWithCheck().
Create a cell as a copy of a cell from an alternative mesh. Each node of the new cell will be created with duplication check and reused if there is already a node withing the tolerance distance tol. tol=-1 disables this duplication check.
References createCell(), and createNodeWithCheck().
Referenced by GIMLI::addTriangleBoundary().
void GIMLI::Mesh::create2DGrid | ( | const RVector & | x, |
const RVector & | y, | ||
int | markerType = 0 , |
||
bool | worldBoundaryMarker = false |
||
) |
Default boundary marker are -x[1], +x[2], +z[3], -z[4]. If worldBoundaryMarker is set it becomes +z[-1], else[-2].
References clear(), createCell(), createNeighbourInfos(), GIMLI::Boundary::leftCell(), setDimension(), and GIMLI::unique().
Referenced by addVTUPiece_(), and GIMLI::createMesh2D().
void GIMLI::Mesh::create3DGrid | ( | const RVector & | x, |
const RVector & | y, | ||
const RVector & | z, | ||
int | markerType = 0 , |
||
bool | worldBoundaryMarker = false |
||
) |
Default boundary marker are -x[1], +x[2], +z[3], -z[4], -y[5], +z[6]. If worldBoundaryMarker is set it becomes +z[-1], else[-2]. You can exportBoundaryVTU and take a look with Paraview.
References clear(), createCell(), createNeighbourInfos(), GIMLI::Boundary::leftCell(), setDimension(), and GIMLI::unique().
Referenced by GIMLI::createMesh3D(), and exportVTK().
Boundary * GIMLI::Mesh::createBoundary | ( | const IndexArray & | nodes, |
int | marker = 0 , |
||
bool | check = true |
||
) |
Create a boundary from the given node indieces
References GIMLI::BaseEntity::id().
Cell * GIMLI::Mesh::createCell | ( | int | marker = 0 | ) |
Create empty cell without a node or a shape.
Referenced by copyCell(), create2DGrid(), create3DGrid(), createCell(), GIMLI::createMesh2D(), GIMLI::createMesh3D(), createMeshByCells(), createNeighbourInfosCell_(), createRefined_(), exportVTK(), GIMLI::interpolateSurface(), loadBinary(), loadBinaryV2(), and operator=().
Cell * GIMLI::Mesh::createCell | ( | const IndexArray & | nodes, |
int | marker = 0 |
||
) |
Create a cell from the given node indieces
References createCell(), and GIMLI::BaseEntity::id().
|
inline |
Create one dimensional grid. Boundary on the domain border will get marker: 1,2 for: left, right.
Referenced by GIMLI::createGrid().
|
inline |
Create two dimensional grid. Boundary on the domain border will get the marker: 1,2,3,4 for: left, right, top, bottom
|
inline |
Create three dimensional grid. Boundary on the domain border will get the marker: 1,2,3,4,5,6 for: left, right, top, bottom, front, back
Mesh GIMLI::Mesh::createH2 | ( | ) | const |
Create and copy global H2 mesh of this mesh.
References cellMarkers(), createRefined_(), dimension(), and setCellAttributes().
Referenced by GIMLI::ModellingBase::startModel().
Mesh GIMLI::Mesh::createHull | ( | ) | const |
void GIMLI::Mesh::createMeshByBoundaries | ( | const Mesh & | mesh, |
const std::vector< Boundary * > & | bounds | ||
) |
Create a partly mesh without cells from mesh, based on a vector of ptrs to boundaries
Create new boundaries
References clear(), dim(), GIMLI::BaseEntity::id(), and setDimension().
Referenced by createSubMesh(), and exportBoundaryVTU().
Mesh GIMLI::Mesh::createMeshByCellIdx | ( | const IndexArray & | idxList | ) |
Create a new mesh that is a part from this mesh, based on cell-ids
References createMeshByCellIdx(), and dimension().
Referenced by createMeshByCellIdx(), createMeshByMarker(), and findNearestNode().
void GIMLI::Mesh::createMeshByCellIdx | ( | const Mesh & | mesh, |
const IndexArray & | idxList | ||
) |
Create a partly mesh from mesh, based on cell-ids
References clear(), createMeshByCells(), dim(), setDimension(), GIMLI::str(), and GIMLI::unique().
Create a partly mesh without cells from mesh, based on a vector of ptrs to boundaries
Create new cells
copy all boundary with marker != 0
check that all nodes have a common cell
Copy data if available
Create all remaining boundaries
References clear(), createCell(), createNeighbourInfos(), dataMap(), dim(), GIMLI::BaseEntity::id(), and setDimension().
Referenced by createMeshByCellIdx(), and createSubMesh().
void GIMLI::Mesh::createMeshByMarker | ( | const Mesh & | mesh, |
int | from, | ||
int | to = -1 |
||
) |
Create a partly mesh from mesh, based on meshs attributes. For a single attribute set to to 0, for unlimited set to to -1
References createMeshByCellIdx().
void GIMLI::Mesh::createNeighbourInfos | ( | bool | force = false | ) |
Search and set to each boundary the corresponding left and right cell.
References cleanNeighbourInfos(), and createNeighbourInfosCell_().
Referenced by GIMLI::addTriangleBoundary(), boundarySizedNormals(), create2DGrid(), create3DGrid(), createMeshByCells(), createNeighbourInfosCell_(), fillEmptyCells(), GIMLI::interpolateSurface(), load(), operator=(), GIMLI::setAllNeumannBoundaryConditions(), GIMLI::setDefaultBERTBoundaryConditions(), and smooth().
void GIMLI::Mesh::createNeighbourInfosCell_ | ( | Cell * | c | ) |
Create and store boundaries and neighboring information for this cell.
References GIMLI::Cell::boundaryNodes(), clear(), createCell(), createNeighbourInfos(), GIMLI::Cell::findNeighbourCell(), GIMLI::Boundary::leftCell(), setDimension(), and GIMLI::unique().
Referenced by createNeighbourInfos(), and findNearestNode().
Ensure is geometry check
References createNodeWithCheck(), GIMLI::x(), GIMLI::y(), and GIMLI::z().
Node * GIMLI::Mesh::createNodeWithCheck | ( | const RVector3 & | pos, |
double | tol = 1e-6 , |
||
bool | warn = false , |
||
bool | edgeCheck = false |
||
) |
Create new Node with duplication checks. Returns the already existing node when its within a tolerance distance to pos. If edgeCheck is set, any 2d (p1) boundary edges will be checked for any intersection with pos and splitted if necessary.
References dim(), GIMLI::KDTreeWrapper::insert(), GIMLI::KDTreeWrapper::nearest(), GIMLI::Boundary::rtti(), and GIMLI::str().
Referenced by copyBoundary(), copyCell(), createNodeGC_(), and importSTL().
Mesh GIMLI::Mesh::createP2 | ( | ) | const |
Create and copy global P2 mesh of this mesh.
References createRefined_(), dimension(), and GIMLI::BaseEntity::id().
Referenced by GIMLI::DCSRMultiElectrodeModelling::checkPrimpotentials_(), and GIMLI::ModellingBase::startModel().
|
protected |
Copy data if available
copy original marker into the new mesh
References cellMarkers(), clear(), createCell(), dataMap(), GIMLI::BaseEntity::id(), GIMLI::Cell::rtti(), GIMLI::Boundary::rtti(), setCellMarkers(), and GIMLI::Shape::xyz().
Referenced by createH2(), and createP2().
Create a secondary node, which is stored in an aditional list for additional use. If tolerance tol set to a value > 0, then it will be checked if there is already a node at this position and return a ptr to the existing node instead of creating a new.
References GIMLI::KDTreeWrapper::insert(), and GIMLI::KDTreeWrapper::nearest().
Referenced by operator=().
Syntactic sugar to extract a part of the mesh based on cells.
References createMeshByCells(), and dimension().
Syntactic sugar to extract a part of the mesh based on boundaries.
References createMeshByBoundaries(), and dimension().
Syntactic sugar to extract a part of the mesh based on nodes with associated cells and boundaries.
References dimension().
|
inline |
Return the data with a given name. If there is no such data an exception is thrown.
Referenced by addExportData(), exportVTK(), exportVTU(), GIMLI::getComplexResistivities(), and loadBinaryV2().
void GIMLI::Mesh::dataInfo | ( | ) | const |
Print data map info.
References GIMLI::str().
|
inline |
Return the full data map read only.
Referenced by createMeshByCells(), createRefined_(), and exportVTK().
void GIMLI::Mesh::deleteCells | ( | const std::vector< Cell * > & | cells | ) |
Delete all given cells from the given mesh. Warning will be really deleted.
References staticGeometry_.
|
inline |
Shortcut for dimension.
References GIMLI::x(), GIMLI::y(), and GIMLI::z().
Referenced by GIMLI::DCMultiElectrodeModelling::assembleStiffnessMatrixDCFEMByPass(), boundarySizedNormals(), GIMLI::TravelTimeDijkstraModelling::createDefaultStartModel(), createHull(), createMeshByBoundaries(), createMeshByCellIdx(), createMeshByCells(), createNodeWithCheck(), fillEmptyCells(), GIMLI::interpolate(), load(), loadBinary(), loadBinaryV2(), operator=(), prolongateEmptyCellsValues(), saveBinary(), and GIMLI::DCMultiElectrodeModelling::searchElectrodes_().
|
inline |
Return the dimension of this mesh.
Referenced by createH2(), GIMLI::DCMultiElectrodeModelling::createJacobian(), createMeshByCellIdx(), createP2(), createSubMesh(), saveBinary(), saveBinaryV2(), GIMLI::setAllNeumannBoundaryConditions(), and GIMLI::setDefaultBERTBoundaryConditions().
Return the divergence for each cell of a given vector field for each boundary. The divergence is calculated by simple 1 point boundary integration over each cell. Higher order integration needs to be implemented. Contact the author if you need this.
References boundarySizedNormals(), cellSizes(), GIMLI::BaseEntity::id(), and GIMLI::Boundary::leftCell().
void GIMLI::Mesh::exportBoundaryVTU | ( | const std::string & | fbody, |
bool | binary = false |
||
) | const |
Export the boundary of this mesh in vtu format: Visualization Toolkit Unstructured Points Data (http://www.vtk.org) Set Binary to true writes the datacontent in binary format. The file suffix .vtu will be added or substituted if .vtu or .vtk is found.
References addVTUPiece_(), boundaries(), and createMeshByBoundaries().
RVector GIMLI::Mesh::exportData | ( | const std::string & | name | ) | const |
DEPRECATED Return the data with a given name. If there is no such data an exception is thrown.
Referenced by GIMLI::interpolate().
|
inline |
DEPRECATED Return the full data map read only.
Referenced by GIMLI::interpolate(), operator=(), and GIMLI::prepExportSensitivityData().
int GIMLI::Mesh::exportMidCellValue | ( | const std::string & | fileName, |
const RVector & | data1, | ||
const RVector & | data2 | ||
) | const |
Very simple export filter. Write to file fileName: x_0 y_0 [z_0] data1_0 [data2]_0
x_n y_n [z_n] data1_n [data2]_n
with n = number of cells, and x y [z] == cell center
References GIMLI::MeshEntity::center(), GIMLI::BaseEntity::id(), GIMLI::Matrix< ValueType >::push_back(), and GIMLI::saveMatrixCol().
void GIMLI::Mesh::exportVTK | ( | const std::string & | fbody, |
bool | writeCells = true |
||
) | const |
Export mesh and whole exportData map
void GIMLI::Mesh::exportVTK | ( | const std::string & | fbody, |
const std::vector< RVector3 > & | vec, | ||
bool | writeCells = true |
||
) | const |
Export mesh and whole exportData map and vector data in vec
References dataMap().
void GIMLI::Mesh::exportVTK | ( | const std::string & | fbody, |
const RVector & | arr | ||
) | const |
Export mesh with one additional array that will called 'arr'
References cellAttributes(), clear(), create3DGrid(), createCell(), data(), dataMap(), GIMLI::BaseEntity::id(), GIMLI::nonZero(), positions(), GIMLI::Cell::rtti(), GIMLI::Boundary::rtti(), GIMLI::str(), GIMLI::x(), GIMLI::y(), and GIMLI::z().
void GIMLI::Mesh::exportVTU | ( | const std::string & | filename, |
bool | binary = false |
||
) | const |
Export the mesh in filename using vtu format: Visualization Toolkit Unstructured Points Data (http://www.vtk.org) Set binary to true writes the data content in binary format. The file suffix .vtu will be added or substituted if .vtu or .vtk is found. exportData, cell.marker and cell.attribute will be exported as data.
References addVTUPiece_(), cellAttributes(), and data().
void GIMLI::Mesh::fillEmptyCells | ( | const std::vector< Cell * > & | emptyList, |
double | background = -1.0 |
||
) |
DEPRECATED Prolongate the attribute of each cell in emptyList by the attribute from neighbouring cells. The attributes have to be lower than TOLERANCE. This function is called recursively until all zero-attribute-cells in emptyList are filled with an attribute greater than Zero.
References cellAttributes(), createNeighbourInfos(), dim(), GIMLI::Boundary::norm(), and smooth().
std::vector< Boundary *> GIMLI::Mesh::findBoundaryByMarker | ( | int | marker | ) | const |
Return an index list of all nodes that match the marker
Return a vector of boundary ptrs with the boundary marker equal marker.
Referenced by GIMLI::addTriangleBoundary().
std::vector< Boundary *> GIMLI::Mesh::findBoundaryByMarker | ( | int | from, |
int | to | ||
) | const |
Return a vector of boundary ptrs with the boundary marker between [from and to).
for to equal open end set to = MAX_INT
Return ptr to the cell that match position pos, counter holds amount of touch tests. Searching is done first by nearest-neighbour-kd-tree search, followed by slope-search if extensive is set. Return NULL if no cell can be found.
** *sigh, no luck with simple kd-tree search, try more expensive full slope search
References GIMLI::Vector< ValueType >::clear(), GIMLI::KDTreeWrapper::nearest(), and GIMLI::BaseEntity::untag().
Referenced by findCellsAlongRay(), GIMLI::interpolate(), interpolationMatrix(), and GIMLI::DCMultiElectrodeModelling::searchElectrodes_().
Shortcut for findCell(const RVector3 & pos, size_t & counter, bool extensive)
std::vector< Cell *> GIMLI::Mesh::findCellByAttribute | ( | double | from, |
double | to = 0.0 |
||
) | const |
Return vector of cell ptrs with attribute match the range [from .. to).
For single attribute match to is set to 0.0, for open end set to = -1.0
std::vector< Cell *> GIMLI::Mesh::findCellByMarker | ( | int | from, |
int | to = 0 |
||
) | const |
Return vector of cell ptrs with marker match the range [from .. to).
For single marker match to is set to 0, for open end set to = -1
Referenced by GIMLI::CreateSensitivityColMT< ValueType >::calc1(), and GIMLI::Region::resize().
std::vector< Cell *> GIMLI::Mesh::findCellsAlongRay | ( | const RVector3 & | start, |
const RVector3 & | dir, | ||
R3Vector & | pos | ||
) | const |
Return vector of cells that are intersected with a given ray from start to end. Intersecting positions, i.e., the travel path are stored in pos. Note this will not yet check if the ray lies completely along a boundary. This will probably fail and need to be implemented.
References GIMLI::Cell::boundary(), GIMLI::Vector< ValueType >::clean(), findCell(), GIMLI::Shape::intersectRay(), and GIMLI::KDTreeWrapper::nearest().
Index GIMLI::Mesh::findNearestNode | ( | const RVector3 & | pos | ) |
Return the index to the node of this mesh with the smallest distance to pos.
References createMeshByCellIdx(), createNeighbourInfosCell_(), GIMLI::BaseEntity::id(), GIMLI::Shape::isInside(), GIMLI::KDTreeWrapper::nearest(), GIMLI::BaseEntity::tag(), and GIMLI::BaseEntity::tagged().
IndexArray GIMLI::Mesh::findNodesIdxByMarker | ( | int | marker | ) | const |
Return an index vector of all nodes that match the marker
References GIMLI::find(), and nodeMarkers().
Referenced by GIMLI::DCMultiElectrodeModelling::assembleStiffnessMatrixDCFEMByPass(), and GIMLI::DCMultiElectrodeModelling::searchElectrodes_().
|
inline |
Return True if date with such a name exists.
Referenced by GIMLI::DCMultiElectrodeModelling::assembleStiffnessMatrixDCFEMByPass(), and GIMLI::getComplexResistivities().
|
inline |
Return read only reference for all defined hole regions.
Referenced by operator=(), and GIMLI::TriangleWrapper::transformMeshToTriangle_().
void GIMLI::Mesh::importSTL | ( | const std::string & | fileName, |
bool | isBinary = false |
||
) |
Import Ascii STL, and save triangles as Boundary(ies).
References createNodeWithCheck(), and GIMLI::toStr().
RSparseMapMatrix GIMLI::Mesh::interpolationMatrix | ( | const R3Vector & | q | ) |
Return the interpolation matrix I from all node positions to the query points q. I is a (len(q) x nodeCount()) SparseMapMatrix. To perform the interpolation just calculate the matrix vector product. uInterpolated = I.mult(uPerNode) or uInterpolated = I * uPerNode
void GIMLI::Mesh::interpolationMatrix | ( | const R3Vector & | q, |
RSparseMapMatrix & | I | ||
) |
|
inline |
Return if the mesh is a geometry definition.
Referenced by operator=().
void GIMLI::Mesh::load | ( | const std::string & | fileName, |
bool | createNeighbours = true , |
||
IOFormat | format = Binary |
||
) |
Load Mesh from file and try to import fileformat regarding file suffix. If createNeighbourInfos is set, the mesh is checked for consistency and missing boundaries will be created.
References createNeighbourInfos(), dim(), GIMLI::BaseEntity::id(), loadBinary(), loadBinaryV2(), and saveBinary().
Referenced by Mesh().
void GIMLI::Mesh::loadBinary | ( | const std::string & | fileName | ) |
Be carefull with interchanging binary meshs between 32-64bit architecture. Atm we save fixed int for counter and idx. We have to replace and test it with uint32 or uint16
References clear(), createCell(), dim(), nodeMarker(), setDimension(), setGeometry(), GIMLI::str(), and GIMLI::toStr().
Referenced by load().
void GIMLI::Mesh::loadBinaryV2 | ( | const std::string & | fbody | ) |
Load mesh in binary format v.2.0. should be possible to interchange on all little endian platforms. Format see saveBinaryV2. If something goes wrong while reading, an exception is thrown.
References clear(), createCell(), data(), dim(), setDimension(), GIMLI::str(), and GIMLI::toStr().
Referenced by load().
void GIMLI::Mesh::mapBoundaryMarker | ( | const std::map< int, int > & | aMap | ) |
Change all boundary marker that match bMap.first to bMap.second.
|
inline |
Return true if createNeighbourInfos is called once
Referenced by operator=().
R3Vector GIMLI::Mesh::nodeCenters | ( | ) | const |
DEPRECATED Return all node positions.
IVector GIMLI::Mesh::nodeMarker | ( | ) | const |
IVector GIMLI::Mesh::nodeMarkers | ( | ) | const |
Return a vector of all node marker
Referenced by findNodesIdxByMarker(), and nodeMarker().
std::vector< Node *> GIMLI::Mesh::nodes | ( | const IndexArray & | ids | ) | const |
Return vector of node from index list
Return a vector of nodes ptrs matching BVector b.
References GIMLI::find().
Copy assignment operator.
References addHoleMarker(), addRegionMarker(), cellAttributes(), clear(), createCell(), createNeighbourInfos(), createSecondaryNode(), dim(), exportDataMap(), holeMarker(), isGeometry(), neighboursKnown(), setCellAttributes(), setExportDataMap(), setGeometry(), setStaticGeometry(), and staticGeometry().
R3Vector GIMLI::Mesh::positions | ( | bool | withSecNodes = false | ) | const |
Return a vector of all node positions
Referenced by exportVTK(), and GIMLI::interpolate().
R3Vector GIMLI::Mesh::positions | ( | const IndexArray & | idx | ) | const |
Return a vector of node positions for an index vector
void GIMLI::Mesh::prolongateEmptyCellsValues | ( | RVector & | vals, |
double | background = -1.0 |
||
) | const |
Prolongate the empty (lower than TOLERANCE.) cell values in vals from its neighbouring cells. This function is called recursively until all zero-attribute-values in vals are filled with an attribute greater than Zero. RVector vals need to be of size cellCount(). If Background is unequal -1.0 all empty values will be set to background.
References dim(), GIMLI::find(), GIMLI::BaseEntity::id(), GIMLI::Boundary::norm(), and smooth().
Rotate mesh the with RVector3 r, r in radian, If you want to rotate in degree, use degToRad(const RVector3 & deg). And return a reference to the mesh (no copy)
int GIMLI::Mesh::saveBinary | ( | const std::string & | fileName | ) | const |
Be carefull with interchanging binary meshs between 32-64bit architecture. Atm we save fixed int for counter and idx. We have to replace and test it with uint32 or uint16
write vertex dummy-infos;
write cell dummy-infos
References dim(), dimension(), GIMLI::BaseEntity::id(), and GIMLI::Boundary::leftCell().
Referenced by load().
void GIMLI::Mesh::saveBinaryV2 | ( | const std::string & | fbody | ) | const |
Save mesh in binary format v.2.0. should be possible to interchange on all little endian platforms. Contains all export data, marker an attribute. If something goes wrong while writing, an exception is thrown. sizes uint8 (1 byte), uint16 (2 byte), uint32 (4 byte), uint64 (8 byte) Format:
References dimension(), GIMLI::BaseEntity::id(), and GIMLI::Boundary::leftCell().
Scale the mesh with RVector3 s. And return a reference to the mesh (no copy)
|
inline |
Scale the mesh with s. Shortcut for scale(RVector3(s,s,s))
|
inline |
DEPRECATED
void GIMLI::Mesh::setBoundaryMarkers | ( | const IndexArray & | ids, |
int | marker | ||
) |
Set the marker to all boundaries in index array.
void GIMLI::Mesh::setBoundaryMarkers | ( | const IVector & | marker | ) |
Set all cell marker the values in attribute.
void GIMLI::Mesh::setCellAttributes | ( | const RVector & | attribute | ) |
Set all cell attributes to the valaues in vector attribute.
References GIMLI::toStr().
Referenced by createH2(), GIMLI::DCMultiElectrodeModelling::createJacobian(), GIMLI::DCMultiElectrodeModelling::mapERTModel(), and operator=().
void GIMLI::Mesh::setCellAttributes | ( | double | attribute | ) |
Set all cell attributes to the scalar value: attribute
void GIMLI::Mesh::setCellMarkers | ( | const IndexArray & | ids, |
int | marker | ||
) |
Set the cell marker of all indices in ids to marker.
Referenced by GIMLI::createGrid(), and createRefined_().
void GIMLI::Mesh::setCellMarkers | ( | const IVector & | marker | ) |
Set all cell marker the values in attribute.
void GIMLI::Mesh::setCellMarkers | ( | const RVector & | attribute | ) |
Set all cell marker the values in attribute (casting to int)
References GIMLI::str().
|
inline |
Set the comment for VTK Ascii export headline.
|
inline |
Replace the datamap by m
|
inline |
Set the dimension of the mesh. [Default = 2]
Referenced by create2DGrid(), create3DGrid(), createMeshByBoundaries(), createMeshByCellIdx(), createMeshByCells(), createNeighbourInfosCell_(), loadBinary(), and loadBinaryV2().
|
inline |
DEPRECATED Set the full data map.
Referenced by operator=().
void GIMLI::Mesh::setGeometry | ( | bool | b | ) |
Mesh is marked as geometry definition or PLC so createNode will allways with check.
Referenced by loadBinary(), and operator=().
void GIMLI::Mesh::setStaticGeometry | ( | bool | stat | ) |
If the mesh is static in geometry and shape some useful informations are cached. (cell sizes, boundary sizes, ...) For dynamic meshes, i.e., node positions can be moved, you have to set staticGeometry to false to avoid any caching.
References staticGeometry_.
Referenced by operator=().
void GIMLI::Mesh::smooth | ( | bool | nodeMoving = true , |
bool | edgeSwapping = true , |
||
uint | smoothFunction = 1 , |
||
uint | smoothIteration = 10 |
||
) |
Smooth the mesh via moving all free nodes into the average of all neighboring nodes. Repeat this smoothIteration times. There is currently only this smoothFunction. EdgeSwapping is deactivated.
References addRegionMarker(), createNeighbourInfos(), GIMLI::KDTreeWrapper::insert(), GIMLI::KDTreeWrapper::size(), GIMLI::toStr(), and GIMLI::KDTreeWrapper::tree().
Referenced by fillEmptyCells(), and prolongateEmptyCellsValues().
|
inline |
Return true if this mesh have static geometry. [Default=True]
Referenced by operator=().
void GIMLI::Mesh::swapCoordinates | ( | Index | i, |
Index | j | ||
) |
Swap coordinate i with j for i and j lower then dimension of the mesh
apply a 4x4 transformation matrix to the whole mesh
References GIMLI::load().
Translate the mesh with RVector3 t. And return a reference to the mesh (no copy)
|
protected |
A static geometry mesh caches geometry informations.
Referenced by boundarySizedNormals(), boundarySizes(), cellSizes(), cellToBoundaryInterpolation(), deleteCells(), and setStaticGeometry().