-using namespace std;
/*
File Mesh.cxx
$Header$
#include <list>
#include <map>
+#include <sstream>
+#include "MEDMEM_DriversDef.hxx"
#include "MEDMEM_Field.hxx"
#include "MEDMEM_Mesh.hxx"
+#include "MEDMEM_Support.hxx"
#include "MEDMEM_Family.hxx"
#include "MEDMEM_Group.hxx"
+#include "MEDMEM_Coordinate.hxx"
#include "MEDMEM_Connectivity.hxx"
#include "MEDMEM_CellModel.hxx"
+#include "MEDMEM_DriverFactory.hxx"
+
+using namespace std;
+using namespace MEDMEM;
+using namespace MED_EN;
+
+//#include "MEDMEM_Grid.hxx" this inclision should have never be here !!!
+
//update Families with content list
//int family_count(int family_number, int count, int * entities_number, int * entities_list) ;
// ------- Drivers Management Part
-// MESH::INSTANCE_DE<MED_MESH_RDONLY_DRIVER> MESH::inst_med_rdonly ;
-//const MESH::INSTANCE * const MESH::instances[] = { &MESH::inst_med_rdonly , &MESH::inst_med_rdwr } ;
-
-MESH::INSTANCE_DE<MED_MESH_RDWR_DRIVER> MESH::inst_med ;
-const MESH::INSTANCE * const MESH::instances[] = { &MESH::inst_med } ;
-
-/*! Add a MESH driver of type (MED_DRIVER, ....) associated with file <fileName>. The meshname used in the file
- is <driverName>. addDriver returns an int handler. */
+/*! Add a %MESH driver of type %driverTypes (MED_DRIVER, ....) associated with file fileName. The meshname used in the file
+ is driverName. addDriver returns an integer handler. */
int MESH::addDriver(driverTypes driverType,
- const string & fileName/*="Default File Name.med"*/,const string & driverName/*="Default Mesh Name"*/) {
+ const string & fileName/*="Default File Name.med"*/,
+ const string & driverName/*="Default Mesh Name"*/,
+ med_mode_acces access) {
- const char * LOC = "MESH::addDriver(driverTypes driverType, const string & fileName=\"Default File Name.med\",const string & driverName=\"Default Mesh Name\") : ";
+ const char * LOC = "MESH::addDriver(driverTypes driverType, const string & fileName=\"Default File Name.med\",const string & driverName=\"Default Mesh Name\",MED_EN::med_mode_acces access) : ";
GENDRIVER * driver;
- int current;
BEGIN_OF(LOC);
- driver = instances[driverType]->run(fileName, this) ;
+ SCRUTE(driverType);
+
+ driver = DRIVERFACTORY::buildDriverForMesh(driverType,fileName,this,
+ driverName,access) ;
+
_drivers.push_back(driver);
- current = _drivers.size()-1;
+
+ int current = _drivers.size()-1;
_drivers[current]->setMeshName(driverName);
- return current;
END_OF(LOC);
+ return current;
}
/*! Add an existing MESH driver. */
-int MESH::addDriver(MED_MESH_DRIVER & driver) {
+int MESH::addDriver(GENDRIVER & driver) {
const char * LOC = "MESH::addDriver(GENDRIVER &) : ";
BEGIN_OF(LOC);
// A faire : Vérifier que le driver est de type MESH.
- _drivers.push_back(&driver);
+ GENDRIVER * newDriver = driver.copy() ;
+
+ _drivers.push_back(newDriver);
return _drivers.size()-1;
-
+
END_OF(LOC);
}
void MESH::init() {
- string _name = "NOT DEFINED"; // A POSITIONNER EN FCT DES IOS ?
+ const char * LOC = "MESH::init(): ";
- _numberOfMEDNodeFamily = MED_INVALID;
- _MEDArrayNodeFamily = (int * ) NULL; // SOLUTION TEMPORAIRE
- _numberOfMEDCellFamily = (int * ) NULL;
- _numberOfMEDFaceFamily = (int * ) NULL;
- _numberOfMEDEdgeFamily = (int * ) NULL;
- _MEDArrayCellFamily = (int **) NULL; // SOLUTION TEMPORAIRE
- _MEDArrayFaceFamily = (int **) NULL; // SOLUTION TEMPORAIRE
- _MEDArrayEdgeFamily = (int **) NULL; // SOLUTION TEMPORAIRE
-
- COORDINATE * _coordinate = (COORDINATE *) NULL;
- CONNECTIVITY * _connectivity = (CONNECTIVITY *) NULL;
+ BEGIN_OF(LOC);
+
+ _name = "NOT DEFINED"; // A POSITIONNER EN FCT DES IOS ?
+
+ _coordinate = (COORDINATE *) NULL;
+ _connectivity = (CONNECTIVITY *) NULL;
_spaceDimension = MED_INVALID; // 0 ?!?
_meshDimension = MED_INVALID;
_numberOfNodes = MED_INVALID;
-
- _numberOfNodesFamilies = 0; // MED_INVALID ?!?
- _numberOfCellsFamilies = 0;
- _numberOfFacesFamilies = 0;
- _numberOfEdgesFamilies = 0;
- _numberOfNodesGroups = 0; // MED_INVALID ?!?
- _numberOfCellsGroups = 0;
- _numberOfFacesGroups = 0;
- _numberOfEdgesGroups = 0;
+ _isAGrid = false;
+
+ _arePresentOptionnalNodesNumbers = 0;
+ END_OF(LOC);
};
/*! Create an empty MESH. */
-MESH::MESH() {
+MESH::MESH():_coordinate(NULL),_connectivity(NULL), _isAGrid(false) {
init();
};
-MESH::MESH(const MESH &m)
+MESH::MESH(MESH &m)
{
- // VERIFIER QUE TS LES OPERATEURS DE RECOPIE DES ATTRIBUTS
- // SONT CORRECTS.
- *this = m;
+ _name=m._name;
+ _isAGrid = m._isAGrid;
+
+ if (m._coordinate != NULL)
+ _coordinate = new COORDINATE(* m._coordinate);
+ else
+ _coordinate = (COORDINATE *) NULL;
+
+ if (m._connectivity != NULL)
+ _connectivity = new CONNECTIVITY(* m._connectivity);
+ else
+ _connectivity = (CONNECTIVITY *) NULL;
+
+ _spaceDimension = m._spaceDimension;
+ _meshDimension = m._meshDimension;
+ _numberOfNodes = m._numberOfNodes;
+
+ _familyNode = m._familyNode;
+ for (int i=0; i<m._familyNode.size(); i++)
+ {
+ _familyNode[i] = new FAMILY(* m._familyNode[i]);
+ _familyNode[i]->setMesh(this);
+ }
+
+ _familyCell = m._familyCell;
+ for (int i=0; i<m._familyCell.size(); i++)
+ {
+ _familyCell[i] = new FAMILY(* m._familyCell[i]);
+ _familyCell[i]->setMesh(this);
+ }
+
+ _familyFace = m._familyFace;
+ for (int i=0; i<m._familyFace.size(); i++)
+ {
+ _familyFace[i] = new FAMILY(* m._familyFace[i]);
+ _familyFace[i]->setMesh(this);
+ }
+
+ _familyEdge = m._familyEdge;
+ for (int i=0; i<m._familyEdge.size(); i++)
+ {
+ _familyEdge[i] = new FAMILY(* m._familyEdge[i]);
+ _familyEdge[i]->setMesh(this);
+ }
+
+ _groupNode = m._groupNode;
+ for (int i=0; i<m._groupNode.size(); i++)
+ {
+ _groupNode[i] = new GROUP(* m._groupNode[i]);
+ _groupNode[i]->setMesh(this);
+ }
+
+ _groupCell = m._groupCell;
+ for (int i=0; i<m._groupCell.size(); i++)
+ {
+ _groupCell[i] = new GROUP(* m._groupCell[i]);
+ _groupCell[i]->setMesh(this);
+ }
+
+ _groupFace = m._groupFace;
+ for (int i=0; i<m._groupFace.size(); i++)
+ {
+ _groupFace[i] = new GROUP(* m._groupFace[i]);
+ _groupFace[i]->setMesh(this);
+ }
+
+ _groupEdge = m._groupEdge;
+ for (int i=0; i<m._groupEdge.size(); i++)
+ {
+ _groupEdge[i] = new GROUP(* m._groupEdge[i]);
+ _groupEdge[i]->setMesh(this);
+ }
+
+ //_drivers = m._drivers; //Recopie des drivers?
}
MESH::~MESH() {
- if ( _MEDArrayNodeFamily != (int * ) NULL) delete [] _MEDArrayNodeFamily; // SOLUTION TEMPORAIRE
- if ( _numberOfMEDCellFamily != (int * ) NULL) delete [] _numberOfMEDCellFamily;
- if ( _numberOfMEDFaceFamily != (int * ) NULL) delete [] _numberOfMEDFaceFamily;
- if ( _numberOfMEDEdgeFamily != (int * ) NULL) delete [] _numberOfMEDEdgeFamily;
- // IL FAUT FAIRE UNE BOUCLE DE DESALLOCATION
- if ( _MEDArrayCellFamily != (int **) NULL) delete [] _MEDArrayCellFamily; // SOLUTION TEMPORAIRE
- if ( _MEDArrayFaceFamily != (int **) NULL) delete [] _MEDArrayFaceFamily; // SOLUTION TEMPORAIRE
- if ( _MEDArrayEdgeFamily != (int **) NULL) delete [] _MEDArrayEdgeFamily; // SOLUTION TEMPORAIRE
-
- if (_coordinate != NULL) delete _coordinate ;
- if (_connectivity != NULL) delete _connectivity ;
+ MESSAGE("MESH::~MESH() : Destroying the Mesh");
+ if (_coordinate != ((COORDINATE *) NULL)) delete _coordinate ;
+ if (_connectivity != ((CONNECTIVITY *) NULL)) delete _connectivity ;
int size ;
size = _familyNode.size() ;
for (int i=0;i<size;i++)
for (int i=0;i<size;i++)
delete _groupEdge[i] ;
+ MESSAGE("In this object MESH there is(are) " << _drivers.size() << " driver(s)");
+
+ for (unsigned int index=0; index < _drivers.size(); index++ )
+ {
+ SCRUTE(_drivers[index]);
+ if ( _drivers[index] != NULL) delete _drivers[index];
+ }
+
}
MESH & MESH::operator=(const MESH &m)
{
+ const char * LOC = "MESH & MESH::operator=(const MESH &m) : ";
+ BEGIN_OF(LOC);
+ MESSAGE(LOC <<"Not yet implemented, operating on the object " << m);
// A FAIRE.........
// ATTENTION CET OPERATEUR DE RECOPIE EST DANGEREUX POUR LES
// reverse_nodal_connectivity = m.reverse_nodal_connectivity;
// reverse_nodal_connectivity_index = m.reverse_nodal_connectivity_index ;
// }
+ END_OF(LOC);
+
return *this;
}
-/*! Create a MESH object using a MESH driver of type (MED_DRIVER, ....) associated with file <fileName>.
- The meshname <driverName> must already exists in the file.*/
-MESH::MESH(driverTypes driverType, const string & fileName/*=""*/, const string & driverName/*=""*/) {
+/*! Create a %MESH object using a %MESH driver of type %driverTypes (MED_DRIVER, ....) associated with file fileName.
+ The meshname driverName must already exists in the file.*/
+MESH::MESH(driverTypes driverType, const string & fileName/*=""*/, const string & driverName/*=""*/) throw (MEDEXCEPTION)
+{
const char * LOC ="MESH::MESH(driverTypes driverType, const string & fileName="", const string & driverName="") : ";
-
int current;
BEGIN_OF(LOC);
init();
-
- current = addDriver(driverType,fileName,driverName);
- switch(_drivers[current]->getAccessMode() ) {
- case MED_RDONLY : {
- MESSAGE("MESH::MESH(driverTypes driverType, .....) : driverType must have a MED_RDWR accessMode");
- rmDriver(current);
- break;}
- default : {
- }
- }
+ GENDRIVER *myDriver=DRIVERFACTORY::buildDriverForMesh(driverType,fileName,this,driverName,MED_LECT);
+ current = addDriver(*myDriver);
+ delete myDriver;
_drivers[current]->open();
_drivers[current]->read();
_drivers[current]->close();
+
+// if (_isAGrid)
+// ((GRID *) this)->fillMeshAfterRead();
+
END_OF(LOC);
};
// return N;
// }
-ostream & operator<<(ostream &os, MESH &myMesh)
+ostream & MEDMEM::operator<<(ostream &os, const MESH &myMesh)
{
-// int space_dimension = myMesh.get_space_dimension();
-// os << "Space dimension "<< space_dimension << endl ;
-// int nodes_count = myMesh.get_nodes_count() ;
-// os << "Nodes count : "<< nodes_count << endl ;
-// double * coord = myMesh.get_coordinates();
-// os << "Coordinates :" << endl ;
-// string * coordinate_name = myMesh.get_coordinates_name() ;
-// string * coordinate_unit = myMesh.get_coordinates_unit() ;
-
-// for(int i=0;i<space_dimension;i++) {
-// os << " - name "<< i << " : "<< coordinate_name[i] << endl;
-// os << " - unit "<< i << " : "<< coordinate_unit[i] << endl ;
-// }
-
-// for(int j=0;j<nodes_count;j++) {
-// for(int i=0;i<space_dimension;i++)
-// os << " coord["<< i <<","<< j << "] : "<< coord[j+i*nodes_count] <<" ," ;
-// os << endl;
-// }
-
-// int cells_types_count = myMesh.get_cells_types_count() ;
-// os << "cells_types_count : " << cells_types_count << endl ;
-// for(int i=1;i<cells_types_count;i++) {
-// os << "cell type : " << myMesh.get_cells_type(i) << endl ;
-// os << " - cells count : " << myMesh.get_cells_count(i) << endl ;
-// int *conectivity = myMesh.get_nodal_connectivity(i) ;
-// int nodes_cells_count = myMesh.get_cells_type(i).get_nodes_count() ;
-// for(int j=0;j<myMesh.get_cells_count(i);j++) {
-// os << " " ;
-// for(int k=0;k<nodes_cells_count;k++)
-// os <<conectivity[j*nodes_cells_count+k] << " " ;
-// os << endl ;
-// }
-// }
-
-// int nodes_families_count = myMesh.get_nodes_families_count() ;
-// os << "nodes_families_count : " << nodes_families_count << endl ;
-// vector<FamilyNode*> nodes_Families = myMesh.get_nodes_Families() ;
-// for(int i=0;i<nodes_families_count;i++)
-// os << " - Famile "<<i<<" : "<< (FamilyNode &) (*nodes_Families[i]) << endl ;
-// os << endl ;
-
-// int cells_families_count = myMesh.get_cells_families_count() ;
-// os << "cells_families_count : " << cells_families_count << endl ;
-// vector<FamilyCell*> cells_Families = myMesh.get_cells_Families() ;
-// for(int i=0;i<cells_families_count;i++)
-// os << " - Famile "<<i<<" : "<< (*cells_Families[i]) << endl ;
-// os << endl ;
-
- return os ;
+ int spacedimension = myMesh.getSpaceDimension();
+ int meshdimension = myMesh.getMeshDimension();
+ int numberofnodes = myMesh.getNumberOfNodes();
+
+ os << "Space Dimension : " << spacedimension << endl << endl;
+
+ os << "Mesh Dimension : " << meshdimension << endl << endl;
+
+ const double * coordinates = myMesh.getCoordinates(MED_FULL_INTERLACE);
+ os << "SHOW NODES COORDINATES : " << endl;
+
+ os << "Name :" << endl;
+ const string * coordinatesnames = myMesh.getCoordinatesNames();
+ for(int i=0; i<spacedimension ; i++)
+ {
+ os << " - " << coordinatesnames[i] << endl;
+ }
+ os << "Unit :" << endl;
+ const string * coordinatesunits = myMesh.getCoordinatesUnits();
+ for(int i=0; i<spacedimension ; i++)
+ {
+ os << " - " << coordinatesunits[i] << endl;
+ }
+ for(int i=0; i<numberofnodes ; i++)
+ {
+ os << "Nodes " << i+1 << " : ";
+ for (int j=0; j<spacedimension ; j++)
+ os << coordinates[i*spacedimension+j] << " ";
+ os << endl;
+ }
+
+ os << endl << "SHOW CONNECTIVITY :" << endl;
+ os << *myMesh._connectivity << endl;
+
+ medEntityMesh entity;
+ os << endl << "SHOW FAMILIES :" << endl << endl;
+ for (int k=1; k<=4; k++)
+ {
+ if (k==1) entity = MED_NODE;
+ if (k==2) entity = MED_CELL;
+ if (k==3) entity = MED_FACE;
+ if (k==4) entity = MED_EDGE;
+ int numberoffamilies = myMesh.getNumberOfFamilies(entity);
+ os << "NumberOfFamilies on "<< entNames[entity] <<" : "<<numberoffamilies<<endl;
+ using namespace MED_EN;
+ for (int i=1; i<numberoffamilies+1;i++)
+ {
+ os << * myMesh.getFamily(entity,i) << endl;
+ }
+ }
+
+ os << endl << "SHOW GROUPS :" << endl << endl;
+ for (int k=1; k<=4; k++)
+ {
+ if (k==1) entity = MED_NODE;
+ if (k==2) entity = MED_CELL;
+ if (k==3) entity = MED_FACE;
+ if (k==4) entity = MED_EDGE;
+ int numberofgroups = myMesh.getNumberOfGroups(entity);
+ os << "NumberOfGroups on "<< entNames[entity] <<" : "<<numberofgroups<<endl;
+ using namespace MED_EN;
+ for (int i=1; i<numberofgroups+1;i++)
+ {
+ os << * myMesh.getGroup(entity,i) << endl;
+ }
+ }
+
+ return os;
}
/*!
Return -1 if not found.
*/
-int MESH::getElementNumber(medConnectivity ConnectivityType, medEntityMesh Entity, medGeometryElement Type, int * connectivity)
+int MESH::getElementNumber(medConnectivity ConnectivityType, medEntityMesh Entity, medGeometryElement Type, int * connectivity) const
{
const char* LOC="MESH::getElementNumber " ;
BEGIN_OF(LOC) ;
else
numberOfValue = myType.getNumberOfNodes() ; // nodes
- int * myReverseConnectivityValue = getReverseConnectivity(ConnectivityType,Entity) ;
- int * myReverseConnectivityIndex = getReverseConnectivityIndex(ConnectivityType,Entity) ;
+ const int * myReverseConnectivityValue = getReverseConnectivity(ConnectivityType,Entity) ;
+ const int * myReverseConnectivityIndex = getReverseConnectivityIndex(ConnectivityType,Entity) ;
// First node or face/edge
int indexBegin = myReverseConnectivityIndex[connectivity[0]-1] ;
For instance, we could get only face in 3D and edge in 2D.
*/
-SUPPORT * MESH::getBoundaryElements(medEntityMesh Entity) throw (MEDEXCEPTION)
+SUPPORT * MESH::getBoundaryElements(medEntityMesh Entity)
+ throw (MEDEXCEPTION)
{
const char * LOC = "MESH::getBoundaryElements : " ;
BEGIN_OF(LOC) ;
mySupport->setAll(false);
- int * myConnectivityValue = getReverseConnectivity(MED_DESCENDING) ;
- int * myConnectivityIndex = getReverseConnectivityIndex(MED_DESCENDING) ;
+ const int * myConnectivityValue = getReverseConnectivity(MED_DESCENDING) ;
+ const int * myConnectivityIndex = getReverseConnectivityIndex(MED_DESCENDING) ;
int numberOf = getNumberOfElements(Entity,MED_ALL_ELEMENTS) ;
list<int> myElementsList ;
int size = 0 ;
- SCRUTE(numberOf) ;
for (int i=0 ; i<numberOf; i++)
if (myConnectivityValue[myConnectivityIndex[i]] == 0) {
- SCRUTE(i+1) ;
myElementsList.push_back(i+1) ;
size++ ;
}
- SCRUTE(size) ;
// Well, we must know how many geometric type we have found
int * myListArray = new int[size] ;
int id = 0 ;
list<int>::iterator myElementsListIt ;
for (myElementsListIt=myElementsList.begin();myElementsListIt!=myElementsList.end();myElementsListIt++) {
myListArray[id]=(*myElementsListIt) ;
- SCRUTE(id);
- SCRUTE(myListArray[id]);
id ++ ;
}
medGeometryElement* geometricType ;
int * numberOfGaussPoint ;
int * geometricTypeNumber ;
- int * numberOfEntities ;
- MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY() ;
+ int * numberOfElements ;
+ //MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY() ;
int * mySkyLineArrayIndex ;
int numberOfType = getNumberOfTypes(Entity) ;
if (numberOfType == 1) { // wonderfull : it's easy !
numberOfGeometricType = 1 ;
geometricType = new medGeometryElement[1] ;
- medGeometryElement * allType = getTypes(Entity);
+ const medGeometryElement * allType = getTypes(Entity);
geometricType[0] = allType[0] ;
numberOfGaussPoint = new int[1] ;
numberOfGaussPoint[0] = 1 ;
geometricTypeNumber = new int[1] ; // not use, but initialized to nothing
geometricTypeNumber[0] = 0 ;
- numberOfEntities = new int[1] ;
- numberOfEntities[0] = size ;
+ numberOfElements = new int[1] ;
+ numberOfElements[0] = size ;
mySkyLineArrayIndex = new int[2] ;
mySkyLineArrayIndex[0]=1 ;
mySkyLineArrayIndex[1]=1+size ;
}
numberOfGeometricType = theType.size() ;
geometricType = new medGeometryElement[numberOfGeometricType] ;
- medGeometryElement * allType = getTypes(Entity);
+ //const medGeometryElement * allType = getTypes(Entity); !! UNUSZED VARIABLE !!
numberOfGaussPoint = new int[numberOfGeometricType] ;
geometricTypeNumber = new int[numberOfGeometricType] ; // not use, but initialized to nothing
- numberOfEntities = new int[numberOfGeometricType] ;
+ numberOfElements = new int[numberOfGeometricType] ;
mySkyLineArrayIndex = new int[numberOfGeometricType+1] ;
int index = 0 ;
mySkyLineArrayIndex[0]=1 ;
geometricType[index] = (*theTypeIt).first ;
numberOfGaussPoint[index] = 1 ;
geometricTypeNumber[index] = 0 ;
- numberOfEntities[index] = (*theTypeIt).second ;
- mySkyLineArrayIndex[index+1]=mySkyLineArrayIndex[index]+numberOfEntities[index] ;
+ numberOfElements[index] = (*theTypeIt).second ;
+ mySkyLineArrayIndex[index+1]=mySkyLineArrayIndex[index]+numberOfElements[index] ;
index++ ;
}
}
- mySkyLineArray->setMEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
+ //mySkyLineArray->setMEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
+ MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
mySupport->setNumberOfGeometricType(numberOfGeometricType) ;
mySupport->setGeometricType(geometricType) ;
mySupport->setNumberOfGaussPoint(numberOfGaussPoint) ;
+ mySupport->setNumberOfElements(numberOfElements) ;
+ mySupport->setTotalNumberOfElements(size) ;
// mySupport->setGeometricTypeNumber(geometricTypeNumber) ;
- mySupport->setNumberOfEntities(numberOfEntities) ;
- mySupport->setTotalNumberOfEntities(size) ;
mySupport->setNumber(mySkyLineArray) ;
+ delete[] numberOfElements;
+ delete[] geometricTypeNumber;
+ delete[] numberOfGaussPoint;
+ delete[] geometricType;
+ delete[] mySkyLineArrayIndex;
+ delete[] myListArray;
+// delete mySkyLineArray;
+
END_OF(LOC) ;
return mySupport ;
}
-FIELD<double>* MESH::getVolume(const SUPPORT * Support) throw (MEDEXCEPTION)
+FIELD<double>* MESH::getVolume(const SUPPORT * Support) const throw (MEDEXCEPTION)
{
+ const char * LOC = "MESH::getVolume(SUPPORT*) : ";
+ BEGIN_OF(LOC);
+
// Support must be on 3D elements
- MESSAGE("MESH::getVolume(SUPPORT*)");
// Make sure that the MESH class is the same as the MESH class attribut
// in the class Support
MESH* myMesh = Support->getMesh();
if (this != myMesh)
- throw MEDEXCEPTION("MESH::getVolume(SUPPORT*) no compatibility between *this and SUPPORT::_mesh !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
int dim_space = getSpaceDimension();
medEntityMesh support_entity = Support->getEntity();
bool onAll = Support->isOnAllElements();
int nb_type, length_values;
- medGeometryElement* types;
+ const medGeometryElement* types;
int nb_entity_type;
// !!!! WARNING : use of nodal global numbering in the mesh !!!!
- int* global_connectivity;
+ const int* global_connectivity;
- if (onAll)
- {
- nb_type = myMesh->getNumberOfTypes(support_entity);
- length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
- types = getTypes(support_entity);
- }
- else
+// if (onAll)
+// {
+// nb_type = myMesh->getNumberOfTypes(support_entity);
+// length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
+// types = getTypes(support_entity);
+// }
+// else
{
nb_type = Support->getNumberOfTypes();
length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
}
int index;
- FIELD<double>* Volume;
-
- Volume = new FIELD<double>::FIELD(Support,1);
+ FIELD<double>* Volume = new FIELD<double>(Support,1);
// double *volume = new double [length_values];
Volume->setName("VOLUME");
Volume->setDescription("cells volume");
Volume->setOrderNumber(0);
Volume->setTime(0.0);
- double *volume = Volume->getValue(MED_FULL_INTERLACE);
+ //const double *volume = Volume->getValue(MED_FULL_INTERLACE);
+ MEDARRAY<double> *volume = Volume->getvalue();
- index = 0;
+ index = 1;
const double * coord = getCoordinates(MED_FULL_INTERLACE);
for (int i=0;i<nb_type;i++)
}
else
{
- nb_entity_type = Support->getNumberOfElements(type);
-
- int * supp_number = Support->getNumber(type);
- int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
- int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
- global_connectivity = new int[(type%100)*nb_entity_type];
-
- for (int k_type = 0; k_type<nb_entity_type; k_type++) {
- for (int j_ent = 0; j_ent<(type%100); j_ent++) {
- global_connectivity[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
- }
- }
+ // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all"));
+
+ nb_entity_type = Support->getNumberOfElements(type);
+
+ const int * supp_number = Support->getNumber(type);
+ const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
+ const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
+ int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
+
+ for (int k_type = 0; k_type<nb_entity_type; k_type++) {
+ for (int j_ent = 0; j_ent<(type%100); j_ent++) {
+ global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
+ }
+ }
+ global_connectivity = global_connectivity_tmp ;
}
switch (type)
(x2-x1)*((y3-y1)*(z4-z1) - (z3-z1)*(y4-y1)) +
(x4-x1)*((y3-y1)*(z2-z1) - (z3-z1)*(y2-y1)))/6.0;
- volume[index] = xvolume ;
+ //volume[index] = xvolume ;
+ volume->setIJ(index,1,xvolume) ;
index++;
}
break;
(x5-x1)*((y4-y1)*(z3-z1) - (z4-z1)*(y3-y1)))
)/6.0;
- volume[index] = xvolume ;
+ //volume[index] = xvolume ;
+ volume->setIJ(index,1,xvolume) ;
index = index++;
}
break;
xvolume = -2.0*(2.0*(A + B + D + E + G + H) + C + F + P)/9.0;
- volume[index] = xvolume ;
+ //volume[index] = xvolume ;
+ volume->setIJ(index,1,xvolume) ;
index++;
}
break;
V + W) + 2.0*(I + R + U + X + Y + Z) +
AA)/27.0;
- volume[index] = xvolume ;
+ //volume[index] = xvolume ;
+ volume->setIJ(index,1,xvolume) ;
index++;
}
break;
}
default :
- throw MEDEXCEPTION("MESH::getVolume(SUPPORT*) Bad Support to get Volumes on it !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Volumes on it !"));
break;
}
+
+ if (!onAll) delete [] global_connectivity ;
}
return Volume;
}
-FIELD<double>* MESH::getArea(const SUPPORT * Support) throw (MEDEXCEPTION)
+FIELD<double>* MESH::getArea(const SUPPORT * Support) const throw (MEDEXCEPTION)
{
+ const char * LOC = "MESH::getArea(SUPPORT*) : ";
+ BEGIN_OF(LOC);
+
// Support must be on 2D elements
- MESSAGE("MESH::getArea(SUPPORT*)");
// Make sure that the MESH class is the same as the MESH class attribut
// in the class Support
MESH* myMesh = Support->getMesh();
if (this != myMesh)
- throw MEDEXCEPTION("MESH::getArea(SUPPORT*) no compatibility between *this and SUPPORT::_mesh !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
int dim_space = getSpaceDimension();
medEntityMesh support_entity = Support->getEntity();
bool onAll = Support->isOnAllElements();
int nb_type, length_values;
- medGeometryElement* types;
+ const medGeometryElement* types;
int nb_entity_type;
// !!!! WARNING : use of nodal global numbering in the mesh !!!!
- int* global_connectivity;
+ const int* global_connectivity;
- if (onAll)
- {
- nb_type = myMesh->getNumberOfTypes(support_entity);
- length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
- types = getTypes(support_entity);
- }
- else
+// if (onAll)
+// {
+// nb_type = myMesh->getNumberOfTypes(support_entity);
+// length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
+// types = getTypes(support_entity);
+// }
+// else
{
nb_type = Support->getNumberOfTypes();
length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
int index;
FIELD<double>* Area;
- Area = new FIELD<double>::FIELD(Support,1);
+ Area = new FIELD<double>(Support,1);
Area->setName("AREA");
Area->setDescription("cells or faces area");
}
else
{
- nb_entity_type = Support->getNumberOfElements(type);
+ // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all !"));
+
+ nb_entity_type = Support->getNumberOfElements(type);
+
+ const int * supp_number = Support->getNumber(type);
+ const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
+ const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
+ int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
+
+ for (int k_type = 0; k_type<nb_entity_type; k_type++) {
+ for (int j_ent = 0; j_ent<(type%100); j_ent++) {
+ global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
+ }
+ }
+
+ global_connectivity = global_connectivity_tmp ;
- int * supp_number = Support->getNumber(type);
- int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
- int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
- global_connectivity = new int[(type%100)*nb_entity_type];
-
- for (int k_type = 0; k_type<nb_entity_type; k_type++) {
- for (int j_ent = 0; j_ent<(type%100); j_ent++) {
- global_connectivity[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
- }
- }
}
switch (type)
break;
}
default :
- throw MEDEXCEPTION("MESH::getArea(SUPPORT*) Bad Support to get Areas on it !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Areas on it !"));
break;
}
+
+ if (!onAll) delete [] global_connectivity ;
}
Area->setValue(MED_FULL_INTERLACE,area);
-
+ delete[] area;
return Area;
}
-FIELD<double>* MESH::getLength(const SUPPORT * Support) throw (MEDEXCEPTION)
+FIELD<double>* MESH::getLength(const SUPPORT * Support) const throw (MEDEXCEPTION)
{
+ const char * LOC = "MESH::getLength(SUPPORT*) : ";
+ BEGIN_OF(LOC);
+
// Support must be on 1D elements
- MESSAGE("MESH::getLength(SUPPORT*)");
// Make sure that the MESH class is the same as the MESH class attribut
// in the class Support
MESH* myMesh = Support->getMesh();
if (this != myMesh)
- throw MEDEXCEPTION("MESH::getLength(SUPPORT*) no compatibility between *this and SUPPORT::_mesh !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
int dim_space = getSpaceDimension();
medEntityMesh support_entity = Support->getEntity();
bool onAll = Support->isOnAllElements();
int nb_type, length_values;
- medGeometryElement* types;
+ const medGeometryElement* types;
int nb_entity_type;
// !!!! WARNING : use of nodal global numbering in the mesh !!!!
- int* global_connectivity;
+ const int* global_connectivity;
- if (onAll)
- {
- nb_type = myMesh->getNumberOfTypes(support_entity);
- length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
- types = getTypes(support_entity);
- }
- else
+// if (onAll)
+// {
+// nb_type = myMesh->getNumberOfTypes(support_entity);
+// length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
+// types = getTypes(support_entity);
+// }
+// else
{
nb_type = Support->getNumberOfTypes();
length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
int index;
FIELD<double>* Length;
- Length = new FIELD<double>::FIELD(Support,1);
+ Length = new FIELD<double>(Support,1);
// double *length = new double [length_values];
Length->setValueType(MED_REEL64);
- double *length = Length->getValue(MED_FULL_INTERLACE);
+ //const double *length = Length->getValue(MED_FULL_INTERLACE);
+ MEDARRAY<double> * length = Length->getvalue();
const double * coord = getCoordinates(MED_FULL_INTERLACE);
- index = 0;
+ index = 1;
for (int i=0;i<nb_type;i++)
{
}
else
{
- nb_entity_type = Support->getNumberOfElements(type);
+ //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all !"));
- int * supp_number = Support->getNumber(type);
- int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
- int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
- global_connectivity = new int[(type%100)*nb_entity_type];
+ nb_entity_type = Support->getNumberOfElements(type);
- for (int k_type = 0; k_type<nb_entity_type; k_type++) {
- for (int j_ent = 0; j_ent<(type%100); j_ent++) {
- global_connectivity[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
- }
- }
+ const int * supp_number = Support->getNumber(type);
+ const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
+ const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
+ int* global_connectivity_tmp = new int[(type%100)*nb_entity_type];
+
+ for (int k_type = 0; k_type<nb_entity_type; k_type++) {
+ for (int j_ent = 0; j_ent<(type%100); j_ent++) {
+ global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
+ }
+ }
+ global_connectivity = global_connectivity_tmp ;
}
switch (type)
xlength = sqrt((x1 - x2)*(x1 - x2) + (y1 - y2)*(y1 - y2) +
(z1 - z2)*(z1 - z2));
- length[index] = xlength ;
+ length->setIJ(index,1,xlength) ;
index++;
}
break;
}
default :
- throw MEDEXCEPTION("MESH::getLength(SUPPORT*) Bad Support to get Lengths on it !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Lengths on it !"));
break;
}
+
+ if (!onAll) delete [] global_connectivity ;
}
return Length;
}
-FIELD<double>* MESH::getNormal(const SUPPORT * Support) throw (MEDEXCEPTION)
+FIELD<double>* MESH::getNormal(const SUPPORT * Support) const throw (MEDEXCEPTION)
{
+ const char * LOC = "MESH::getNormal(SUPPORT*) : ";
+ BEGIN_OF(LOC);
+
// Support must be on 2D or 1D elements
- MESSAGE("MESH::getNormal(SUPPORT*)");
// Make sure that the MESH class is the same as the MESH class attribut
// in the class Support
MESH* myMesh = Support->getMesh();
if (this != myMesh)
- throw MEDEXCEPTION("MESH::getNormal(SUPPORT*) no compatibility between *this and SUPPORT::_mesh : pointeur problem !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh : pointeur problem !"));
int dim_space = getSpaceDimension();
medEntityMesh support_entity = Support->getEntity();
bool onAll = Support->isOnAllElements();
int nb_type, length_values;
- medGeometryElement* types;
+ const medGeometryElement* types;
int nb_entity_type;
// !!!! WARNING : use of nodal global numbering in the mesh !!!!
- int* global_connectivity;
+ const int* global_connectivity;
- if (onAll)
- {
- nb_type = myMesh->getNumberOfTypes(support_entity);
- length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
- types = getTypes(support_entity);
- }
- else
+// if (onAll)
+// {
+// nb_type = myMesh->getNumberOfTypes(support_entity);
+// length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
+// types = getTypes(support_entity);
+// }
+// else
{
nb_type = Support->getNumberOfTypes();
length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
int index;
- FIELD<double>* Normal = new FIELD<double>::FIELD(Support,dim_space);
+ FIELD<double>* Normal = new FIELD<double>(Support,dim_space);
Normal->setName("NORMAL");
Normal->setDescription("cells or faces normal");
- for (int k=0;k<dim_space;k++) {
- int kp1 = k + 1;
- Normal->setComponentName(kp1,"normal");
- Normal->setComponentDescription(kp1,"desc-comp");
- Normal->setMEDComponentUnit(kp1,"unit");
+ for (int k=1;k<=dim_space;k++) {
+ Normal->setComponentName(k,"normal");
+ Normal->setComponentDescription(k,"desc-comp");
+ Normal->setMEDComponentUnit(k,"unit");
}
Normal->setValueType(MED_REEL64);
- Normal->setIterationNumber(0);
- Normal->setOrderNumber(0);
+ Normal->setIterationNumber(MED_NOPDT);
+ Normal->setOrderNumber(MED_NONOR);
Normal->setTime(0.0);
double * normal = new double [dim_space*length_values];
(type==MED_QUAD4) || (type==MED_QUAD8)) &&
(dim_space != 3)) || (((type==MED_SEG2) || (type==MED_SEG3)) &&
(dim_space != 2)) )
- throw MEDEXCEPTION("MESH::getNormal(SUPPORT*) no compatibility between *this and SUPPORT::_mesh : dimension problem !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh : dimension problem !"));
double xnormal1, xnormal2, xnormal3 ;
}
else
{
- nb_entity_type = Support->getNumberOfElements(type);
+ // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all for instance !"));
+ nb_entity_type = Support->getNumberOfElements(type);
+
+ const int * supp_number = Support->getNumber(type);
+ const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
+ const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
+ int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
+
+ for (int k_type = 0; k_type<nb_entity_type; k_type++) {
+ for (int j_ent = 0; j_ent<(type%100); j_ent++) {
+ global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
+ }
+ }
+
+ global_connectivity = global_connectivity_tmp ;
- int * supp_number = Support->getNumber(type);
- int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
- int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
- global_connectivity = new int[(type%100)*nb_entity_type];
-
- for (int k_type = 0; k_type<nb_entity_type; k_type++) {
- for (int j_ent = 0; j_ent<(type%100); j_ent++) {
- global_connectivity[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
- }
- }
}
switch (type)
int N2 = global_connectivity[tria_index+1]-1;
int N3 = global_connectivity[tria_index+2]-1;
- double xarea;
+ //double xarea; !! UNUSED VARIABLE !!
double x1 = coord[dim_space*N1];
double x2 = coord[dim_space*N2];
double x3 = coord[dim_space*N3];
((x4-x1)*(z2-z1) - (x2-x1)*(z4-z1)) +
((x2-x1)*(y4-y1) - (x4-x1)*(y2-y1))*
((x2-x1)*(y4-y1) - (x4-x1)*(y2-y1))) +
- sqrt(((y3-y2)*(z4-z3) - (y4-y3)*(z3-z2))*
- ((y3-y2)*(z4-z3) - (y4-y3)*(z3-z2)) +
- ((x4-x3)*(z2-z2) - (x3-x2)*(z4-z3))*
- ((x4-x3)*(z4-z2) - (x3-x2)*(z4-z3)) +
- ((x3-x2)*(y4-y3) - (x4-x3)*(y3-y2))*
- ((x3-x2)*(y4-y3) - (x4-x3)*(y3-y2))))/2.0;
+ sqrt(((y4-y3)*(z2-z3) - (y2-y3)*(z4-z3))*
+ ((y4-y3)*(z2-z3) - (y2-y3)*(z4-z3)) +
+ ((x2-x3)*(z4-z3) - (x4-x3)*(z2-z3))*
+ ((x2-x3)*(z4-z3) - (x4-x3)*(z2-z3)) +
+ ((x4-x3)*(y2-y3) - (x2-x3)*(y4-y3))*
+ ((x4-x3)*(y2-y3) - (x2-x3)*(y4-y3))))/2.0;
+
+// sqrt(((y3-y2)*(z4-z3) - (y4-y3)*(z3-z2))*
+// ((y3-y2)*(z4-z3) - (y4-y3)*(z3-z2)) +
+// ((x4-x3)*(z2-z2) - (x3-x2)*(z4-z3))*
+// ((x4-x3)*(z4-z2) - (x3-x2)*(z4-z3)) +
+// ((x3-x2)*(y4-y3) - (x4-x3)*(y3-y2))*
+// ((x3-x2)*(y4-y3) - (x4-x3)*(y3-y2))))/2.0;
xnormal1 = xnormal1*xarea;
xnormal2 = xnormal2*xarea;
break;
}
default :
- throw MEDEXCEPTION("MESH::getNormal(SUPPORT*) Bad Support to get Normals on it !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Normals on it !"));
break;
}
+
+ if (!onAll) delete [] global_connectivity ;
}
Normal->setValue(MED_FULL_INTERLACE,normal);
+ delete[] normal ;
+
+ END_OF(LOC);
return Normal;
}
-FIELD<double>* MESH::getBarycenter(const SUPPORT * Support) throw (MEDEXCEPTION)
+FIELD<double>* MESH::getBarycenter(const SUPPORT * Support) const throw (MEDEXCEPTION)
{
- MESSAGE("MESH::getBarycenter(SUPPORT*)");
+ const char * LOC = "MESH::getBarycenter(SUPPORT*) : ";
+ BEGIN_OF(LOC);
// Make sure that the MESH class is the same as the MESH class attribut
// in the class Support
MESH* myMesh = Support->getMesh();
if (this != myMesh)
- throw MEDEXCEPTION("MESH::getBarycenter(SUPPORT*) no compatibility between *this and SUPPORT::_mesh !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
int dim_space = getSpaceDimension();
medEntityMesh support_entity = Support->getEntity();
bool onAll = Support->isOnAllElements();
int nb_type, length_values;
- medGeometryElement* types;
+ const medGeometryElement* types;
int nb_entity_type;
// !!!! WARNING : use of nodal global numbering in the mesh !!!!
- int* global_connectivity;
+ const int* global_connectivity;
- if (onAll)
- {
- nb_type = myMesh->getNumberOfTypes(support_entity);
- length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
- types = getTypes(support_entity);
- }
- else
+// if (onAll)
+// {
+// nb_type = myMesh->getNumberOfTypes(support_entity);
+// length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
+// types = getTypes(support_entity);
+// }
+// else
{
nb_type = Support->getNumberOfTypes();
length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
int index;
FIELD<double>* Barycenter;
- Barycenter = new FIELD<double>::FIELD(Support,dim_space);
+ Barycenter = new FIELD<double>(Support,dim_space);
Barycenter->setName("BARYCENTER");
Barycenter->setDescription("cells or faces barycenter");
}
else
{
- nb_entity_type = Support->getNumberOfElements(type);
-
- int * supp_number = Support->getNumber(type);
- int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
- int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
- global_connectivity = new int[(type%100)*nb_entity_type];
-
- for (int k_type = 0; k_type<nb_entity_type; k_type++) {
- for (int j_ent = 0; j_ent<(type%100); j_ent++) {
- global_connectivity[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
- }
- }
+ // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all !"));
+ nb_entity_type = Support->getNumberOfElements(type);
+
+ const int * supp_number = Support->getNumber(type);
+ const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
+ const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
+ int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
+
+ for (int k_type = 0; k_type<nb_entity_type; k_type++) {
+ for (int j_ent = 0; j_ent<(type%100); j_ent++) {
+ global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
+ }
+ }
+ global_connectivity = global_connectivity_tmp;
}
switch (type)
break;
}
default :
- throw MEDEXCEPTION("MESH::getBarycenter(SUPPORT*) Bad Support to get a barycenter on it (in fact unknown type) !");
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get a barycenter on it (in fact unknown type) !"));
break;
}
+
+ if (!onAll) delete [] global_connectivity ;
}
Barycenter->setValue(MED_FULL_INTERLACE,barycenter);
+ delete[] barycenter ;
+
+ END_OF(LOC);
+
return Barycenter;
}
+
+//=======================================================================
+//function : checkGridFillCoords
+//purpose : if this->_isAGrid, assure that _coordinate is filled
+//=======================================================================
+
+// inline void MESH::checkGridFillCoords() const
+// {
+// if (_isAGrid)
+// ((GRID *) this)->fillCoordinates();
+// }
+
+//=======================================================================
+//function : checkGridFillConnectivity
+//purpose : if this->_isAGrid, assure that _connectivity is filled
+//=======================================================================
+
+// inline void MESH::checkGridFillConnectivity() const
+// {
+// if (_isAGrid)
+// ((GRID *) this)->fillConnectivity();
+// }
+
+bool MESH::isEmpty() const
+{
+ bool notempty = _name != "NOT DEFINED" || _coordinate != NULL || _connectivity != NULL ||
+ _spaceDimension !=MED_INVALID || _meshDimension !=MED_INVALID ||
+ _numberOfNodes !=MED_INVALID || _groupNode.size() != 0 ||
+ _familyNode.size() != 0 || _groupCell.size() != 0 ||
+ _familyCell.size() != 0 || _groupFace.size() != 0 ||
+ _familyFace.size() != 0 || _groupEdge.size() != 0 ||
+ _familyEdge.size() != 0 || _isAGrid != 0 ;
+ return !notempty;
+}
+
+void MESH::read(int index)
+{
+ const char * LOC = "MESH::read(int index=0) : ";
+ BEGIN_OF(LOC);
+
+ if (_drivers[index]) {
+ _drivers[index]->open();
+ _drivers[index]->read();
+ _drivers[index]->close();
+ }
+ else
+ throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
+ << "The index given is invalid, index must be between 0 and |"
+ << _drivers.size()
+ )
+ );
+// if (_isAGrid)
+// ((GRID *) this)->fillMeshAfterRead();
+
+ END_OF(LOC);
+}
+//=======================================================================
+//function : getSkin
+//purpose :
+//=======================================================================
+
+SUPPORT * MESH::getSkin(const SUPPORT * Support3D) throw (MEDEXCEPTION)
+{
+ const char * LOC = "MESH::getSkin : " ;
+ BEGIN_OF(LOC) ;
+ // some test :
+ if (this != Support3D->getMesh())
+ throw MEDEXCEPTION(STRING(LOC) << "no compatibility between *this and SUPPORT::_mesh !");
+ if (_meshDimension != 3 || Support3D->getEntity() != MED_CELL)
+ throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Defined on 3D cells only"));
+
+ // well, all rigth !
+ SUPPORT * mySupport = new SUPPORT(this,"Skin",MED_FACE);
+ mySupport->setAll(false);
+
+ list<int> myElementsList ;
+ int i,j, size = 0 ;
+
+ calculateConnectivity(MED_FULL_INTERLACE, MED_DESCENDING, MED_CELL);
+ if (Support3D->isOnAllElements())
+ {
+ int * myConnectivityValue = const_cast <int*> (getReverseConnectivity(MED_DESCENDING)) ;
+ int nbFaces = getNumberOfElements(MED_FACE,MED_ALL_ELEMENTS);
+ for (i=0, j=1 ; j<=nbFaces; ++j, i += 2)
+ {
+ int cellNb1 = myConnectivityValue [i];
+ int cellNb2 = myConnectivityValue [i+1];
+ //MESSAGE( " FACE # " << j << " -- Cells: " << cellNb1 << ", " << cellNb2 );
+ if ((cellNb1 == 0 || cellNb2 == 0) && (cellNb1 + cellNb2 > 0))
+ {
+ myElementsList.push_back( j ) ;
+ size++ ;
+ }
+ }
+ }
+ else
+ {
+ map<int,int> FaceNbEncounterNb;
+
+ int * myConnectivityValue = const_cast <int *>
+ (getConnectivity(MED_FULL_INTERLACE, MED_DESCENDING,
+ MED_CELL, MED_ALL_ELEMENTS));
+ int * myConnectivityIndex = const_cast <int *> (getConnectivityIndex(MED_DESCENDING, MED_CELL));
+ int * myCellNbs = const_cast <int *> (Support3D->getnumber()->getValue());
+ int nbCells = Support3D->getnumber()->getLength();
+ for (i=0; i<nbCells; ++i)
+ {
+ int cellNb = myCellNbs [ i ];
+ int faceFirst = myConnectivityIndex[ cellNb-1 ];
+ int faceLast = myConnectivityIndex[ cellNb ];
+ for (j = faceFirst; j < faceLast; ++j)
+ {
+ int faceNb = abs( myConnectivityValue [ j-1 ] );
+ //MESSAGE( "Cell # " << i << " -- Face: " << faceNb);
+ if (FaceNbEncounterNb.find( faceNb ) == FaceNbEncounterNb.end())
+ FaceNbEncounterNb[ faceNb ] = 1;
+ else
+ FaceNbEncounterNb[ faceNb ] ++;
+ }
+ }
+ map<int,int>::iterator FaceNbEncounterNbItr;
+ for (FaceNbEncounterNbItr = FaceNbEncounterNb.begin();
+ FaceNbEncounterNbItr != FaceNbEncounterNb.end();
+ FaceNbEncounterNbItr ++)
+ if ((*FaceNbEncounterNbItr).second == 1)
+ {
+ myElementsList.push_back( (*FaceNbEncounterNbItr).first) ;
+ size++ ;
+ }
+ }
+ // Well, we must know how many geometric type we have found
+ int * myListArray = new int[size] ;
+ int id = 0 ;
+ list<int>::iterator myElementsListIt ;
+ for (myElementsListIt=myElementsList.begin();myElementsListIt!=myElementsList.end();myElementsListIt++) {
+ myListArray[id]=(*myElementsListIt) ;
+ id ++ ;
+ }
+
+ int numberOfGeometricType ;
+ medGeometryElement* geometricType ;
+ int * numberOfGaussPoint ;
+ int * geometricTypeNumber ;
+ int * numberOfEntities ;
+ // MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY() ;
+ int * mySkyLineArrayIndex ;
+
+ int numberOfType = getNumberOfTypes(MED_FACE) ;
+ if (numberOfType == 1) { // wonderfull : it's easy !
+ numberOfGeometricType = 1 ;
+ geometricType = new medGeometryElement[1] ;
+ const medGeometryElement * allType = getTypes(MED_FACE);
+ geometricType[0] = allType[0] ;
+ numberOfGaussPoint = new int[1] ;
+ numberOfGaussPoint[0] = 1 ;
+ geometricTypeNumber = new int[1] ; // not use, but initialized to nothing
+ geometricTypeNumber[0] = 0 ;
+ numberOfEntities = new int[1] ;
+ numberOfEntities[0] = size ;
+ mySkyLineArrayIndex = new int[2] ;
+ mySkyLineArrayIndex[0]=1 ;
+ mySkyLineArrayIndex[1]=1+size ;
+ }
+ else {// hemmm
+ map<medGeometryElement,int> theType ;
+ for (myElementsListIt=myElementsList.begin();myElementsListIt!=myElementsList.end();myElementsListIt++) {
+ medGeometryElement myType = getElementType(MED_FACE,*myElementsListIt) ;
+ if (theType.find(myType) != theType.end() )
+ theType[myType]+=1 ;
+ else
+ theType[myType]=1 ;
+ }
+ numberOfGeometricType = theType.size() ;
+ geometricType = new medGeometryElement[numberOfGeometricType] ;
+ //const medGeometryElement * allType = getTypes(MED_FACE); !! UNUSED VARIABLE !!
+ numberOfGaussPoint = new int[numberOfGeometricType] ;
+ geometricTypeNumber = new int[numberOfGeometricType] ; // not use, but initialized to nothing
+ numberOfEntities = new int[numberOfGeometricType] ;
+ mySkyLineArrayIndex = new int[numberOfGeometricType+1] ;
+ int index = 0 ;
+ mySkyLineArrayIndex[0]=1 ;
+ map<medGeometryElement,int>::iterator theTypeIt ;
+ for (theTypeIt=theType.begin();theTypeIt!=theType.end();theTypeIt++) {
+ geometricType[index] = (*theTypeIt).first ;
+ numberOfGaussPoint[index] = 1 ;
+ geometricTypeNumber[index] = 0 ;
+ numberOfEntities[index] = (*theTypeIt).second ;
+ mySkyLineArrayIndex[index+1]=mySkyLineArrayIndex[index]+numberOfEntities[index] ;
+ index++ ;
+ }
+ }
+ // mySkyLineArray->setMEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
+ MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
+
+ mySupport->setNumberOfGeometricType(numberOfGeometricType) ;
+ mySupport->setGeometricType(geometricType) ;
+ mySupport->setNumberOfGaussPoint(numberOfGaussPoint) ;
+ // mySupport->setGeometricTypeNumber(geometricTypeNumber) ;
+ mySupport->setNumberOfElements(numberOfEntities) ;
+ mySupport->setTotalNumberOfElements(size) ;
+ mySupport->setNumber(mySkyLineArray) ;
+
+ delete[] numberOfEntities;
+ delete[] geometricTypeNumber;
+ delete[] numberOfGaussPoint;
+ delete[] geometricType;
+ delete[] mySkyLineArrayIndex;
+ delete[] myListArray;
+// delete mySkyLineArray;
+
+ END_OF(LOC) ;
+ return mySupport ;
+
+}
+
+/*!
+ return a SUPPORT pointer on the union of all SUPPORTs in Supports.
+ You should delete this pointer after use to avoid memory leaks.
+*/
+SUPPORT * MESH::mergeSupports(const vector<SUPPORT *> Supports) const throw (MEDEXCEPTION)
+{
+ const char * LOC = "MESH:::mergeSupports(const vector<SUPPORT *> ) : " ;
+ BEGIN_OF(LOC) ;
+
+ SUPPORT * returnedSupport;
+ string returnedSupportName;
+ string returnedSupportDescription;
+ char * returnedSupportNameChar;
+ char * returnedSupportDescriptionChar;
+ int size = Supports.size();
+
+ if (size == 1)
+ {
+ MESSAGE(LOC <<" there is only one SUPPORT in the argument list, the method return a copy of this object !");
+ SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
+
+ returnedSupport = new SUPPORT(*obj);
+
+ int lenName = strlen((Supports[0]->getName()).c_str()) + 8 + 1;
+ int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 8 + 1;
+
+ returnedSupportNameChar = new char[lenName];
+ returnedSupportDescriptionChar = new char[lenDescription];
+
+ returnedSupportNameChar = strcpy(returnedSupportNameChar,"Copy of ");
+ returnedSupportNameChar = strcat(returnedSupportNameChar,(Supports[0]->getName()).c_str());
+ returnedSupportDescriptionChar = strcpy(returnedSupportDescriptionChar,"Copy of ");
+ returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
+ (Supports[0]->getDescription()).c_str());
+
+ returnedSupportName = string(returnedSupportNameChar);
+ returnedSupportDescription = string(returnedSupportDescriptionChar);
+
+ returnedSupport->setName(returnedSupportName);
+ returnedSupport->setDescription(returnedSupportDescription);
+
+ delete [] returnedSupportNameChar;
+ delete [] returnedSupportDescriptionChar;
+ }
+ else
+ {
+ SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
+ returnedSupport = new SUPPORT(*obj);
+
+ int lenName = strlen((Supports[0]->getName()).c_str()) + 9 + 1;
+ int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 9 + 1;
+
+ for(int i = 1;i<size;i++)
+ {
+ obj = const_cast <SUPPORT *> (Supports[i]);
+ returnedSupport->blending(obj);
+
+ if (i == (size-1))
+ {
+ lenName = lenName + 5 + strlen((Supports[i]->getName()).c_str());
+ lenDescription = lenDescription + 5 +
+ strlen((Supports[i]->getDescription()).c_str());
+ }
+ else
+ {
+ lenName = lenName + 2 + strlen((Supports[i]->getName()).c_str());
+ lenDescription = lenDescription + 2 +
+ strlen((Supports[i]->getDescription()).c_str());
+ }
+ }
+
+ returnedSupportNameChar = new char[lenName];
+ returnedSupportDescriptionChar = new char[lenDescription];
+
+ returnedSupportNameChar = strcpy(returnedSupportNameChar,"Merge of ");
+ returnedSupportDescriptionChar = strcpy(returnedSupportDescriptionChar,"Merge of ");
+
+ returnedSupportNameChar = strcat(returnedSupportNameChar,(Supports[0]->getName()).c_str());
+ returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
+ (Supports[0]->getDescription()).c_str());
+
+ for(int i = 1;i<size;i++)
+ {
+ if (i == (size-1))
+ {
+ returnedSupportNameChar = strcat(returnedSupportNameChar," and ");
+ returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar," and ");
+
+ returnedSupportNameChar = strcat(returnedSupportNameChar,
+ (Supports[i]->getName()).c_str());
+ returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
+ (Supports[i]->getDescription()).c_str());
+ }
+ else
+ {
+ returnedSupportNameChar = strcat(returnedSupportNameChar,", ");
+ returnedSupportNameChar = strcat(returnedSupportNameChar,
+ (Supports[i]->getName()).c_str());
+
+ returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,", ");
+ returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
+ (Supports[i]->getDescription()).c_str());
+ }
+ }
+
+ returnedSupportName = string(returnedSupportNameChar);
+ returnedSupport->setName(returnedSupportName);
+
+ returnedSupportDescription = string(returnedSupportDescriptionChar);
+ returnedSupport->setDescription(returnedSupportDescription);
+
+ delete [] returnedSupportNameChar;
+ delete [] returnedSupportDescriptionChar;
+ }
+
+ END_OF(LOC) ;
+ return returnedSupport;
+}
+
+/*!
+ return a SUPPORT pointer on the intersection of all SUPPORTs in Supports.
+ The (SUPPORT *) NULL pointer is returned if the intersection is empty.
+ You should delete this pointer after use to avois memory leaks.
+*/
+SUPPORT * MESH::intersectSupports(const vector<SUPPORT *> Supports) const throw (MEDEXCEPTION)
+{
+ const char * LOC = "MESH:::intersectSupports(const vector<SUPPORT *> ) : " ;
+ BEGIN_OF(LOC) ;
+
+ SUPPORT * returnedSupport;
+ string returnedSupportName;
+ string returnedSupportDescription;
+ char * returnedSupportNameChar;
+ char * returnedSupportDescriptionChar;
+ int size = Supports.size();
+
+ if (size == 1)
+ {
+ MESSAGE(LOC <<" there is only one SUPPORT in the argument list, the method return a copy of this object !");
+ SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
+
+ returnedSupport = new SUPPORT(*obj);
+
+ int lenName = strlen((Supports[0]->getName()).c_str()) + 8 + 1;
+ int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 8 + 1;
+
+ returnedSupportNameChar = new char[lenName];
+ returnedSupportDescriptionChar = new char[lenDescription];
+
+ returnedSupportNameChar = strcpy(returnedSupportNameChar,"Copy of ");
+ returnedSupportNameChar = strcat(returnedSupportNameChar,(Supports[0]->getName()).c_str());
+ returnedSupportDescriptionChar = strcpy(returnedSupportDescriptionChar,"Copy of ");
+ returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
+ (Supports[0]->getDescription()).c_str());
+
+ returnedSupportName = string(returnedSupportNameChar);
+ returnedSupportDescription = string(returnedSupportDescriptionChar);
+
+ returnedSupport->setName(returnedSupportName);
+ returnedSupport->setDescription(returnedSupportDescription);
+
+ delete [] returnedSupportNameChar;
+ delete [] returnedSupportDescriptionChar;
+ }
+ else
+ {
+ SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
+ returnedSupport = new SUPPORT(*obj);
+
+ int lenName = strlen((Supports[0]->getName()).c_str()) + 16 + 1;
+ int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 16 + 1;
+
+ for(int i = 1;i<size;i++)
+ {
+ obj = const_cast <SUPPORT *> (Supports[i]);
+ returnedSupport->intersecting(obj);
+
+ if (i == (size-1))
+ {
+ lenName = lenName + 5 + strlen((Supports[i]->getName()).c_str());
+ lenDescription = lenDescription + 5 +
+ strlen((Supports[i]->getDescription()).c_str());
+ }
+ else
+ {
+ lenName = lenName + 2 + strlen((Supports[i]->getName()).c_str());
+ lenDescription = lenDescription + 2 +
+ strlen((Supports[i]->getDescription()).c_str());
+ }
+ }
+
+ if(returnedSupport != (SUPPORT *) NULL)
+ {
+ returnedSupportNameChar = new char[lenName];
+ returnedSupportDescriptionChar = new char[lenDescription];
+
+ returnedSupportNameChar = strcpy(returnedSupportNameChar,
+ "Intersection of ");
+ returnedSupportDescriptionChar =
+ strcpy(returnedSupportDescriptionChar,"Intersection of ");
+
+ returnedSupportNameChar = strcat(returnedSupportNameChar,
+ (Supports[0]->getName()).c_str());
+ returnedSupportDescriptionChar =
+ strcat(returnedSupportDescriptionChar,
+ (Supports[0]->getDescription()).c_str());
+
+ for(int i = 1;i<size;i++)
+ {
+ if (i == (size-1))
+ {
+ returnedSupportNameChar = strcat(returnedSupportNameChar,
+ " and ");
+ returnedSupportDescriptionChar =
+ strcat(returnedSupportDescriptionChar," and ");
+
+ returnedSupportNameChar =
+ strcat(returnedSupportNameChar,
+ (Supports[i]->getName()).c_str());
+ returnedSupportDescriptionChar =
+ strcat(returnedSupportDescriptionChar,
+ (Supports[i]->getDescription()).c_str());
+ }
+ else
+ {
+ returnedSupportNameChar = strcat(returnedSupportNameChar,
+ ", ");
+ returnedSupportNameChar =
+ strcat(returnedSupportNameChar,
+ (Supports[i]->getName()).c_str());
+
+ returnedSupportDescriptionChar =
+ strcat(returnedSupportDescriptionChar,", ");
+ returnedSupportDescriptionChar =
+ strcat(returnedSupportDescriptionChar,
+ (Supports[i]->getDescription()).c_str());
+ }
+ }
+
+ returnedSupportName = string(returnedSupportNameChar);
+ returnedSupport->setName(returnedSupportName);
+
+ returnedSupportDescription = string(returnedSupportDescriptionChar);
+ returnedSupport->setDescription(returnedSupportDescription);
+
+ delete [] returnedSupportNameChar;
+ delete [] returnedSupportDescriptionChar;
+ }
+ }
+
+ END_OF(LOC) ;
+ return returnedSupport;
+}
+
+
+// internal helper type
+struct _cell
+{
+ //int number;
+ std::vector<int> groups;
+ MED_EN::medGeometryElement geometricType;
+};
+
+// Create families from groups
+void MESH::createFamilies()
+{
+ int nbFam=0; // count the families we create, used as identifier ???????????
+
+ int idFamNode = 0; // identifier for node families
+ int idFamElement = 0; // identifier for cell, face or edge families
+
+ // Main loop on mesh's entities
+ for (medEntityMesh entity=MED_CELL; entity!=MED_ALL_ENTITIES; ++entity)
+ {
+ int numberofgroups = getNumberOfGroups(entity);
+ if(!numberofgroups)
+ continue; // no groups for this entity
+
+ vector< vector<FAMILY*> > whichFamilyInGroup(numberofgroups); // this container is used to update groups at the end
+
+ // make myFamilies points to the member corresponding to entity
+ vector<FAMILY*>* myFamilies;
+ switch ( entity )
+ {
+ case MED_CELL :
+ myFamilies = & _familyCell;
+ break;
+ case MED_FACE :
+ myFamilies = & _familyFace;
+ break;
+ case MED_EDGE :
+ myFamilies = & _familyEdge;
+ break;
+ case MED_NODE :
+ myFamilies = & _familyNode;
+ break;
+ }
+
+ vector<GROUP*> myGroups=getGroups(entity); // get a copy of the groups ptr for the entity
+ // get a copy of the (old) family ptrs before clearing
+ vector<FAMILY*> myOldFamilies=getFamilies(entity);
+ myFamilies->clear();
+
+
+ // 1 - Create a vector containing for each cell (of the entity) an information structure
+ // giving geometric type and the groups it belong to
+
+ med_int numberOfTypes=getNumberOfTypes(entity);
+ const int * index=getGlobalNumberingIndex(entity);
+ const medGeometryElement* geometricTypes=_connectivity->getGeometricTypes(entity); // pb avec entity=MED_NODE???
+ med_int numberOfCells=index[numberOfTypes]-1; // total number of cells for that entity
+ SCRUTE(numberOfTypes);
+ SCRUTE(numberOfCells);
+ vector< _cell > tab_cell(numberOfCells);
+ for(med_int t=0; t!=numberOfTypes; ++t)
+ for(int n=index[t]-1; n!=index[t+1]-1; ++n)
+ tab_cell[n].geometricType=geometricTypes[t];
+
+
+ // 2 - Scan cells in groups and update in tab_cell the container of groups a cell belong to
+
+ for (unsigned g=0; g!=myGroups.size(); ++g)
+ {
+ // scan cells that belongs to the group
+ const int* groupCells=myGroups[g]->getnumber()->getValue();
+ int nbCells=myGroups[g]->getnumber()->getLength();
+ for(int c=0; c!=nbCells; ++c)
+ tab_cell[groupCells[c]-1].groups.push_back(g);
+ }
+
+
+ // 3 - Scan the cells vector, genarate family name, and create a map associating the family names
+ // whith the vector of contained cells
+
+ map< string,vector<int> > tab_families;
+ map< string,vector<int> >::iterator fam;
+ for(int n=0; n!=numberOfCells; ++n)
+ {
+ ostringstream key; // to generate the name of the family
+ key << "FAM";
+ if(tab_cell[n].groups.empty()) // this cell don't belong to any group
+ key << "_NONE" << entity;
+
+ for(vector<int>::const_iterator it=tab_cell[n].groups.begin(); it!=tab_cell[n].groups.end(); ++it)
+ {
+ string groupName=myGroups[*it]->getName();
+ if(groupName.empty())
+ key << "_G" << *it;
+ else
+ key << "_" << groupName;
+ }
+
+ tab_families[key.str()].push_back(n+1); // fill the vector of contained cells associated whith the family
+ /* fam = tab_families.find(key.str());
+ if( fam != tab_families.end())
+ fam->second.push_back(n+1); // +1 : convention Fortran de MED
+ else
+ {
+ vector<int> newfamily;
+ newfamily.push_back(n+1); // +1 : convention Fortran de MED
+ tab_families.insert(make_pair(key.str(),newfamily));
+ }*/
+ }
+
+
+ // 4 - Scan the family map, create MED Families, check if it already exist.
+
+ for( fam=tab_families.begin(); fam!=tab_families.end(); ++fam)
+ {
+ vector<medGeometryElement> tab_types_geometriques;
+ medGeometryElement geometrictype=MED_NONE;
+ vector<int> tab_index_types_geometriques;
+ vector<int> tab_nombres_elements;
+
+ // scan family cells and fill the tab that are needed by the create a MED FAMILY
+ for( int i=0; i!=fam->second.size(); ++i)
+ {
+ int ncell=fam->second[i]-1;
+ if(tab_cell[ncell].geometricType != geometrictype)
+ {
+ // new geometric type -> we store it and complete index tabs
+ if(!tab_index_types_geometriques.empty())
+ tab_nombres_elements.push_back(i+1-tab_index_types_geometriques.back());
+ tab_types_geometriques.push_back( (geometrictype=tab_cell[ncell].geometricType));
+ tab_index_types_geometriques.push_back(i+1);
+ }
+ }
+ // store and complete index tabs for the last geometric type
+ tab_nombres_elements.push_back(fam->second.size()+1-tab_index_types_geometriques.back());
+ tab_index_types_geometriques.push_back(fam->second.size()+1);
+
+ // create a empty MED FAMILY and fill it with the tabs we constructed
+ FAMILY* newFam = new FAMILY();
+ newFam->setTotalNumberOfElements(fam->second.size());
+ newFam->setName(fam->first);
+ newFam->setMesh(this);
+ newFam->setNumberOfGeometricType(tab_types_geometriques.size());
+ newFam->setGeometricType(&tab_types_geometriques[0]); // we know the tab is not empy
+ newFam->setNumberOfElements(&tab_nombres_elements[0]);
+ newFam->setNumber(&tab_index_types_geometriques[0],&fam->second[0]);
+ newFam->setEntity(entity);
+ newFam->setAll(false);
+
+ int idFam = 0;
+
+ switch ( entity )
+ {
+ case MED_NODE :
+ ++idFamNode;
+ idFam = idFamNode;
+ break;
+ case MED_CELL:
+ ++idFamElement;
+ idFam = -idFamElement;
+ break;
+ case MED_FACE :
+ ++idFamElement;
+ idFam = -idFamElement;
+ break;
+ case MED_EDGE :
+ ++idFamElement;
+ idFam = -idFamElement;
+ break;
+ }
+
+ newFam->setIdentifier(idFam);
+
+ // Update links between families and groups
+
+ int ncell1=fam->second[0]-1; // number of first cell in family
+ int numberOfGroups=tab_cell[ncell1].groups.size(); // number of groups in the family
+ if(numberOfGroups)
+ {
+ newFam->setNumberOfGroups(numberOfGroups);
+ string * groupNames=new string[numberOfGroups];
+
+ // iterate on the groups the family belongs to
+ vector<int>::const_iterator it=tab_cell[ncell1].groups.begin();
+ for(int ng=0 ; it!=tab_cell[ncell1].groups.end(); ++it, ++ng)
+ {
+ whichFamilyInGroup[*it].push_back(newFam);
+ groupNames[ng]=myGroups[*it]->getName();
+ }
+ newFam->setGroupsNames(groupNames);
+ }
+
+ int sizeOfFamVect = myFamilies->size();
+
+ MESSAGE(" MESH::createFamilies() entity " << entity << " size " << sizeOfFamVect);
+
+ myFamilies->push_back(newFam);
+ }
+
+ // delete old families
+ for (int i=0;i<myOldFamilies.size();i++)
+ delete myOldFamilies[i] ;
+
+ // update references in groups
+ for (int i=0;i<myGroups.size();i++)
+ {
+ myGroups[i]->setNumberOfFamilies(whichFamilyInGroup[i].size());
+ myGroups[i]->setFamilies(whichFamilyInGroup[i]);
+ }
+
+ // re-scan the cells vector, and fill the family vector with cells.
+ // creation of support, check if it already exist.
+ }
+}
