/**
* IJKMesh class is defined by
- * - case 1: file name of mesh med file (MEDCouplingIMesh)
- * - case 2: 1D cartesian, xmin and xmax and number of cells
- * - case 3: 2D cartesian, xmin, xmax, ymin and ymax and numbers of cells in x direction and y direction
- * - case 4: 3D cartesian, xmin, xmax, ymin, ymax, zmin and zmax and numbers of cells in x direction, y direction and z direction
+ * - case 1: file name of mesh med file (MEDCouplingCMesh)
+ * - case 2: 1D cartesian, xmin and xmax and number of cells (MEDCouplingIMesh)
+ * - case 3: 2D cartesian, xmin, xmax, ymin and ymax and numbers of cells in x direction and y direction (MEDCouplingIMesh)
+ * - case 4: 3D cartesian, xmin, xmax, ymin, ymax, zmin and zmax and numbers of cells in x direction, y direction and z direction (MEDCouplingIMesh)
*
- * Mesh structure
- * nx , ny , nz cell structure
- * nx+1, ny+1, nz+1 node structure
+ * Mesh cell and node structures are stored in a single MEDCouplingStructuredMesh _mesh
+ * nx * ny * nz cell structure
+ * (nx+1) * (ny+1) * (nz+1) node structure
*
* The face structure is more tricky because it depends on the dimension
- * if dim=1, a single grid : nx+1
- * if dim=2, union of two grids : nx,ny+1 and nx,ny+1
- * if dim=3, union of three grids : nx,ny,nz+1, nx,ny+1,nz and nx+1,ny,nz
+ * if dim=1, a single grid :
+ * nx+1 nodes
+ * if dim=2, union of two grids :
+ * (nx+1)*ny nodes (faces orthogonal to x-axis), origin (0,dy/2)
+ * nx*(ny+1) nodes (faces orthogonal to x-axis), origin (dx/2,0)
+ * if dim=3, union of three grids :
+ * nx*ny*(nz+1) (faces orthogonal to x-axis), origin (0,dy/2,dz/2)
+ * nx*(ny+1)*nz (faces orthogonal to y-axis), origin (dx/2,0,dz/2)
+ * (nx+1)*ny*nz (faces orthogonal to z-axis), origin (dx/2,dy/2,0)
+ *
+ * Mesh face structures are stored in a vector of meshes _faceMeshes of size meshDim
+ * The face centers are located on the nodes of the meshes in _faceMeshes
+ * The origins of the meshes in _faceMeshes are shifted from the origin of _mesh by dx/2 on x-axis, dy/2 on y-axis and dz/2 on z-axis
*
* - number of nodes surounding each cell : known constant : 2*_Ndim
* - number of faces surounding each cell : known constant : 2 _Ndim
* - normal vectors surrounding each cell
* - measure of each cell : known constant : dx*dy*dz
+ * - measures of faces : known constants : dx*dy, dx*dz and dy*dz
*/
namespace MEDCoupling
class IJKCell;
class IJKFace;
-typedef enum
+enum EntityType
{
CELLS = 0,
NODES = 1,
FACES = 2,
- } EntityType;
+ };
#include <vector>
#include <string>
public: //----------------------------------------------------------------
/**
- * default constructor
+ * \brief default constructor
*/
Mesh ( void ) ;
/**
- * constructor with data to load a structured MEDCouplingIMesh
+ * \brief constructor with data to load a structured MEDCouplingIMesh
* @param filename name of structured mesh file
* @param meshLevel : relative mesh dimension : 0->cells, 1->Faces etc
*/
- Mesh ( const std::string filename, int meshLevel=0 ) ;
+ Mesh ( const std::string filename, const std::string & meshName="" , int meshLevel=0 ) ;
/**
- * constructor with data for a regular 1D grid
+ * \brief constructor with data for a regular 1D grid
* @param xmin : minimum x
* @param xmax : maximum x
* @param nx : Number of cells in x direction
Mesh( double xmin, double xmax, int nx, std::string meshName="MESH1D_Regular_Grid" ) ;
/**
- * constructor with data for a regular 2D grid
+ * \brief constructor with data for a regular 2D grid
* @param xmin : minimum x
* @param xmax : maximum x
* @param ymin : minimum y
Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, std::string meshName="MESH2D_Regular_Rectangle_Grid") ;
/**
- * constructor with data for a regular 3D grid
+ * \brief constructor with data for a regular 3D grid
* @param xmin : minimum x
* @param xmax : maximum x
* @param ymin : minimum y
*/
Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz, std::string meshName="MESH3D_Regular_Cuboid_Grid") ;
- Mesh( const MEDCoupling::MEDCouplingIMesh* mesh ) ;
+ Mesh( const MEDCoupling::MEDCouplingMesh* mesh ) ;
+ /**
+ * \brief Computes the global cell number from its IJK position
+ * @param int i : cell index along x-axis
+ * @param int j : cell index along y-axis
+ * @param int k : cell index along z-axis
+ * @return global cell number
+ */
+ int getCellNumber(int i, int j, int k) const { return _mesh->getCellIdFromPos ( i, j, k); };
+ /**
+ * \brief Computes the global node number from its IJK position
+ * @param int i : node index along x-axis
+ * @param int j : node index along y-axis
+ * @param int k : node index along z-axis
+ * @return global node number
+ */
+ int getNodeNumber(int i, int j, int k) const { return _mesh->getNodeIdFromPos ( i, j, k); };
+ /**
+ * \brief Computes the global face number from its IJK position
+ * @param int face grid number corresponding to the direction of the face normal : 0->x, 1->y, 2->z
+ * @param int i : node index along x-axis
+ * @param int j : node index along y-axis
+ * @param int k : node index along z-axis
+ * @return global node number
+ */
+ int getFaceNumber(int face_grid_number, int i, int j, int k) const { return _faceMeshes[face_grid_number]->getNodeIdFromPos ( i, j, k); };
+ /**
+ * \brief Computes the IJK position of a cell from its index number
+ * @param int cellId : global cell index
+ * @return vector of i,j,k indices of the cell
+ */
+ std::vector< int > getIJKCellCoordinates(int cellId) const { return _mesh->getLocationFromCellId(cellId); };
+ /**
+ * \brief Computes the IJK position of a node from its index number
+ * @param int nodeId : global node index
+ * @return vector of i,j,k indices of the node
+ */
+ std::vector< int > getIJKNodeCoordinates(int nodeId) const { return _mesh->getLocationFromNodeId(nodeId); };
+ /**
+ * \brief Computes the IJK position of a face from its index number
+ * @param int face grid number corresponding to the direction of the face normal : 0->x, 1->y, 2->z
+ * @param int faceId : global face index
+ * @return vector of i,j,k indices of the node
+ */
+ std::vector< int > getIJKFaceCoordinates( int face_grid_number, int faceId) const { return _faceMeshes[face_grid_number]->getLocationFromNodeId(faceId); };
+ /**
+ * \brief Computes the indices of nodes surrounding a given cell
+ * @param int cellId : global cell index
+ * @return vector of indices of the nodes surrounding the cell
+ */
+ std::vector< mcIdType > getNodeIdsOfCell(mcIdType cellId) const { std::vector< mcIdType > conn; _mesh_>getNodeIdsOfCell(mcIdType cellId, conn) ; return conn; };
+ /**
+ * \brief Computes the coordinates of a node
+ * @param int nodeId : global node index
+ * @return vector of components of the node coordinates
+ */
+ std::vector< double > getNodeCoordinates (mcIdType nodeId) const { std::vector< double > coo; _mesh_>getCoordinatesOfNode(mcIdType nodeId, coo) ; return coo; };
+ /**
+ * \brief Computes the coordinates of a face center of mass
+ * @param int face grid number corresponding to the direction of the face normal : 0->x, 1->y, 2->z
+ * @param int faceId : global face index
+ * @return vector of components of the face coordinates
+ */
+ std::vector< double > getFaceCoordinates ( int face_grid_number, mcIdType faceId) const { std::vector< double > coo; _faceMeshes[face_grid_number]>getCoordinatesOfNode(mcIdType faceId, coo) ; return coo; };
+ /**
+ * \brief Computes the isobarycenter of a cell
+ * @param int cellId : cell number
+ */
+ std::vector< double > getCellCenterCoordinates (mcIdType cellId) const ;
+
+ double getMeasureOfAnyCell () const;
+
/**
- * constructor with data
+ * \brief constructor with data
* @param filename : file name of structured mesh med file
* @param meshLevel : relative mesh dimension : 0->cells, 1->Faces etc
*/
- void readMeshMed( const std::string filename, int meshLevel=0 ) ;
+ void readMeshMed( const std::string filename, const std::string & meshName="" , int meshLevel=0 ) ;
/**
- * constructor by copy
+ * \brief constructor by copy
* @param mesh : The Mesh object to be copied
*/
Mesh ( const IJKMesh & mesh ) ;
/**
- * destructor
+ * \brief destructor
*/
~Mesh( void ) ;
/**
- * return mesh name
+ * \brief return mesh name
* @return _name
*/
- std::string getName( void ) const ;
+ std::string getName( void ) const { return _mesh->getName (); };
/**
- * return Space dimension
- * @return _spaceDim
+ * \brief return Space dimension
+ * @return spaceDim
*/
- int getSpaceDimension( void ) const ;
+ int getSpaceDimension( void ) const { return _mesh->getSpaceDimension (); };
/**
- * return Mesh dimension
- * @return _meshDim
+ * \brief return Mesh dimension
+ * @return meshDim
*/
- int getMeshDimension( void ) const ;
+ int getMeshDimension( void ) const { return _mesh->getMeshDimension (); };
/**
- * return the number of nodes in this mesh
+ * \brief return the number of nodes in this mesh
* @return _numberOfNodes
*/
int getNumberOfNodes ( void ) const ;
/**
- * return the number of faces in this mesh
+ * \brief return the number of faces in this mesh
* @return _numberOfFaces
*/
int getNumberOfFaces ( void ) const ;
/**
- * return the number of cells in this mesh
+ * \brief return the number of cells in this mesh
* @return _numberOfCells
*/
int getNumberOfCells ( void ) const ;
/**
- * return The cell i in this mesh
- * @return _cells[i]
- */
- IJKCell& getCell ( int i ) const ;
-
- /**
- * return The face i in this mesh
- * @return _faces[i]
- */
- IJKFace& getFace ( int i ) const ;
-
- /**
- * return The node i in this mesh
- * @return _nodes[i]
+ * \brief return the number of edges in this mesh
+ * @return _numberOfEdges
*/
- IJKNode& getNode ( int i ) const ;
+ int getNumberOfEdges ( void ) const ;
/**
- * return number of cell in x direction (structured mesh)
- * return _nX
+ * \brief return number of cell in x direction (structured mesh)
+ * @return _nX
*/
int getNx( void ) const ;
/**
- * return number of cell in y direction (structured mesh)
- * return _nY
+ * \brief return number of cell in y direction (structured mesh)
+ * @return _nY
*/
int getNy( void ) const ;
/**
- * return number of cell in z direction (structured mesh)
- * return _nZ
+ * \brief return number of cell in z direction (structured mesh)
+ * @return _nZ
*/
int getNz( void ) const ;
- double getXMin( void ) const ;
-
- double getXMax( void ) const ;
-
- double getYMin( void ) const ;
-
- double getYMax( void ) const ;
-
- double getZMin( void ) const ;
-
- double getZMax( void ) const ;
-
- std::vector<double> getDXYZ() const ;
-
std::vector<int> getCellGridStructure() const;
std::vector<int> getNodeGridStructure() const;
+ std::vector< vector< int > > getFaceGridStructures() const;
+
/**
- * surcharge operator =
+ * \brief overload operator =
* @param mesh : The Mesh object to be copied
*/
const Mesh& operator= ( const Mesh& mesh ) ;
/**
- * return the mesh MEDCoupling
- * return _mesh
+ * \brief return the MEDCouplingStructuredMesh
+ * @return _mesh
*/
- MEDCoupling::MCAuto<MEDCoupling::MEDCouplingIMesh> getMEDCouplingIMesh ( void ) const ;
+ MEDCoupling::MCAuto<MEDCoupling::MEDCouplingStructuredMesh> getMEDCouplingStructuredMesh ( void ) const ;
/**
- * return the list of face group names
- * return _faaceGroupNames
+ * \brief return the MEDCouplingStructuredMesh
+ * @return _faceMesh
+ */
+ std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDCouplingStructuredMesh> > getFaceMeshes const;
+
+ /**
+ * \brief return the list of face group names
+ * @return _faceGroupNames
*/
std::vector<std::string> getNameOfFaceGroups( void ) const ;
/**
- * return the list of node group names
- * return _nodeGroupNames
+ * \brief return the list of node group names
+ * @return _nodeGroupNames
*/
std::vector<std::string> getNameOfNodeGroups( void ) const ;
/**
- * write mesh in the VTK format
+ * \brief write the cell mesh in the VTK format
*/
void writeVTK ( const std::string fileName ) const ;
/**
- * write mesh in the MED format
+ * \brief write the cell and face meshes in the VTK format
*/
- void writeMED ( const std::string fileName ) const ;
+ void writeVTKAll meshes ( const std::string fileName ) const ;
+
+ /**
+ * \brief write the cell mesh in the MED format
+ */
+ void writeMED ( const std::string fileName, bool fromScratch = true ) const ;
+
+ /**
+ * \brief write the cell and face meshes in the MED format
+ */
+ void writeMEDAllMeshes ( const std::string fileName, bool fromScratch = true ) const ;
/*
* Functions to manage periodic boundary condition in square/cubic geometries
*/
- void setPeriodicFaces() ;
+ void setPeriodicFaces(bool check_groups= false, bool use_central_inversion=false) ;
int getIndexFacePeriodic(int indexFace, bool check_groups= false, bool use_central_inversion=false);
- void setBoundaryNodes();
+ void setBoundaryNodesFromFaces();
+ std::map<int,int> getIndexFacePeriodic( void ) const;
bool isIndexFacePeriodicSet() const ;
bool isBorderNode(int nodeid) const ;
bool isQuadrangular() const ;
bool isHexahedral() const ;
bool isStructured() const ;
- std::string getElementTypes() const ;
-
+
+ // epsilon used in mesh comparisons
+ double getComparisonEpsilon() const {return _epsilon;};
+ void setComparisonEpsilon(double epsilon){ _epsilon=epsilon;};
+ // Quick comparison of two meshes to see if they are identical with high probability (three cells are compared)
+ void checkFastEquivalWith( Mesh m) const { return _mesh()->checkFastEquivalWith(m.getMEDCouplingStructuredMesh(),_epsilon);};
+ // Deep comparison of two meshes to see if they are identical Except for their names and units
+ bool isEqualWithoutConsideringStr( Mesh m) const { return _mesh->isEqualWithoutConsideringStr(m.getMEDCouplingStructuredMesh(),_epsilon);};
+
+ std::vector< std::string > getElementTypesNames() const ;
/**
- * Compute the minimum value over all cells of the ratio cell perimeter/cell vaolume
+ * \brief Compute the minimum value over all cells of the ratio cell perimeter/cell volume
*/
- double minRatioVolSurf();
+ double minRatioVolSurf() const{ return _cellMeasure / *max_element(begin(_faceMeasures), end(_faceMeasures));};
/**
- * Return the maximum number of neighbours around an element (cells around a cell or nodes around a node)
+ * \brief Compute the maximum number of neighbours around an element (cells around a cell or nodes around a node)
*/
int getMaxNbNeighbours(EntityType type) const;
/**
- * return the measure of all cells (length in 1D, surface in 2D or volume in 3D)
- * @return _measureOfCells
+ * return the measure of a cell (length in 1D, surface in 2D or volume in 3D)
+ * @return _cellMeasure
*/
- double getCellsMeasure ( void ) const ;
+ double getCellMeasure ( ) const { return _cellMeasure;};
+
+ /**
+ * return the measure of a cell (length in 1D, surface in 2D or volume in 3D)
+ * @return _faceMeasures
+ */
+ std::vector< double > getFaceMeasures ( ) const { return _faceMeasures;};
/**
* return normal vectors around each cell
*/
- std::vector< Vector > getNormalVectors (void) const ;
+ std::vector< Vector > getNormalVectors (void) const { return _faceNormals;};
private: //----------------------------------------------------------------
- std::string _name;
-
- /**
- * Space dimension
- */
- int _spaceDim ;
-
- /**
- * Mesh dimension
+ /**
+ * \brief The cell mesh
+ * Question : can _mesh be const since no buildUnstructured is applied?
*/
- int _meshDim ;
-
- /*
- * Structured mesh parameters
- */
-
- double _xMin;
-
- double _xMax;
-
- double _yMin;
-
- double _yMax;
-
- double _zMin;
-
- double _zMax;
-
- std::vector<int> _nxyz;//Number of cells in each direction
-
- std::vector<double> _dxyz;//lenght depth and height of each cell
-
- /*
- * The number of nodes in this mesh.
+ MEDCoupling::MCAuto<MEDCoupling::MEDCouplingStructuredMesh> _mesh;//This is either a MEDCouplingIMesh (creation from scratch) or a MEDCouplingCMesh (loaded from a file)
+ /**
+ * \brief The face meshes
*/
- int _numberOfNodes;
-
- /*
- * The numbers of faces in this mesh.
+ std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDCouplingStructuredMesh> > _faceMeshes;//These are MEDCouplingIMesh
+ /**
+ * \brief Generate the face meshes from the cell mesh
*/
- int _numberOfFaces;
+ void setFaceMeshes();
- /*
- * The number of cells in this mesh.
- */
- int _numberOfCells;
+ double _cellMeasure;
+ std::vector< double > _faceMeasures;
+ std::vector< std::vector< double > > _faceNormals;
- /*
- * The names of face groups.
+ /* Boundary data */
+ /**
+ * \brief The names of face groups.
*/
std::vector<std::string> _faceGroupNames;
- /*
- * The names of node groups.
+ /**
+ * \brief The names of node groups.
*/
std::vector<std::string> _nodeGroupNames;
- /*
- * The list of face groups.
+ /**
+ * \brief The list of face groups.
*/
std::vector<MEDCoupling::MEDCouplingUMesh *> _faceGroups;
- /*
- * The list of node groups.
+ /**
+ * \brief The list of node groups.
*/
std::vector<MEDCoupling::DataArrayIdType *> _nodeGroups;
- /*
- * The mesh MEDCouplingIMesh
- */
- MEDCoupling::MCAuto<MEDCoupling::MEDCouplingIMesh> _mesh;
- std::vector< INTERP_KERNEL::NormalizedCellType > _eltsTypes;//List of cell types contained in the mesh
+
+ /**
+ * \brief List of boundary faces
+ */
+ std::vector< int > _boundaryFaceIds;
+ /**
+ * \brief List of boundary nodes
+ */
+ std::vector< int > _boundaryNodeIds;
- /*
- * Tools to manage periodic boundary conditions in square/cube geometries
+ /**
+ * \brief Tools to manage periodic boundary conditions in square/cube geometries
*/
bool _indexFacePeriodicSet;
std::map<int,int> _indexFacePeriodicMap;
- /* List of boundary faces*/
- std::vector< int > _boundaryFaceIds;
- /* List of boundary nodes*/
- std::vector< int > _boundaryNodeIds;
-
- /*
- * The number of nodes surrounding each cell.
- */
- int _numberOfNodesSurroundingCells ;
-
- /*
- * The number of faces surounding each cell.
- */
- int _numberOfFacesSurroundingCells ;
-
- /*
- * The length or surface or volume of each cell.
- */
- double _measureOfCells ;
-
- /*
- * The normal vectors surrounding each cell.
- */
- std::vector< Vector > _normalVectorsAroundCells ;
+ double _epsilon;
};
#endif /* IJKMESH_HXX_ */
* Authors: CDMATH
*/
+#include <iostream>
+#include <stdio.h>
+#include <stdlib.h>
+#include <iterator>
+#include <algorithm>
+
#include "IJKMesh.hxx"
#include "IJKNode.hxx"
#include "IJKCell.hxx"
#include "IJKFace.hxx"
+#include "MEDCouplingIMesh.hxx"
+#include "MEDCouplingUMesh.hxx"
#include "MEDFileCMesh.hxx"
#include "MEDLoader.hxx"
-#include "MEDCouplingIMesh.hxx"
#include "CdmathException.hxx"
-#include <iostream>
-#include <stdio.h>
-#include <stdlib.h>
-#include <iterator>
-#include <algorithm>
-
using namespace MEDCoupling;
using namespace std;
//----------------------------------------------------------------------
-Mesh::Mesh( void )
+IJKMesh::IJKMesh( void )
//----------------------------------------------------------------------
{
_mesh=NULL;
- _spaceDim = 0 ;
- _meshDim = 0 ;
- _numberOfNodes = 0;
- _numberOfFaces = 0;
- _numberOfCells = 0;
- _xMin=0.;
- _xMax=0.;
- _yMin=0.;
- _yMax=0.;
- _zMin=0.;
- _zMax=0.;
- _nxyz.resize(0);
- _dxyz.resize(0.);
+ _faceMeshes=std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDCouplingStructuredMesh> >(0)
+ _measureField=NULL;
_faceGroupNames.resize(0);
_faceGroups.resize(0);
_nodeGroupNames.resize(0);
_nodeGroups.resize(0);
_indexFacePeriodicSet=false;
_name="";
+ _epsilon=1e-6;
}
//----------------------------------------------------------------------
-Mesh::~Mesh( void )
+IJKMesh::~IJKMesh( void )
//----------------------------------------------------------------------
{
+ _measureField->decrRef();
}
std::string
-Mesh::getName( void ) const
+IJKMesh::getName( void ) const
{
return _name;
}
-Mesh::Mesh( const MEDCoupling::MEDCouplingIMesh* mesh )
+IJKMesh::IJKMesh( const MEDCoupling::MEDCouplingStructuredMesh* mesh )
{
- _spaceDim=mesh->getSpaceDimension();
- _meshDim=mesh->getMeshDimension();
- _dxyz=mesh->getDXYZ();
- _nxyz=mesh->getCellGridStructure();
- double* Box0=new double[2*_spaceDim];
- mesh->getBoundingBox(Box0);
- _name=mesh->getName();
+ _epsilon=1e-6;
_indexFacePeriodicSet=false;
-
- _xMin=Box0[0];
- _xMax=Box0[1];
- if (_spaceDim>=2)
- {
- _yMin=Box0[2];
- _yMax=Box0[3];
- }
- if (_spaceDim>=3)
- {
- _zMin=Box0[4];
- _zMax=Box0[5];
- }
+ _measureField = mesh->getMeasureField(true);
+
+ //_mesh=mesh if _mesh is declared const
+ _mesh=mesh->clone(false);//No deep copy : it is assumed node coordinates and cell connectivity will not change
+ setFaceMeshes();
+}
+
+void
+IJKMesh::setFaceMeshes()
+{
+ int meshDim = _mesh->getMeshDimension();
+
+ _faceMeshes=std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDCouplingStructuredMesh> >( meshDim )
+
+ std:vector< int > nodeStr = _mesh->getNodeGridStructure();
+ std:vector< int > dxyz = _mesh->getDXYZ();
+ std:vector< double > origin = _mesh->getOrigin();
- double *originPtr = new double[_spaceDim];
- double *dxyzPtr = new double[_spaceDim];
- int *nodeStrctPtr = new int[_spaceDim];
+ double originPtr[meshDim];
+ double dxyzPtr[meshDim];
+ mcIdType nodeStrctPtr[meshDim];
- for(int i=0;i<_spaceDim;i++)
+ /* Prepare the creation of face meshes, and the filling of face normals, face measures and cell measure */
+ for(int i=0; i<meshDim; i++)
{
- originPtr[i]=Box0[2*i];
- nodeStrctPtr[i]=_nxyz[i]+1;
- dxyzPtr[i]=_dxyz[i];
+ originPtr[i]=origin[i];
+ nodeStrctPtr[i]=nodeStr[i];
+ dxyzPtr[i]=dxyz[i];
+ }
+ _cellMeasure=1;
+ _faceNormals=std:vector< std:vector< double > > (meshDim, std:vector< double >(meshDim,0));
+ /* Creation of face meshes, and filling of face normals, face measures and cell measure */
+ for(int i=0; i<meshDim; i++)
+ {
+ _cellMeasure*=dxyz[i];
+ _faceMeasures[i]=1;
+ _faceNormals[i][i]=1;
+ for(int j=0; j<meshDim; j++)
+ if(j != i)
+ {
+ nodeStrctPtr[j]-=1;
+ originPtr[j]+=dxyz[j]/2;
+ _faceMeasures[i]*=dxyz[j];
+ }
+ _faceMesh[i]=MEDCouplingIMesh::New(_mesh->getName()+"_faces_"+std::to_string(i),
+ _mesh->getMeshDimension(),
+ nodeStrctPtr,
+ nodeStrctPtr+meshDim,
+ originPtr,
+ originPtr+meshDim,
+ dxyzPtr,
+ dxyzPtr+meshDim);
+ for(int j=0; j<meshDim; j++)
+ if(j != i)
+ {
+ nodeStrctPtr[j]+=1;
+ originPtr[j]-=dxyz[j]/2;
+ }
}
- _mesh=MEDCouplingIMesh::New(_name,
- _spaceDim,
- nodeStrctPtr,
- nodeStrctPtr+_spaceDim,
- originPtr,
- originPtr+_spaceDim,
- dxyzPtr,
- dxyzPtr+_spaceDim);
- delete [] originPtr;
- delete [] dxyzPtr;
- delete [] nodeStrctPtr;
- delete [] Box0 ;
}
//----------------------------------------------------------------------
-Mesh::Mesh( const IJKMesh& m )
+IJKMesh::IJKMesh( const IJKMesh& m )
//----------------------------------------------------------------------
{
- _spaceDim = m.getSpaceDimension() ;
- _meshDim = m.getMeshDimension() ;
- _name=m.getName();
- _xMax=m.getXMax();
- _yMin=m.getYMin();
- _yMax=m.getYMax();
- _zMin=m.getZMin();
- _zMax=m.getZMax();
- _nxyz = m.getCellGridStructure() ;
- _dxyz=m.getDXYZ();
-
- _numberOfNodes = m.getNumberOfNodes();
- _numberOfFaces = m.getNumberOfFaces();
- _numberOfCells = m.getNumberOfCells();
+ _measureField = m->getMeasureField(true);
+ _faceMeshes=m.getFaceMeshes();
_faceGroupNames = m.getNameOfFaceGroups() ;
_faceGroups = m.getFaceGroups() ;
_nodeGroupNames = m.getNameOfNodeGroups() ;
if(_indexFacePeriodicSet)
_indexFacePeriodicMap=m.getIndexFacePeriodic();
- MCAuto<MEDCouplingIMesh> m1=m.getMEDCouplingIMesh()->deepCopy();
- _mesh=m1;
+ //_mesh=m if _mesh is declared const
+ _mesh=m.getMEDCouplingIMesh()->clone(false);
}
//----------------------------------------------------------------------
-Mesh::Mesh( const std::string filename, int meshLevel )
+IJKMesh::IJKMesh( const std::string filename, const std::string & meshName, int meshLevel )
//----------------------------------------------------------------------
{
- readMeshMed(filename, meshLevel);
+ readMeshMed(filename, meshName, meshLevel);
}
//----------------------------------------------------------------------
void
-Mesh::readMeshMed( const std::string filename, const int meshLevel)
+IJKMesh::readMeshMed( const std::string filename, const std::string & meshName , int meshLevel )
//----------------------------------------------------------------------
{
- MEDFileCMesh *m=MEDFileCMesh::New(filename.c_str());//reads the first mesh encountered in the file, otherwise call New (const char *fileName, const char *mName, int dt=-1, int it=-1)
+ MEDFileCMesh * m;//Here we would like a MEDFileStructuredMesh but that class does not exist
+
+ if( meshName == "" )
+ m=MEDFileCMesh::New(filename.c_str());//reads the first mesh encountered in the file, otherwise call New (const char *fileName, const char *mName, int dt=-1, int it=-1)
+ else
+ m=MEDFileCMesh::New(filename.c_str(), meshName.c_str());//seeks the mesh named meshName in the file
+
_mesh=m->getMeshAtLevel(meshLevel);
_mesh->checkConsistencyLight();
_mesh->setName(_mesh->getName());
- _meshDim=_mesh->getMeshDimension();
- _spaceDim=_mesh->getSpaceDimension();
- _name=_mesh->getName();
+ _epsilon=1e-6;
_indexFacePeriodicSet=false;
- MEDCoupling::MEDCouplingIMesh* structuredMesh = dynamic_cast<MEDCoupling::MEDCouplingIMesh*> (_mesh.retn());
- if(structuredMesh)
- {
- _dxyz=structuredMesh->getDXYZ();
- _nxyz=structuredMesh->getCellGridStructure();
- double* Box0=new double[2*_spaceDim];
- structuredMesh->getBoundingBox(Box0);
-
- _xMin=Box0[0];
- _xMax=Box0[1];
- std::cout<<"nx= "<<_nxyz[0];
- if (_spaceDim>=2)
- {
- _yMin=Box0[2];
- _yMax=Box0[3];
- std::cout<<", "<<"ny= "<<_nxyz[1];
- }
- if (_spaceDim>=3)
- {
- _zMin=Box0[4];
- _zMax=Box0[5];
- std::cout<<", "<<"nz= "<<_nxyz[2];
- }
- }
- else
- throw CdmathException("Mesh::readMeshMed med file does not contain a structured MedcouplingIMesh mesh");
+ _measureField = mesh->getMeasureField(true);
+
+ setFaceMeshes();
cout<<endl<< "Loaded file "<< filename<<endl;
- cout<<"Structured Mesh name= "<<m->getName()<<", mesh dim="<< _meshDim<< ", space dim="<< _spaceDim<< ", nb cells= "<<getNumberOfCells()<< ", nb nodes= "<<getNumberOfNodes()<<endl;
+ cout<<"Structured Mesh name= "<<_mesh->getName()<<", mesh dim="<< _mesh->getMeshDimension()<< ", space dim="<< _mesh->getSpaceDimension()<< ", nb cells= "<<_mesh->getNumberOfCells()<< ", nb nodes= "<<_mesh->getNumberOfNodes()<<endl;
m->decrRef();
}
-void
-Mesh::setPeriodicFaces()
-{
- _indexFacePeriodicSet=true;
-}
-
-bool
-Mesh::isBorderNode(int nodeid) const
-{
- return getNode(nodeid).isBorder();
-}
-
-bool
-Mesh::isTriangular() const
-{
- return false;
-}
-
-bool
-Mesh::isQuadrangular() const
-{
- return _meshDim==2;
-}
-bool
-Mesh::isHexahedral() const
-{
- return _meshDim==3;
-}
-bool
-Mesh::isStructured() const
-{
- return true;
-}
-
-std::string
-Mesh::getElementTypes() const
-{
- if( _meshDim==1 )
- return "Segments ";
- else if( _meshDim==2 )
- return "Quadrangles ";
- else if( _meshDim==3 )
- return "Hexahedra ";
- else
- {
- cout<< "Mesh " + _name + " does not have acceptable dimension. Dimension is " << _meshDim<<endl;
- throw CdmathException("Mesh::getElementTypes : wrong dimension");
- }
-}
-
//----------------------------------------------------------------------
-Mesh::Mesh( double xmin, double xmax, int nx, std::string meshName )
+IJKMesh::IJKMesh( double xmin, double xmax, int nx, std::string meshName )
//----------------------------------------------------------------------
{
if(nx<=0)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx) : nx <= 0");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx) : nx <= 0");
if(xmin>=xmax)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx) : xmin >= xmax");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx) : xmin >= xmax");
double dx = (xmax - xmin)/nx ;
- _spaceDim = 1 ;
- _meshDim = 1 ;
- _name=meshName;
+ double meshDim = 1 ;
+ _epsilon=1e-6;
+ _isStructured = true;
_indexFacePeriodicSet=false;
- _xMin=xmin;
- _xMax=xmax;
- _yMin=0.;
- _yMax=0.;
- _zMin=0.;
- _zMax=0.;
-
- _dxyz.resize(_spaceDim);
- _dxyz[0]=dx;
- _nxyz.resize(_spaceDim);
- _nxyz[0]=nx;
-
- double *originPtr = new double[_spaceDim];
- double *dxyzPtr = new double[_spaceDim];
- int *nodeStrctPtr = new int[_spaceDim];
+ double originPtr[meshDim];
+ double dxyzPtr[meshDim];
+ mcIdType nodeStrctPtr[meshDim];
originPtr[0]=xmin;
nodeStrctPtr[0]=nx+1;
dxyzPtr[0]=dx;
- _mesh=MEDCouplingIMesh::New(meshName,
- _spaceDim,
+ _mesh=MEDCouplingIIJKMesh::New(meshName,
+ meshDim,
nodeStrctPtr,
- nodeStrctPtr+_spaceDim,
+ nodeStrctPtr+meshDim,
originPtr,
- originPtr+_spaceDim,
+ originPtr+meshDim,
dxyzPtr,
- dxyzPtr+_spaceDim);
- delete [] originPtr;
- delete [] dxyzPtr;
- delete [] nodeStrctPtr;
-
- _numberOfCells = _mesh->getNumberOfCells() ;
-
- _numberOfNodes = _mesh->getNumberOfNodes() ;
-
- _numberOfFaces = _numberOfNodes;
-
+ dxyzPtr+meshDim);
+ _measureField = _mesh->getMeasureField(true);
+ setFaceMeshes();
}
//----------------------------------------------------------------------
-Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, std::string meshName)
+IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, std::string meshName)
//----------------------------------------------------------------------
{
if(nx<=0 || ny<=0)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny) : nx <= 0 or ny <= 0");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny) : nx <= 0 or ny <= 0");
if(xmin>=xmax)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny) : xmin >= xmax");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny) : xmin >= xmax");
if(ymin>=ymax)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny) : ymin >= ymax");
-
- _xMin=xmin;
- _xMax=xmax;
- _yMin=ymin;
- _yMax=ymax;
- _zMin=0.;
- _zMax=0.;
-
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny) : ymin >= ymax");
double dx = (xmax - xmin)/nx ;
double dy = (ymax - ymin)/ny ;
- _spaceDim = 2 ;
- _meshDim = 2 ;
- _name=meshName;
+ double meshDim = 2 ;
+ _epsilon=1e-6;
_indexFacePeriodicSet=false;
- _nxyz.resize(_spaceDim);
- _nxyz[0]=nx;
- _nxyz[1]=ny;
-
- _dxyz.resize(_spaceDim);
- _dxyz[0]=dx;
- _dxyz[1]=dy;
- double *originPtr = new double[_spaceDim];
- double *dxyzPtr = new double[_spaceDim];
- int *nodeStrctPtr = new int[_spaceDim];
+ double originPtr[meshDim];
+ double dxyzPtr[meshDim];
+ mcIdType nodeStrctPtr[meshDim];
originPtr[0]=xmin;
originPtr[1]=ymin;
dxyzPtr[1]=dy;
_mesh=MEDCouplingIMesh::New(meshName,
- _spaceDim,
+ meshDim,
nodeStrctPtr,
- nodeStrctPtr+_spaceDim,
+ nodeStrctPtr+meshDim,
originPtr,
- originPtr+_spaceDim,
+ originPtr+meshDim,
dxyzPtr,
- dxyzPtr+_spaceDim);
-
- delete [] originPtr;
- delete [] dxyzPtr;
- delete [] nodeStrctPtr;
-
- _numberOfCells = _mesh->getNumberOfCells() ;
-
- _numberOfNodes = _mesh->getNumberOfNodes() ;
-
- _numberOfFaces = nx*(ny+1)+ny*(nx+1);
-
+ dxyzPtr+meshDim);
+ _measureField = _mesh->getMeasureField(true);
+ setFaceMeshes();
}
//----------------------------------------------------------------------
-Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz, std::string meshName)
+IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz, std::string meshName)
//----------------------------------------------------------------------
{
if(nx<=0 || ny<=0 || nz<=0)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : nx <= 0 or ny <= 0 or nz <= 0");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : nx <= 0 or ny <= 0 or nz <= 0");
if(xmin>=xmax)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : xmin >= xmax");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : xmin >= xmax");
if(ymin>=ymax)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : ymin >= ymax");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : ymin >= ymax");
if(zmin>=zmax)
- throw CdmathException("Mesh::Mesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : zmin >= zmax");
+ throw CdmathException("IJKMesh::IJKMesh( double xmin, double xmax, int nx, double ymin, double ymax, int ny, double zmin, double zmax, int nz) : zmin >= zmax");
- _spaceDim = 3;
- _meshDim = 3;
- _name=meshName;
- _indexFacePeriodicSet=false;
- _xMin=xmin;
- _xMax=xmax;
- _yMin=ymin;
- _yMax=ymax;
- _zMin=zmin;
- _zMax=zmax;
+ double meshDim = 3;
+ _epsilon=1e-6;
double dx = (xmax - xmin)/nx ;
double dy = (ymax - ymin)/ny ;
double dz = (zmax - zmin)/nz ;
- _dxyz.resize(_spaceDim);
- _dxyz[0]=dx;
- _dxyz[1]=dy;
- _dxyz[2]=dz;
-
- _nxyz.resize(_spaceDim);
- _nxyz[0]=nx;
- _nxyz[1]=ny;
- _nxyz[2]=nz;
-
- double *originPtr = new double[_spaceDim];
- double *dxyzPtr = new double[_spaceDim];
- int *nodeStrctPtr = new int[_spaceDim];
+ double originPtr[meshDim];
+ double dxyzPtr[meshDim];
+ mcIdType nodeStrctPtr[meshDim];
originPtr[0]=xmin;
originPtr[1]=ymin;
dxyzPtr[2]=dz;
_mesh=MEDCouplingIMesh::New(meshName,
- _spaceDim,
+ meshDim,
nodeStrctPtr,
- nodeStrctPtr+_spaceDim,
+ nodeStrctPtr+meshDim,
originPtr,
- originPtr+_spaceDim,
+ originPtr+meshDim,
dxyzPtr,
- dxyzPtr+_spaceDim);
-
- delete [] originPtr;
- delete [] dxyzPtr;
- delete [] nodeStrctPtr;
-
- _numberOfCells = _mesh->getNumberOfCells() ;
-
- _numberOfNodes = _mesh->getNumberOfNodes() ;
-
- _numberOfFaces = nx*ny*(nz+1)+nx*nz*(ny+1)+ny*nz*(nx+1);
-
+ dxyzPtr+meshDim);
+ _measureField = _mesh->getMeasureField(true);
+ setFaceMeshes();
}
-//----------------------------------------------------------------------
-int
-Mesh::getSpaceDimension( void ) const
-//----------------------------------------------------------------------
+void
+IJKMesh::setPeriodicFaces()
{
- return _spaceDim ;
+ _indexFacePeriodicSet=true;
}
-//----------------------------------------------------------------------
-int
-Mesh::getMeshDimension( void ) const
-//----------------------------------------------------------------------
+bool
+IJKMesh::isBorderNode(int nodeid) const
{
- return _meshDim ;
+ return getNode(nodeid).isBorder();
}
-std::vector<double>
-Mesh::getDXYZ() const
+bool
+IJKMesh::isTriangular() const
{
- return _dxyz;
+ return false;
}
-std::vector<int>
-Mesh::getCellGridStructure() const
+bool
+IJKMesh::isQuadrangular() const
{
- return _nxyz;
+ return _mesh->getMeshDimension()==2;
}
-
-//----------------------------------------------------------------------
-int
-Mesh::getNx( void ) const
-//----------------------------------------------------------------------
+bool
+IJKMesh::isHexahedral() const
{
- return _nxyz[0];
+ return _mesh->getMeshDimension()==3;
}
-
-//----------------------------------------------------------------------
-int
-Mesh::getNy( void ) const
-//----------------------------------------------------------------------
+bool
+IJKMesh::isStructured() const
{
- if(_spaceDim < 2)
- throw CdmathException("int double& Field::operator[ielem] : Ny is not defined in dimension < 2!");
-
- return _nxyz[1];
+ return true;
}
-//----------------------------------------------------------------------
-int
-Mesh::getNz( void ) const
-//----------------------------------------------------------------------
+std::string
+IJKMesh::getElementTypes() const
{
- if(_spaceDim < 3)
- throw CdmathException("int Mesh::getNz( void ) : Nz is not defined in dimension < 3!");
-
- return _nxyz[2];
+ if( _mesh->getMeshDimension()==1 )
+ return "Segments ";
+ else if( _mesh->getMeshDimension()==2 )
+ return "Quadrangles ";
+ else if( _mesh->getMeshDimension()==3 )
+ return "Hexahedra ";
+ else
+ {
+ cout<< "Mesh " + _name + " does not have acceptable dimension. Mesh dimension is " << meshDim<<endl;
+ throw CdmathException("IJKMesh::getElementTypes : wrong dimension");
+ }
}
//----------------------------------------------------------------------
double
-Mesh::getXMin( void ) const
+IJKMesh::getXMin( void ) const
//----------------------------------------------------------------------
-{
- return _xMin ;
+{
+ double Box0[2*_mesh->getMeshDimension()];
+ _mesh->getBoundingBox(Box0);
+
+ return Box0[0] ;
}
//----------------------------------------------------------------------
double
-Mesh::getXMax( void ) const
+IJKMesh::getXMax( void ) const
//----------------------------------------------------------------------
{
- return _xMax ;
+ double Box0[2*_mesh->getMeshDimension()];
+ _mesh->getBoundingBox(Box0);
+
+ return Box0[1] ;
}
//----------------------------------------------------------------------
double
-Mesh::getYMin( void ) const
+IJKMesh::getYMin( void ) const
//----------------------------------------------------------------------
{
- return _yMin ;
+ if(_mesh->getMeshDimension()<2)
+ throw CdmathException("IJKMesh::getYMin : dimension should be >=2");
+
+ double Box0[2*_mesh->getMeshDimension()];
+ _mesh->getBoundingBox(Box0);
+
+ return Box0[2] ;
}
//----------------------------------------------------------------------
double
-Mesh::getYMax( void ) const
+IJKMesh::getYMax( void ) const
//----------------------------------------------------------------------
{
- return _yMax ;
+ if(_mesh->getMeshDimension()<2)
+ throw CdmathException("IJKMesh::getYMax : dimension should be >=2");
+
+ double Box0[2*_mesh->getMeshDimension()];
+ _mesh->getBoundingBox(Box0);
+
+ return Box0[3] ;
}
//----------------------------------------------------------------------
double
-Mesh::getZMin( void ) const
+IJKMesh::getZMin( void ) const
//----------------------------------------------------------------------
{
- return _zMin ;
+ if(_mesh->getMeshDimension()<3)
+ throw CdmathException("IJKMesh::getZMin : dimension should be 3");
+
+ double Box0[2*_mesh->getMeshDimension()];
+ _mesh->getBoundingBox(Box0);
+
+ return Box0[4] ;
}
//----------------------------------------------------------------------
double
-Mesh::getZMax( void ) const
+IJKMesh::getZMax( void ) const
//----------------------------------------------------------------------
{
- return _zMax ;
+ if(_mesh->getMeshDimension()<3)
+ throw CdmathException("IJKMesh::getZMax : dimension should be 3");
+
+ double Box0[2*_mesh->getMeshDimension()];
+ _mesh->getBoundingBox(Box0);
+
+ return Box0[5] ;
}
//----------------------------------------------------------------------
-MCAuto<MEDCouplingIMesh>
-Mesh::getMEDCouplingIMesh( void ) const
+MCAuto<MEDCouplingStructuredMesh>
+IJKMesh::getMEDCouplingStructuredMesh( void ) const
//----------------------------------------------------------------------
{
return _mesh ;
}
+std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDCouplingStructuredMesh> >
+IJKMesh::getFaceMeshes() const
+{
+ return _faceMeshes;
+}
//----------------------------------------------------------------------
int
-Mesh::getNumberOfNodes ( void ) const
+IJKMesh::getNumberOfNodes ( void ) const
//----------------------------------------------------------------------
{
- return _numberOfNodes ;
+ return _mesh->getNumberOfNodes ();
}
//----------------------------------------------------------------------
int
-Mesh::getNumberOfCells ( void ) const
+IJKMesh::getNumberOfCells ( void ) const
//----------------------------------------------------------------------
{
- return _numberOfCells ;
+ return _mesh->getNumberOfCells ();
}
-//----------------------------------------------------------------------
-int
-Mesh::getNumberOfFaces ( void ) const
-//----------------------------------------------------------------------
+std::vector<int>
+IJKMesh::getCellGridStructure() const
{
- return _numberOfFaces ;
+ return _mesh->getCellGridStructure();
}
-//----------------------------------------------------------------------
-Cell&
-Mesh::getCell ( int i ) const
-//----------------------------------------------------------------------
+std::vector<int>
+IJKMesh::getNodeGridStructure() const
{
- return _cells[i] ;
+ return _mesh->getNodeGridStructure();
}
-//----------------------------------------------------------------------
-Face&
-Mesh::getFace ( int i ) const
-//----------------------------------------------------------------------
+std::vector< vector<int> >
+IJKMesh::getFaceGridStructures() const
{
- return _faces[i] ;
+ std::vector< vector<int> > result(_mesh->getmeshDimension());
+ for(int i = 0; i<_mesh->getmeshDimension(); i++)
+ result[i] = _faceMeshes[i]->getNodeGridStructure();
+
+ return result;
}
//----------------------------------------------------------------------
-Node&
-Mesh::getNode ( int i ) const
+int
+IJKMesh::getNumberOfFaces ( void ) const
//----------------------------------------------------------------------
{
- return _nodes[i] ;
+ std::vector<int> NxNyNz = _mesh->getCellGridStructure();
+
+ switch( _mesh->getMeshDimension () )
+ {
+ case 1:
+ return NxNyNz[0]+1;
+ case 2:
+ return NxNyNz[0]*(1+NxNyNz[1]) + NxNyNz[1]*(1 + NxNyNz[0]);
+ case 2:
+ return NxNyNz[0]*NxNyNz[1]*(1+NxNyNz[2]) + NxNyNz[0]*NxNyNz[2]*(1+NxNyNz[1]) + NxNyNz[1]*NxNyNz[2]*(1+NxNyNz[0]);
+ default
+ throw CdmathException("IJKMesh::getNumberOfFaces space dimension must be between 1 and 3");
+ }
}
+int
+IJKMesh::getNumberOfEdges ( void ) const
+{
+ std::vector<int> NxNyNz = _mesh->getCellGridStructure();
+
+ switch( _mesh->getMeshDimension () )
+ {
+ case 1:
+ return NxNyNz[0];
+ case 2:
+ return NxNyNz[0]*(1+NxNyNz[1]) + NxNyNz[1]*(1 + NxNyNz[0]);
+ case 2:
+ return NxNyNz[0]*(1+NxNyNz[1])*(1+NxNyNz[2]) + NxNyNz[1]*(1 + NxNyNz[0])*(1+NxNyNz[2]) + NxNyNz[2]*(1 + NxNyNz[0])*(1+NxNyNz[1]);
+ default
+ throw CdmathException("IJKMesh::getNumberOfEdges space dimension must be between 1 and 3");
+ }
+};
+
vector<string>
-Mesh::getNameOfFaceGroups( void ) const
+IJKMesh::getNameOfFaceGroups( void ) const
{
return _faceGroupNames;
}
vector<MEDCoupling::MEDCouplingIMesh *>
-Mesh::getFaceGroups( void ) const
+IJKMesh::getFaceGroups( void ) const
{
return _faceGroups;
}
vector<string>
-Mesh::getNameOfNodeGroups( void ) const
+IJKMesh::getNameOfNodeGroups( void ) const
{
return _nodeGroupNames;
}
vector<MEDCoupling::DataArrayIdType *>
-Mesh::getNodeGroups( void ) const
+IJKMesh::getNodeGroups( void ) const
{
return _nodeGroups;
}
//----------------------------------------------------------------------
const IJKMesh&
-Mesh::operator= ( const IJKMesh& mesh )
+IJKMesh::operator= ( const IJKMesh& mesh )
//----------------------------------------------------------------------
{
- _spaceDim = mesh.getSpaceDimension() ;
- _meshDim = mesh.getMeshDimension() ;
- _name = mesh.getName();
- _numberOfNodes = mesh.getNumberOfNodes();
- _numberOfFaces = mesh.getNumberOfFaces();
- _numberOfCells = mesh.getNumberOfCells();
+ _epsilon=mesh.getComparisonEpsilon();
_indexFacePeriodicSet= mesh.isIndexFacePeriodicSet();
if(_indexFacePeriodicSet)
_indexFacePeriodicMap=mesh.getIndexFacePeriodic();
- _nxyz = mesh.getCellGridStructure() ;
- _dxyz=mesh.getDXYZ();
- _xMin=mesh.getXMin();
- _xMax=mesh.getXMax();
- _yMin=mesh.getYMin();
- _yMax=mesh.getYMax();
- _zMin=mesh.getZMin();
- _zMax=mesh.getZMax();
-
_faceGroupNames = mesh.getNameOfFaceGroups() ;
_faceGroups = mesh.getFaceGroups() ;
_nodeGroupNames = mesh.getNameOfNodeGroups() ;
_nodeGroups = mesh.getNodeGroups() ;
- MCAuto<MEDCouplingIMesh> m1=mesh.getMEDCouplingIMesh()->deepCopy();
- _mesh=m1;
+ //_mesh=mesh if _mesh is declared const
+ _mesh=mesh.getMEDCouplingStructuredMesh()->clone(false);
+
+ _faceMeshes=mesh.getFaceMeshes());
return *this;
}
-
-bool Mesh::isIndexFacePeriodicSet() const
+
+bool IJKMesh::isIndexFacePeriodicSet() const
{
return _indexFacePeriodicSet;
}
//----------------------------------------------------------------------
double
-Mesh::minRatioVolSurf()
+IJKMesh::minRatioVolSurf()
{
return dx_min;
}
int
-Mesh::getMaxNbNeighbours(EntityType type) const
+IJKMesh::getMaxNbNeighbours(EntityType type) const
{
- return 2*_meshDim;
+ return 2*_mesh->getMeshDimension();
}
//----------------------------------------------------------------------
void
-Mesh::writeVTK ( const std::string fileName ) const
+IJKMesh::writeVTK ( const std::string fileName ) const
//----------------------------------------------------------------------
{
string fname=fileName+".vtu";
_mesh->writeVTK(fname.c_str()) ;
}
+//----------------------------------------------------------------------
+void
+IJKMesh::writeVTKAllMeshes ( const std::string fileName ) const
+//----------------------------------------------------------------------
+{
+ string fname=fileName+".vtu";
+ _mesh->writeVTK(fname.c_str()) ;
+ for(int i=0; i< _mesh->getMeshDimension(); i++)
+ _faceMeshes[i]->writeVTK(fname.c_str()) ;
+}
+
+//----------------------------------------------------------------------
+void
+IJKMesh::writeMED ( const std::string fileName, bool fromScratch ) const
+//----------------------------------------------------------------------
+{
+ string fname=fileName+".med";
+ //Save cell mesh after checking mesh is imesh
+ const MEDCoupling::MEDCouplingIMesh* iMesh = dynamic_cast< const MEDCoupling::MEDCouplingIMesh* > ((const MEDCoupling::MEDCouplingStructuredMesh*) _mesh);
+ if(iMesh)//medcouplingimesh : Use convertToCartesian in order to write mesh
+ MEDCoupling::WriteMesh(fname.c_str(),iMesh->convertToCartesian(), fromScratch);
+ else//medcouplingcmesh : save directly
+ MEDCoupling::WriteMesh(fname.c_str(),_mesh, fromScratch);
+}
//----------------------------------------------------------------------
void
-Mesh::writeMED ( const std::string fileName ) const
+IJKMesh::writeMEDAllMeshes ( const std::string fileName, bool fromScratch ) const
//----------------------------------------------------------------------
{
+ //Save cell mesh after checking mesh is imesh
string fname=fileName+".med";
- MEDCoupling::WriteCMesh(fname.c_str(),_mesh,true);
+ const MEDCoupling::MEDCouplingIMesh* iMesh = dynamic_cast< const MEDCoupling::MEDCouplingIMesh* > ((const MEDCoupling::MEDCouplingStructuredMesh*) _mesh);
+ if(iMesh)//medcouplingimesh : Use convertToCartesian in order to write mesh
+ MEDCoupling::WriteMesh(fname.c_str(),iMesh->convertToCartesian(), fromScratch);
+ else//medcouplingcmesh : save directly
+ MEDCoupling::WriteMesh(fname.c_str(),_mesh, fromScratch);
+
+ //Save face meshes after checking mesh is imesh
+ std::vector< const MEDCoupling::MEDCouplingIMesh* > iMeshes(_mesh->getMeshDimension());
+ for(int i=0; i< _mesh->getMeshDimension(); i++)
+ {
+ const MEDCoupling::MEDCouplingIMesh* iMeshes[i] = dynamic_cast< const MEDCoupling::MEDCouplingIMesh* > ((const MEDCoupling::MEDCouplingStructuredMesh*) _faceMeshes[i]);
+ if(iMesh)//medcouplingimesh : Use convertToCartesian in order to write mesh
+ MEDCoupling::WriteMesh(fname.c_str(),iMeshes[i]->convertToCartesian(), fromScratch);
+ else//medcouplingcmesh : save directly
+ MEDCoupling::WriteMesh(fname.c_str(), _faceMeshes[i], fromScratch);
+ }
+}
+
+std::vector< double >
+IJKMesh::getCellCenterCoordinates (mcIdType cellId) const
+{
+ std::vector< double > result(_mesh->getMeshDimension(),0), coo_node;
+ std::vector< mcIdType > conn=getNodeIdsOfCell(mcIdType cellId);
+ int nbNodes=conn.size();
+
+ for (int i = 0; i< nbNodes; i++)
+ {
+ coo_node = getNodeCoordinates (conn[i]);
+ for(int j=0; j< _mesh->getMeshDimension(); j++)
+ result[j]+=coo_node[j];
+ }
+ for(int j=0; j< _mesh->getMeshDimension(); j++)
+ result[j]/=nbNodes;
+ return result;
+};
+
+//----------------------------------------------------------------------
+int
+IJKMesh::getNx( void ) const
+//----------------------------------------------------------------------
+{
+ return _mesh->getCellGridStructure()[0];
+}
- mu->decrRef();
+//----------------------------------------------------------------------
+int
+IJKMesh::getNy( void ) const
+//----------------------------------------------------------------------
+{
+ if(_mesh->getMeshDimension () < 2)
+ throw CdmathException("int IJKMesh::getNy( void ) : Ny is not defined in dimension < 2!");
+ else
+ return _mesh->getCellGridStructure()[1];
+}
+
+//----------------------------------------------------------------------
+int
+IJKMesh::getNz( void ) const
+//----------------------------------------------------------------------
+{
+ if(_mesh->getMeshDimension () < 3)
+ throw CdmathException("int IJKMesh::getNz( void ) : Nz is not defined in dimension < 3!");
+ else
+ return _mesh->getCellGridStructure()[2];
}
