12 #include "MEDMEM_DriversDef.hxx"
13 #include "MEDMEM_Field.hxx"
14 #include "MEDMEM_Mesh.hxx"
16 #include "MEDMEM_Support.hxx"
17 #include "MEDMEM_Family.hxx"
18 #include "MEDMEM_Group.hxx"
19 #include "MEDMEM_Coordinate.hxx"
20 #include "MEDMEM_Connectivity.hxx"
21 #include "MEDMEM_CellModel.hxx"
22 //#include "MEDMEM_Grid.hxx" this inclision should have never be here !!!
24 //update Families with content list
25 //int family_count(int family_number, int count, int * entities_number, int * entities_list) ;
27 // ------- Drivers Management Part
29 // MESH::INSTANCE_DE<MED_MESH_RDONLY_DRIVER> MESH::inst_med_rdonly ;
30 //const MESH::INSTANCE * const MESH::instances[] = { &MESH::inst_med_rdonly , &MESH::inst_med_rdwr } ;
32 // Add a similar line for your personnal driver (step 3)
34 MESH::INSTANCE_DE<MED_MESH_RDWR_DRIVER> MESH::inst_med ;
35 MESH::INSTANCE_DE<GIBI_MESH_RDWR_DRIVER> MESH::inst_gibi ;
36 MESH::INSTANCE_DE<VTK_MESH_DRIVER> MESH::inst_vtk;
38 // Add your own driver in the driver list (step 4)
39 // Note the list must be coherent with the driver type list defined in MEDMEM_DRIVER.hxx.
40 const MESH::INSTANCE * const MESH::instances[] = { &MESH::inst_med, &MESH::inst_gibi, &MESH::inst_vtk } ;
42 /*! Add a MESH driver of type (MED_DRIVER, ....) associated with file <fileName>. The meshname used in the file
43 is <driverName>. addDriver returns an int handler. */
44 int MESH::addDriver(driverTypes driverType,
45 const string & fileName/*="Default File Name.med"*/,const string & driverName/*="Default Mesh Name"*/) {
47 const char * LOC = "MESH::addDriver(driverTypes driverType, const string & fileName=\"Default File Name.med\",const string & driverName=\"Default Mesh Name\") : ";
51 int itDriver = (int) NO_DRIVER;
57 SCRUTE(instances[driverType]);
62 itDriver = (int) driverType ;
67 itDriver = (int) driverType ;
77 throw MED_EXCEPTION (LOCALIZED(STRING(LOC)<< "NO_DRIVER has been specified to the method which is not allowed"));
81 if (itDriver == ((int) NO_DRIVER))
82 throw MED_EXCEPTION (LOCALIZED(STRING(LOC)<< "othe driver than MED_DRIVER GIBI_DRIVER and VT_DRIVER has been specified to the method which is not allowed"));
84 driver = instances[itDriver]->run(fileName, this) ;
86 _drivers.push_back(driver);
88 current = _drivers.size()-1;
90 _drivers[current]->setMeshName(driverName);
97 /*! Add an existing MESH driver. */
98 int MESH::addDriver(GENDRIVER & driver) {
99 const char * LOC = "MESH::addDriver(GENDRIVER &) : ";
102 // A faire : Vérifier que le driver est de type MESH.
103 GENDRIVER * newDriver = driver.copy() ;
105 _drivers.push_back(newDriver);
106 return _drivers.size()-1;
111 /*! Remove an existing MESH driver. */
112 void MESH::rmDriver (int index/*=0*/) {
113 const char * LOC = "MESH::rmDriver (int index=0): ";
116 if ( _drivers[index] ) {
117 //_drivers.erase(&_drivers[index]);
119 MESSAGE ("detruire");
122 throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
123 << "The index given is invalid, index must be between 0 and |"
132 // ------ End of Drivers Management Part
137 const char * LOC = "MESH::init(): ";
141 string _name = "NOT DEFINED"; // A POSITIONNER EN FCT DES IOS ?
143 _coordinate = (COORDINATE *) NULL;
144 _connectivity = (CONNECTIVITY *) NULL;
146 _spaceDimension = MED_INVALID; // 0 ?!?
147 _meshDimension = MED_INVALID;
148 _numberOfNodes = MED_INVALID;
150 _numberOfNodesFamilies = 0; // MED_INVALID ?!?
151 _numberOfCellsFamilies = 0;
152 _numberOfFacesFamilies = 0;
153 _numberOfEdgesFamilies = 0;
155 _numberOfNodesGroups = 0; // MED_INVALID ?!?
156 _numberOfCellsGroups = 0;
157 _numberOfFacesGroups = 0;
158 _numberOfEdgesGroups = 0;
165 /*! Create an empty MESH. */
166 MESH::MESH():_coordinate(NULL),_connectivity(NULL), _isAGrid(false) {
173 _isAGrid = m._isAGrid;
175 if (m._coordinate != NULL)
176 _coordinate = new COORDINATE(* m._coordinate);
178 _coordinate = (COORDINATE *) NULL;
180 if (m._connectivity != NULL)
181 _connectivity = new CONNECTIVITY(* m._connectivity);
183 _connectivity = (CONNECTIVITY *) NULL;
185 _spaceDimension = m._spaceDimension;
186 _meshDimension = m._meshDimension;
187 _numberOfNodes = m._numberOfNodes;
189 _numberOfNodesFamilies = m._numberOfNodesFamilies;
190 _familyNode = m._familyNode;
191 for (int i=0; i<m._numberOfNodesFamilies; i++)
193 _familyNode[i] = new FAMILY(* m._familyNode[i]);
194 _familyNode[i]->setMesh(this);
197 _numberOfCellsFamilies = m._numberOfCellsFamilies;
198 _familyCell = m._familyCell;
199 for (int i=0; i<m._numberOfCellsFamilies; i++)
201 _familyCell[i] = new FAMILY(* m._familyCell[i]);
202 _familyCell[i]->setMesh(this);
205 _numberOfFacesFamilies = m._numberOfFacesFamilies;
206 _familyFace = m._familyFace;
207 for (int i=0; i<m._numberOfFacesFamilies; i++)
209 _familyFace[i] = new FAMILY(* m._familyFace[i]);
210 _familyFace[i]->setMesh(this);
213 _numberOfEdgesFamilies = m._numberOfEdgesFamilies;
214 _familyEdge = m._familyEdge;
215 for (int i=0; i<m._numberOfEdgesFamilies; i++)
217 _familyEdge[i] = new FAMILY(* m._familyEdge[i]);
218 _familyEdge[i]->setMesh(this);
221 _numberOfNodesGroups = m._numberOfNodesGroups;
222 _groupNode = m._groupNode;
223 for (int i=0; i<m._numberOfNodesGroups; i++)
225 _groupNode[i] = new GROUP(* m._groupNode[i]);
226 _groupNode[i]->setMesh(this);
229 _numberOfCellsGroups = m._numberOfCellsGroups;
230 _groupCell = m._groupCell;
231 for (int i=0; i<m._numberOfCellsGroups; i++)
233 _groupCell[i] = new GROUP(* m._groupCell[i]);
234 _groupCell[i]->setMesh(this);
237 _numberOfFacesGroups = m._numberOfFacesGroups;
238 _groupFace = m._groupFace;
239 for (int i=0; i<m._numberOfFacesGroups; i++)
241 _groupFace[i] = new GROUP(* m._groupFace[i]);
242 _groupFace[i]->setMesh(this);
245 _numberOfEdgesGroups = m._numberOfEdgesGroups;
246 _groupEdge = m._groupEdge;
247 for (int i=0; i<m._numberOfEdgesGroups; i++)
249 _groupEdge[i] = new GROUP(* m._groupEdge[i]);
250 _groupEdge[i]->setMesh(this);
253 //_drivers = m._drivers; //Recopie des drivers?
258 MESSAGE("MESH::~MESH() : Destroying the Mesh");
259 if (_coordinate != ((COORDINATE *) NULL)) delete _coordinate ;
260 if (_connectivity != ((CONNECTIVITY *) NULL)) delete _connectivity ;
262 size = _familyNode.size() ;
263 for (int i=0;i<size;i++)
264 delete _familyNode[i] ;
265 size = _familyCell.size() ;
266 for (int i=0;i<size;i++)
267 delete _familyCell[i] ;
268 size = _familyFace.size() ;
269 for (int i=0;i<size;i++)
270 delete _familyFace[i] ;
271 size = _familyEdge.size() ;
272 for (int i=0;i<size;i++)
273 delete _familyEdge[i] ;
274 size = _groupNode.size() ;
275 for (int i=0;i<size;i++)
276 delete _groupNode[i] ;
277 size = _groupCell.size() ;
278 for (int i=0;i<size;i++)
279 delete _groupCell[i] ;
280 size = _groupFace.size() ;
281 for (int i=0;i<size;i++)
282 delete _groupFace[i] ;
283 size = _groupEdge.size() ;
284 for (int i=0;i<size;i++)
285 delete _groupEdge[i] ;
287 MESSAGE("In this object MESH there is(are) " << _drivers.size() << " driver(s)");
289 for (unsigned int index=0; index < _drivers.size(); index++ )
291 SCRUTE(_drivers[index]);
292 if ( _drivers[index] != NULL) delete _drivers[index];
297 MESH & MESH::operator=(const MESH &m)
299 const char * LOC = "MESH & MESH::operator=(const MESH &m) : ";
302 MESSAGE(LOC <<"Not yet implemented, operating on the object " << m);
305 // ATTENTION CET OPERATEUR DE RECOPIE EST DANGEREUX POUR LES
306 // POINTEURS : ex : nodal_connectivity ???? EXPRES ????
307 // _drivers = m._drivers;
309 // space_dimension=m.space_dimension;
310 // mesh_dimension=m.mesh_dimension;
312 // nodes_count=m.nodes_count;
313 // coordinates=m.coordinates;
314 // coordinates_name=m.coordinates_name ;
315 // coordinates_unit=m.coordinates_unit ;
316 // nodes_numbers=m.nodes_numbers ;
318 // cells_types_count=m.cells_types_count;
319 // cells_type=m.cells_type;
320 // cells_count=m.cells_count;
321 // nodal_connectivity=m.nodal_connectivity;
323 // nodes_families_count=m.nodes_families_count;
324 // nodes_Families=m.nodes_Families;
326 // cells_families_count=m.cells_families_count;
327 // cells_Families=m.cells_Families;
329 // maximum_cell_number_by_node = m.maximum_cell_number_by_node;
330 // if (maximum_cell_number_by_node > 0)
332 // reverse_nodal_connectivity = m.reverse_nodal_connectivity;
333 // reverse_nodal_connectivity_index = m.reverse_nodal_connectivity_index ;
340 /*! Create a MESH object using a MESH driver of type (MED_DRIVER, ....) associated with file <fileName>.
341 The meshname <driverName> must already exists in the file.*/
342 MESH::MESH(driverTypes driverType, const string & fileName/*=""*/, const string & driverName/*=""*/) throw (MEDEXCEPTION)
344 const char * LOC ="MESH::MESH(driverTypes driverType, const string & fileName="", const string & driverName="") : ";
355 MED_MESH_RDONLY_DRIVER myDriver(fileName,this);
356 myDriver.setMeshName(driverName);
357 current = addDriver(myDriver);
362 GIBI_MESH_RDONLY_DRIVER myDriver(fileName,this);
363 current = addDriver(myDriver);
368 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Driver type unknown : can't create driver!"));
373 // current = addDriver(driverType,fileName,driverName);
374 // switch(_drivers[current]->getAccessMode() ) {
375 // case MED_WRONLY : {
376 // MESSAGE("MESH::MESH(driverTypes driverType, .....) : driverType must not be a MED_WRONLY accessMode");
377 // rmDriver(current);
382 _drivers[current]->open();
383 _drivers[current]->read();
384 _drivers[current]->close();
387 // ((GRID *) this)->fillMeshAfterRead();
393 // Node MESH::Donne_Barycentre(const Maille &m) const
395 // Node N=node[m[0]];
396 // for (int i=1;i<m.donne_nbr_sommets();i++) N=N+node[m[i]];
397 // N=N*(1.0/m.donne_nbr_sommets());
401 ostream & operator<<(ostream &os, const MESH &myMesh)
403 int spacedimension = myMesh.getSpaceDimension();
404 int meshdimension = myMesh.getMeshDimension();
405 int numberofnodes = myMesh.getNumberOfNodes();
407 os << "Space Dimension : " << spacedimension << endl << endl;
409 os << "Mesh Dimension : " << meshdimension << endl << endl;
411 const double * coordinates = myMesh.getCoordinates(MED_FULL_INTERLACE);
412 os << "SHOW NODES COORDINATES : " << endl;
414 os << "Name :" << endl;
415 const string * coordinatesnames = myMesh.getCoordinatesNames();
416 for(int i=0; i<spacedimension ; i++)
418 os << " - " << coordinatesnames[i] << endl;
420 os << "Unit :" << endl;
421 const string * coordinatesunits = myMesh.getCoordinatesUnits();
422 for(int i=0; i<spacedimension ; i++)
424 os << " - " << coordinatesunits[i] << endl;
426 for(int i=0; i<numberofnodes ; i++)
428 os << "Nodes " << i+1 << " : ";
429 for (int j=0; j<spacedimension ; j++)
430 os << coordinates[i*spacedimension+j] << " ";
434 os << endl << "SHOW CONNECTIVITY :" << endl;
435 os << *myMesh._connectivity << endl;
437 medEntityMesh entity;
438 os << endl << "SHOW FAMILIES :" << endl << endl;
439 for (int k=1; k<=4; k++)
441 if (k==1) entity = MED_NODE;
442 if (k==2) entity = MED_CELL;
443 if (k==3) entity = MED_FACE;
444 if (k==4) entity = MED_EDGE;
445 int numberoffamilies = myMesh.getNumberOfFamilies(entity);
446 using namespace MED_FR;
447 os << "NumberOfFamilies on "<< entNames[(MED_FR::med_entite_maillage) entity] <<" : "<<numberoffamilies<<endl;
448 using namespace MED_EN;
449 for (int i=1; i<numberoffamilies+1;i++)
451 os << * myMesh.getFamily(entity,i) << endl;
455 os << endl << "SHOW GROUPS :" << endl << endl;
456 for (int k=1; k<=4; k++)
458 if (k==1) entity = MED_NODE;
459 if (k==2) entity = MED_CELL;
460 if (k==3) entity = MED_FACE;
461 if (k==4) entity = MED_EDGE;
462 int numberofgroups = myMesh.getNumberOfGroups(entity);
463 using namespace MED_FR;
464 os << "NumberOfGroups on "<< entNames[(MED_FR::med_entite_maillage) entity] <<" : "<<numberofgroups<<endl;
465 using namespace MED_EN;
466 for (int i=1; i<numberofgroups+1;i++)
468 os << * myMesh.getGroup(entity,i) << endl;
476 Get global number of element which have same connectivity than connectivity argument.
478 It do not take care of connectivity order (3,4,7,10 is same as 7,3,10,4).
480 Return -1 if not found.
482 int MESH::getElementNumber(medConnectivity ConnectivityType, medEntityMesh Entity, medGeometryElement Type, int * connectivity) const
484 const char* LOC="MESH::getElementNumber " ;
487 int numberOfValue ; // size of connectivity array
488 CELLMODEL myType(Type) ;
489 if (ConnectivityType==MED_DESCENDING)
490 numberOfValue = myType.getNumberOfConstituents(1) ; // faces (3D) or edges (2D)
492 numberOfValue = myType.getNumberOfNodes() ; // nodes
494 const int * myReverseConnectivityValue = getReverseConnectivity(ConnectivityType,Entity) ;
495 const int * myReverseConnectivityIndex = getReverseConnectivityIndex(ConnectivityType,Entity) ;
497 // First node or face/edge
498 int indexBegin = myReverseConnectivityIndex[connectivity[0]-1] ;
499 int indexEnd = myReverseConnectivityIndex[connectivity[0]] ;
501 // list to put cells numbers
502 list<int> cellsList ;
503 list<int>::iterator itList ;
504 for (int i=indexBegin; i<indexEnd; i++)
505 cellsList.push_back(myReverseConnectivityValue[i-1]) ;
507 // Foreach node or face/edge in cell, we search cells which are in common.
508 // At the end, we must have only one !
510 for (int i=1; i<numberOfValue; i++) {
511 int connectivity_i = connectivity[i] ;
512 for (itList=cellsList.begin();itList!=cellsList.end();itList++) {
514 for (int j=myReverseConnectivityIndex[connectivity_i-1]; j<myReverseConnectivityIndex[connectivity_i]; j++) {
515 if ((*itList)==myReverseConnectivityValue[j-1]) {
521 itList=cellsList.erase(itList);
522 itList--; // well : rigth if itList = cellsList.begin() ??
527 if (cellsList.size()>1) // we have more than one cell
528 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Found more than one element !")) ;
530 if (cellsList.size()==0)
535 return cellsList.front() ;
540 Return a support which reference all elements on the boundary of mesh.
542 For instance, we could get only face in 3D and edge in 2D.
544 SUPPORT * MESH::getBoundaryElements(medEntityMesh Entity)
547 const char * LOC = "MESH::getBoundaryElements : " ;
550 // actually we could only get face (in 3D) and edge (in 2D)
551 if (_spaceDimension == 3)
552 if (Entity != MED_FACE)
553 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Not defined in 3D mesh for entity "<<Entity<<" !"));
554 if (_spaceDimension == 2)
555 if (Entity != MED_EDGE)
556 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Not defined in 2D mesh for entity "<<Entity<<" !"));
559 SUPPORT * mySupport = new SUPPORT(this,"Boundary",Entity);
560 //mySupport.setDescription("boundary of type ...");
561 mySupport->setAll(false);
564 const int * myConnectivityValue = getReverseConnectivity(MED_DESCENDING) ;
565 const int * myConnectivityIndex = getReverseConnectivityIndex(MED_DESCENDING) ;
566 int numberOf = getNumberOfElements(Entity,MED_ALL_ELEMENTS) ;
567 list<int> myElementsList ;
569 for (int i=0 ; i<numberOf; i++)
570 if (myConnectivityValue[myConnectivityIndex[i]] == 0) {
571 myElementsList.push_back(i+1) ;
574 // Well, we must know how many geometric type we have found
575 int * myListArray = new int[size] ;
577 list<int>::iterator myElementsListIt ;
578 for (myElementsListIt=myElementsList.begin();myElementsListIt!=myElementsList.end();myElementsListIt++) {
579 myListArray[id]=(*myElementsListIt) ;
583 int numberOfGeometricType ;
584 medGeometryElement* geometricType ;
585 int * numberOfGaussPoint ;
586 int * geometricTypeNumber ;
587 int * numberOfElements ;
588 //MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY() ;
589 int * mySkyLineArrayIndex ;
591 int numberOfType = getNumberOfTypes(Entity) ;
592 if (numberOfType == 1) { // wonderfull : it's easy !
593 numberOfGeometricType = 1 ;
594 geometricType = new medGeometryElement[1] ;
595 const medGeometryElement * allType = getTypes(Entity);
596 geometricType[0] = allType[0] ;
597 numberOfGaussPoint = new int[1] ;
598 numberOfGaussPoint[0] = 1 ;
599 geometricTypeNumber = new int[1] ; // not use, but initialized to nothing
600 geometricTypeNumber[0] = 0 ;
601 numberOfElements = new int[1] ;
602 numberOfElements[0] = size ;
603 mySkyLineArrayIndex = new int[2] ;
604 mySkyLineArrayIndex[0]=1 ;
605 mySkyLineArrayIndex[1]=1+size ;
608 map<medGeometryElement,int> theType ;
609 for (myElementsListIt=myElementsList.begin();myElementsListIt!=myElementsList.end();myElementsListIt++) {
610 medGeometryElement myType = getElementType(Entity,*myElementsListIt) ;
611 if (theType.find(myType) != theType.end() )
616 numberOfGeometricType = theType.size() ;
617 geometricType = new medGeometryElement[numberOfGeometricType] ;
618 //const medGeometryElement * allType = getTypes(Entity); !! UNUSZED VARIABLE !!
619 numberOfGaussPoint = new int[numberOfGeometricType] ;
620 geometricTypeNumber = new int[numberOfGeometricType] ; // not use, but initialized to nothing
621 numberOfElements = new int[numberOfGeometricType] ;
622 mySkyLineArrayIndex = new int[numberOfGeometricType+1] ;
624 mySkyLineArrayIndex[0]=1 ;
625 map<medGeometryElement,int>::iterator theTypeIt ;
626 for (theTypeIt=theType.begin();theTypeIt!=theType.end();theTypeIt++) {
627 geometricType[index] = (*theTypeIt).first ;
628 numberOfGaussPoint[index] = 1 ;
629 geometricTypeNumber[index] = 0 ;
630 numberOfElements[index] = (*theTypeIt).second ;
631 mySkyLineArrayIndex[index+1]=mySkyLineArrayIndex[index]+numberOfElements[index] ;
635 //mySkyLineArray->setMEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
636 MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
638 mySupport->setNumberOfGeometricType(numberOfGeometricType) ;
639 mySupport->setGeometricType(geometricType) ;
640 mySupport->setNumberOfGaussPoint(numberOfGaussPoint) ;
641 mySupport->setNumberOfElements(numberOfElements) ;
642 mySupport->setTotalNumberOfElements(size) ;
643 // mySupport->setGeometricTypeNumber(geometricTypeNumber) ;
644 mySupport->setNumber(mySkyLineArray) ;
646 delete[] numberOfElements;
647 delete[] geometricTypeNumber;
648 delete[] numberOfGaussPoint;
649 delete[] geometricType;
650 delete[] mySkyLineArrayIndex;
651 delete[] myListArray;
652 // delete mySkyLineArray;
658 FIELD<double>* MESH::getVolume(const SUPPORT * Support) const throw (MEDEXCEPTION)
660 const char * LOC = "MESH::getVolume(SUPPORT*) : ";
663 // Support must be on 3D elements
665 // Make sure that the MESH class is the same as the MESH class attribut
666 // in the class Support
667 MESH* myMesh = Support->getMesh();
669 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
671 int dim_space = getSpaceDimension();
672 medEntityMesh support_entity = Support->getEntity();
673 bool onAll = Support->isOnAllElements();
675 int nb_type, length_values;
676 const medGeometryElement* types;
678 // !!!! WARNING : use of nodal global numbering in the mesh !!!!
679 const int* global_connectivity;
683 // nb_type = myMesh->getNumberOfTypes(support_entity);
684 // length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
685 // types = getTypes(support_entity);
689 nb_type = Support->getNumberOfTypes();
690 length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
691 types = Support->getTypes();
695 FIELD<double>* Volume = new FIELD<double>::FIELD(Support,1);
696 // double *volume = new double [length_values];
697 Volume->setName("VOLUME");
698 Volume->setDescription("cells volume");
699 Volume->setComponentName(1,"volume");
700 Volume->setComponentDescription(1,"desc-comp");
702 /* string MEDComponentUnit(MED_TAILLE_PNOM,' ');*/
704 string MEDComponentUnit = myMesh->getCoordinatesUnits()[0]+"*"+myMesh->getCoordinatesUnits()[1]+"*"+myMesh->getCoordinatesUnits()[2];
706 Volume->setMEDComponentUnit(1,MEDComponentUnit);
708 Volume->setValueType(MED_REEL64);
710 Volume->setIterationNumber(0);
711 Volume->setOrderNumber(0);
712 Volume->setTime(0.0);
714 //const double *volume = Volume->getValue(MED_FULL_INTERLACE);
715 MEDARRAY<double> *volume = Volume->getvalue();
718 const double * coord = getCoordinates(MED_FULL_INTERLACE);
720 for (int i=0;i<nb_type;i++)
722 medGeometryElement type = types[i] ;
727 nb_entity_type = getNumberOfElements(support_entity,type);
728 global_connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,type);
732 // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all"));
734 nb_entity_type = Support->getNumberOfElements(type);
736 const int * supp_number = Support->getNumber(type);
737 const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
738 const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
739 int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
741 for (int k_type = 0; k_type<nb_entity_type; k_type++) {
742 for (int j_ent = 0; j_ent<(type%100); j_ent++) {
743 global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
746 global_connectivity = global_connectivity_tmp ;
751 case MED_TETRA4 : case MED_TETRA10 :
753 for (int tetra=0;tetra<nb_entity_type;tetra++)
755 int tetra_index = (type%100)*tetra;
757 int N1 = global_connectivity[tetra_index]-1;
758 int N2 = global_connectivity[tetra_index+1]-1;
759 int N3 = global_connectivity[tetra_index+2]-1;
760 int N4 = global_connectivity[tetra_index+3]-1;
762 double x1 = coord[dim_space*N1];
763 double x2 = coord[dim_space*N2];
764 double x3 = coord[dim_space*N3];
765 double x4 = coord[dim_space*N4];
767 double y1 = coord[dim_space*N1+1];
768 double y2 = coord[dim_space*N2+1];
769 double y3 = coord[dim_space*N3+1];
770 double y4 = coord[dim_space*N4+1];
772 double z1 = coord[dim_space*N1+2];
773 double z2 = coord[dim_space*N2+2];
774 double z3 = coord[dim_space*N3+2];
775 double z4 = coord[dim_space*N4+2];
777 xvolume = ((x3-x1)*((y2-y1)*(z4-z1) - (z2-z1)*(y4-y1)) -
778 (x2-x1)*((y3-y1)*(z4-z1) - (z3-z1)*(y4-y1)) +
779 (x4-x1)*((y3-y1)*(z2-z1) - (z3-z1)*(y2-y1)))/6.0;
781 //volume[index] = xvolume ;
782 volume->setIJ(index,1,xvolume) ;
787 case MED_PYRA5 : case MED_PYRA13 :
789 for (int pyra=0;pyra<nb_entity_type;pyra++)
791 int pyra_index = (type%100)*pyra;
793 int N1 = global_connectivity[pyra_index]-1;
794 int N2 = global_connectivity[pyra_index+1]-1;
795 int N3 = global_connectivity[pyra_index+2]-1;
796 int N4 = global_connectivity[pyra_index+3]-1;
797 int N5 = global_connectivity[pyra_index+4]-1;
799 double x1 = coord[dim_space*N1];
800 double x2 = coord[dim_space*N2];
801 double x3 = coord[dim_space*N3];
802 double x4 = coord[dim_space*N4];
803 double x5 = coord[dim_space*N5];
805 double y1 = coord[dim_space*N1+1];
806 double y2 = coord[dim_space*N2+1];
807 double y3 = coord[dim_space*N3+1];
808 double y4 = coord[dim_space*N4+1];
809 double y5 = coord[dim_space*N5+1];
811 double z1 = coord[dim_space*N1+2];
812 double z2 = coord[dim_space*N2+2];
813 double z3 = coord[dim_space*N3+2];
814 double z4 = coord[dim_space*N4+2];
815 double z5 = coord[dim_space*N5+2];
817 xvolume = (((x3-x1)*((y2-y1)*(z5-z1) - (z2-z1)*(y5-y1)) -
818 (x2-x1)*((y3-y1)*(z5-z1) - (z3-z1)*(y5-y1)) +
819 (x5-x1)*((y3-y1)*(z2-z1) - (z3-z1)*(y2-y1))) +
820 ((x4-x1)*((y3-y1)*(z5-z1) - (z3-z1)*(y5-y1)) -
821 (x3-x1)*((y4-y1)*(z5-z1) - (z4-z1)*(y5-y1)) +
822 (x5-x1)*((y4-y1)*(z3-z1) - (z4-z1)*(y3-y1)))
825 //volume[index] = xvolume ;
826 volume->setIJ(index,1,xvolume) ;
831 case MED_PENTA6 : case MED_PENTA15 :
833 for (int penta=0;penta<nb_entity_type;penta++)
835 int penta_index = (type%100)*penta;
837 int N1 = global_connectivity[penta_index]-1;
838 int N2 = global_connectivity[penta_index+1]-1;
839 int N3 = global_connectivity[penta_index+2]-1;
840 int N4 = global_connectivity[penta_index+3]-1;
841 int N5 = global_connectivity[penta_index+4]-1;
842 int N6 = global_connectivity[penta_index+5]-1;
844 double x1 = coord[dim_space*N1];
845 double x2 = coord[dim_space*N2];
846 double x3 = coord[dim_space*N3];
847 double x4 = coord[dim_space*N4];
848 double x5 = coord[dim_space*N5];
849 double x6 = coord[dim_space*N6];
851 double y1 = coord[dim_space*N1+1];
852 double y2 = coord[dim_space*N2+1];
853 double y3 = coord[dim_space*N3+1];
854 double y4 = coord[dim_space*N4+1];
855 double y5 = coord[dim_space*N5+1];
856 double y6 = coord[dim_space*N6+1];
858 double z1 = coord[dim_space*N1+2];
859 double z2 = coord[dim_space*N2+2];
860 double z3 = coord[dim_space*N3+2];
861 double z4 = coord[dim_space*N4+2];
862 double z5 = coord[dim_space*N5+2];
863 double z6 = coord[dim_space*N6+2];
865 double a1 = (x2-x3)/2.0, a2 = (y2-y3)/2.0, a3 = (z2-z3)/2.0;
866 double b1 = (x5-x6)/2.0, b2 = (y5-y6)/2.0, b3 = (z5-z6)/2.0;
867 double c1 = (x4-x1)/2.0, c2 = (y4-y1)/2.0, c3 = (z4-z1)/2.0;
868 double d1 = (x5-x2)/2.0, d2 = (y5-y2)/2.0, d3 = (z5-z2)/2.0;
869 double e1 = (x6-x3)/2.0, e2 = (y6-y3)/2.0, e3 = (z6-z3)/2.0;
870 double f1 = (x1-x3)/2.0, f2 = (y1-y3)/2.0, f3 = (z1-z3)/2.0;
871 double h1 = (x4-x6)/2.0, h2 = (y4-y6)/2.0, h3 = (z4-z6)/2.0;
873 double A = a1*c2*f3 - a1*c3*f2 - a2*c1*f3 + a2*c3*f1 +
875 double B = b1*c2*h3 - b1*c3*h2 - b2*c1*h3 + b2*c3*h1 +
877 double C = (a1*c2*h3 + b1*c2*f3) - (a1*c3*h2 + b1*c3*f2) -
878 (a2*c1*h3 + b2*c1*f3) + (a2*c3*h1 + b2*c3*f1) +
879 (a3*c1*h2 + b3*c1*f2) - (a3*c2*h1 + b3*c2*f1);
880 double D = a1*d2*f3 - a1*d3*f2 - a2*d1*f3 + a2*d3*f1 +
882 double E = b1*d2*h3 - b1*d3*h2 - b2*d1*h3 + b2*d3*h1 +
884 double F = (a1*d2*h3 + b1*d2*f3) - (a1*d3*h2 + b1*d3*f2) -
885 (a2*d1*h3 + b2*d1*f3) + (a2*d3*h1 + b2*d3*f1) +
886 (a3*d1*h2 + b3*d1*f2) - (a3*d2*h1 + b3*d2*f1);
887 double G = a1*e2*f3 - a1*e3*f2 - a2*e1*f3 + a2*e3*f1 +
889 double H = b1*e2*h3 - b1*e3*h2 - b2*e1*h3 + b2*e3*h1 +
891 double P = (a1*e2*h3 + b1*e2*f3) - (a1*e3*h2 + b1*e3*f2) -
892 (a2*e1*h3 + b2*e1*f3) + (a2*e3*h1 + b2*e3*f1) +
893 (a3*e1*h2 + b3*e1*f2) - (a3*e2*h1 + b3*e2*f1);
895 xvolume = -2.0*(2.0*(A + B + D + E + G + H) + C + F + P)/9.0;
897 //volume[index] = xvolume ;
898 volume->setIJ(index,1,xvolume) ;
903 case MED_HEXA8 : case MED_HEXA20 :
905 for (int hexa=0;hexa<nb_entity_type;hexa++)
907 int hexa_index = (type%100)*hexa;
909 int N1 = global_connectivity[hexa_index]-1;
910 int N2 = global_connectivity[hexa_index+1]-1;
911 int N3 = global_connectivity[hexa_index+2]-1;
912 int N4 = global_connectivity[hexa_index+3]-1;
913 int N5 = global_connectivity[hexa_index+4]-1;
914 int N6 = global_connectivity[hexa_index+5]-1;
915 int N7 = global_connectivity[hexa_index+6]-1;
916 int N8 = global_connectivity[hexa_index+7]-1;
918 double x1 = coord[dim_space*N1];
919 double x2 = coord[dim_space*N2];
920 double x3 = coord[dim_space*N3];
921 double x4 = coord[dim_space*N4];
922 double x5 = coord[dim_space*N5];
923 double x6 = coord[dim_space*N6];
924 double x7 = coord[dim_space*N7];
925 double x8 = coord[dim_space*N8];
927 double y1 = coord[dim_space*N1+1];
928 double y2 = coord[dim_space*N2+1];
929 double y3 = coord[dim_space*N3+1];
930 double y4 = coord[dim_space*N4+1];
931 double y5 = coord[dim_space*N5+1];
932 double y6 = coord[dim_space*N6+1];
933 double y7 = coord[dim_space*N7+1];
934 double y8 = coord[dim_space*N8+1];
936 double z1 = coord[dim_space*N1+2];
937 double z2 = coord[dim_space*N2+2];
938 double z3 = coord[dim_space*N3+2];
939 double z4 = coord[dim_space*N4+2];
940 double z5 = coord[dim_space*N5+2];
941 double z6 = coord[dim_space*N6+2];
942 double z7 = coord[dim_space*N7+2];
943 double z8 = coord[dim_space*N8+2];
945 double a1 = (x3-x4)/8.0, a2 = (y3-y4)/8.0, a3 = (z3-z4)/8.0;
946 double b1 = (x2-x1)/8.0, b2 = (y2-y1)/8.0, b3 = (z2-z1)/8.0;
947 double c1 = (x7-x8)/8.0, c2 = (y7-y8)/8.0, c3 = (z7-z8)/8.0;
948 double d1 = (x6-x5)/8.0, d2 = (y6-y5)/8.0, d3 = (z6-z5)/8.0;
949 double e1 = (x3-x2)/8.0, e2 = (y3-y2)/8.0, e3 = (z3-z2)/8.0;
950 double f1 = (x4-x1)/8.0, f2 = (y4-y1)/8.0, f3 = (z4-z1)/8.0;
951 double h1 = (x7-x6)/8.0, h2 = (y7-y6)/8.0, h3 = (z7-z6)/8.0;
952 double p1 = (x8-x5)/8.0, p2 = (y8-y5)/8.0, p3 = (z8-z5)/8.0;
953 double q1 = (x3-x7)/8.0, q2 = (y3-y7)/8.0, q3 = (z3-z7)/8.0;
954 double r1 = (x4-x8)/8.0, r2 = (y4-y8)/8.0, r3 = (z4-z8)/8.0;
955 double s1 = (x2-x6)/8.0, s2 = (y2-y6)/8.0, s3 = (z2-z6)/8.0;
956 double t1 = (x1-x5)/8.0, t2 = (y1-y5)/8.0, t3 = (z1-z5)/8.0;
958 double A = a1*e2*q3 - a1*e3*q2 - a2*e1*q3 + a2*e3*q1 +
960 double B = c1*h2*q3 - c1*h3*q2 - c2*h1*q3 + c2*h3*q1 +
962 double C = (a1*h2 + c1*e2)*q3 - (a1*h3 + c1*e3)*q2 -
963 (a2*h1 + c2*e1)*q3 + (a2*h3 + c2*e3)*q1 +
964 (a3*h1 + c3*e1)*q2 - (a3*h2 + c3*e2)*q1;
965 double D = b1*e2*s3 - b1*e3*s2 - b2*e1*s3 + b2*e3*s1 +
967 double E = d1*h2*s3 - d1*h3*s2 - d2*h1*s3 + d2*h3*s1 +
969 double F = (b1*h2 + d1*e2)*s3 - (b1*h3 + d1*e3)*s2 -
970 (b2*h1 + d2*e1)*s3 + (b2*h3 + d2*e3)*s1 +
971 (b3*h1 + d3*e1)*s2 - (b3*h2 + d3*e2)*s1;
972 double G = (a1*e2*s3 + b1*e2*q3) - (a1*e3*s2 + b1*e3*q2) -
973 (a2*e1*s3 + b2*e1*q3) + (a2*e3*s1 + b2*e3*q1) +
974 (a3*e1*s2 + b3*e1*q2) - (a3*e2*s1 + b3*e2*q1);
975 double H = (c1*h2*s3 + d1*h2*q3) - (c1*h3*s2 + d1*h3*q2) -
976 (c2*h1*s3 + d2*h1*q3) + (c2*h3*s1 + d2*h3*q1) +
977 (c3*h1*s2 + d3*h1*q2) - (c3*h2*s1 + d3*h2*q1);
978 double I = ((a1*h2 + c1*e2)*s3 + (b1*h2 + d1*e2)*q3) -
979 ((a1*h3 + c1*e3)*s2 + (b1*h3 + d1*e3)*q2) -
980 ((a2*h1 + c2*e1)*s3 + (b2*h1 + d2*e1)*q3) +
981 ((a2*h3 + c2*e3)*s1 + (b2*h3 + d2*e3)*q1) +
982 ((a3*h1 + c3*e1)*s2 + (b3*h1 + d3*e1)*q2) -
983 ((a3*h2 + c3*e2)*s1 + (b3*h2 + d3*e2)*q1);
984 double J = a1*f2*r3 - a1*f3*r2 - a2*f1*r3 + a2*f3*r1 +
986 double K = c1*p2*r3 - c1*p3*r2 - c2*p1*r3 + c2*p3*r1 +
988 double L = (a1*p2 + c1*f2)*r3 - (a1*p3 + c1*f3)*r2 -
989 (a2*p1 + c2*f1)*r3 + (a2*p3 + c2*f3)*r1 +
990 (a3*p1 + c3*f1)*r2 - (a3*p2 + c3*f2)*r1;
991 double M = b1*f2*t3 - b1*f3*t2 - b2*f1*t3 + b2*f3*t1 +
993 double N = d1*p2*t3 - d1*p3*t2 - d2*p1*t3 + d2*p3*t1 +
995 double O = (b1*p2 + d1*f2)*t3 - (b1*p3 + d1*f3)*t2 -
996 (b2*p1 + d2*f1)*t3 + (b2*p3 + d2*f3)*t1 +
997 (b3*p1 + d3*f1)*t2 - (b3*p2 + d3*f2)*t1;
998 double P = (a1*f2*t3 + b1*f2*r3) - (a1*f3*t2 + b1*f3*r2) -
999 (a2*f1*t3 + b2*f1*r3) + (a2*f3*t1 + b2*f3*r1) +
1000 (a3*f1*t2 + b3*f1*r2) - (a3*f2*t1 + b3*f2*r1);
1001 double Q = (c1*p2*t3 + d1*p2*r3) - (c1*p3*t2 + d1*p3*r2) -
1002 (c2*p1*t3 + d2*p1*r3) + (c2*p3*t1 + d2*p3*r1) +
1003 (c3*p1*t2 + d3*p1*r2) - (c3*p2*t1 + d3*p2*r1);
1004 double R = ((a1*p2 + c1*f2)*t3 + (b1*p2 + d1*f2)*r3) -
1005 ((a1*p3 + c1*f3)*t2 + (b1*p3 + d1*f3)*r2) -
1006 ((a2*p1 + c2*f1)*t3 + (b2*p1 + d2*f1)*r3) +
1007 ((a2*p3 + c2*f3)*t1 + (b2*p3 + d2*f3)*r1) +
1008 ((a3*p1 + c3*f1)*t2 + (b3*p1 + d3*f1)*r2) -
1009 ((a3*p2 + c3*f2)*t1 + (b3*p2 + d3*f2)*r1);
1010 double S = (a1*e2*r3 + a1*f2*q3) - (a1*e3*r2 + a1*f3*q2) -
1011 (a2*e1*r3 + a2*f1*q3) + (a2*e3*r1 + a2*f3*q1) +
1012 (a3*e1*r2 + a3*f1*q2) - (a3*e2*r1 + a3*f2*q1);
1013 double T = (c1*h2*r3 + c1*p2*q3) - (c1*h3*r2 + c1*p3*q2) -
1014 (c2*h1*r3 + c2*p1*q3) + (c2*h3*r1 + c2*p3*q1) +
1015 (c3*h1*r2 + c3*p1*q2) - (c3*h2*r1 + c3*p2*q1);
1016 double U = ((a1*h2 + c1*e2)*r3 + (a1*p2 + c1*f2)*q3) -
1017 ((a1*h3 + c1*e3)*r2 + (a1*p3 + c1*f3)*q2) -
1018 ((a2*h1 + c2*e1)*r3 + (a2*p1 + c2*f1)*q3) +
1019 ((a2*h3 + c2*e3)*r1 + (a2*p3 + c2*f3)*q1) +
1020 ((a3*h1 + c3*e1)*r2 + (a3*p1 + c3*f1)*q2) -
1021 ((a3*h2 + c3*e2)*r1 + (a3*p2 + c3*f2)*q1);
1022 double V = (b1*e2*t3 + b1*f2*s3) - (b1*e3*t2 + b1*f3*s2) -
1023 (b2*e1*t3 + b2*f1*s3) + (b2*e3*t1 + b2*f3*s1) +
1024 (b3*e1*t2 + b3*f1*s2) - (b3*e2*t1 + b3*f2*s1);
1025 double W = (d1*h2*t3 + d1*p2*s3) - (d1*h3*t2 + d1*p3*s2) -
1026 (d2*h1*t3 + d2*p1*s3) + (d2*h3*t1 + d2*p3*s1) +
1027 (d3*h1*t2 + d3*p1*s2) - (d3*h2*t1 + d3*p2*s1);
1028 double X = ((b1*h2 + d1*e2)*t3 + (b1*p2 + d1*f2)*s3) -
1029 ((b1*h3 + d1*e3)*t2 + (b1*p3 + d1*f3)*s2) -
1030 ((b2*h1 + d2*e1)*t3 + (b2*p1 + d2*f1)*s3) +
1031 ((b2*h3 + d2*e3)*t1 + (b2*p3 + d2*f3)*s1) +
1032 ((b3*h1 + d3*e1)*t2 + (b3*p1 + d3*f1)*s2) -
1033 ((b3*h2 + d3*e2)*t1 + (b3*p2 + d3*f2)*s1);
1034 double Y = (a1*e2*t3 + a1*f2*s3 + b1*e2*r3 + b1*f2*q3) -
1035 (a1*e3*t2 + a1*f3*s2 + b1*e3*r2 + b1*f3*q2) -
1036 (a2*e1*t3 + a2*f1*s3 + b2*e1*r3 + b2*f1*q3) +
1037 (a2*e3*t1 + a2*f3*s1 + b2*e3*r1 + b2*f3*q1) +
1038 (a3*e1*t2 + a3*f1*s2 + b3*e1*r2 + b3*f1*q2) -
1039 (a3*e2*t1 + a3*f2*s1 + b3*e2*r1 + b3*f2*q1);
1040 double Z = (c1*h2*t3 + c1*p2*s3 + d1*h2*r3 + d1*p2*q3) -
1041 (c1*h3*t2 + c1*p3*s2 + d1*h3*r2 + d1*p3*q2) -
1042 (c2*h1*t3 + c2*p1*s3 + d2*h1*r3 + d2*p1*q3) +
1043 (c2*h3*t1 + c2*p3*s1 + d2*h3*r1 + d2*p3*q1) +
1044 (c3*h1*t2 + c3*p1*s2 + d3*h1*r2 + d3*p1*q2) -
1045 (c3*h2*t1 + c3*p2*s1 + d3*h2*r1 + d3*p2*q1);
1046 double AA = ((a1*h2 + c1*e2)*t3 + (a1*p2 + c1*f2)*s3 +
1047 (b1*h2 + d1*e2)*r3 + (b1*p2 + d1*f2)*q3) -
1048 ((a1*h3 + c1*e3)*t2 + (a1*p3 + c1*f3)*s2 +
1049 (b1*h3 + d1*e3)*r2 + (b1*p3 + d1*f3)*q2) -
1050 ((a2*h1 + c2*e1)*t3 + (a2*p1 + c2*f1)*s3 +
1051 (b2*h1 + d2*e1)*r3 + (b2*p1 + d2*f1)*q3) +
1052 ((a2*h3 + c2*e3)*t1 + (a2*p3 + c2*f3)*s1 +
1053 (b2*h3 + d2*e3)*r1 + (b2*p3 + d2*f3)*q1) +
1054 ((a3*h1 + c3*e1)*t2 + (a3*p1 + c3*f1)*s2 +
1055 (b3*h1 + d3*e1)*r2 + (b3*p1 + d3*f1)*q2) -
1056 ((a3*h2 + c3*e2)*t1 + (a3*p2 + c3*f2)*s1 +
1057 (b3*h2 + d3*e2)*r1 + (b3*p2 + d3*f2)*q1);
1059 xvolume = 64.0*(8.0*(A + B + D + E + J + K + M + N) +
1060 4.0*(C + F + G + H + L + O + P + Q + S + T +
1061 V + W) + 2.0*(I + R + U + X + Y + Z) +
1064 //volume[index] = xvolume ;
1065 volume->setIJ(index,1,xvolume) ;
1071 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Volumes on it !"));
1075 if (!onAll) delete [] global_connectivity ;
1081 FIELD<double>* MESH::getArea(const SUPPORT * Support) const throw (MEDEXCEPTION)
1083 const char * LOC = "MESH::getArea(SUPPORT*) : ";
1086 // Support must be on 2D elements
1088 // Make sure that the MESH class is the same as the MESH class attribut
1089 // in the class Support
1090 MESH* myMesh = Support->getMesh();
1092 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
1094 int dim_space = getSpaceDimension();
1095 medEntityMesh support_entity = Support->getEntity();
1096 bool onAll = Support->isOnAllElements();
1098 int nb_type, length_values;
1099 const medGeometryElement* types;
1101 // !!!! WARNING : use of nodal global numbering in the mesh !!!!
1102 const int* global_connectivity;
1106 // nb_type = myMesh->getNumberOfTypes(support_entity);
1107 // length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
1108 // types = getTypes(support_entity);
1112 nb_type = Support->getNumberOfTypes();
1113 length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
1114 types = Support->getTypes();
1118 FIELD<double>* Area;
1120 Area = new FIELD<double>::FIELD(Support,1);
1121 Area->setName("AREA");
1122 Area->setDescription("cells or faces area");
1124 Area->setComponentName(1,"area");
1125 Area->setComponentDescription(1,"desc-comp");
1127 /* string MEDComponentUnit = myMesh->getCoordinatesUnits()[0]+"*"+myMesh->getCoordinatesUnits()[1];*/
1129 string MEDComponentUnit="trucmuch";
1131 Area->setMEDComponentUnit(1,MEDComponentUnit);
1133 Area->setValueType(MED_REEL64);
1135 Area->setIterationNumber(0);
1136 Area->setOrderNumber(0);
1139 double *area = new double[length_values];
1140 //double *area = Area->getValue(MED_FULL_INTERLACE);
1142 const double * coord = getCoordinates(MED_FULL_INTERLACE);
1145 for (int i=0;i<nb_type;i++)
1147 medGeometryElement type = types[i] ;
1152 nb_entity_type = getNumberOfElements(support_entity,type);
1153 global_connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,type);
1157 // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all !"));
1159 nb_entity_type = Support->getNumberOfElements(type);
1161 const int * supp_number = Support->getNumber(type);
1162 const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
1163 const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
1164 int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
1166 for (int k_type = 0; k_type<nb_entity_type; k_type++) {
1167 for (int j_ent = 0; j_ent<(type%100); j_ent++) {
1168 global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
1172 global_connectivity = global_connectivity_tmp ;
1178 case MED_TRIA3 : case MED_TRIA6 :
1180 for (int tria=0;tria<nb_entity_type;tria++)
1182 int tria_index = (type%100)*tria;
1184 int N1 = global_connectivity[tria_index]-1;
1185 int N2 = global_connectivity[tria_index+1]-1;
1186 int N3 = global_connectivity[tria_index+2]-1;
1188 double x1 = coord[dim_space*N1];
1189 double x2 = coord[dim_space*N2];
1190 double x3 = coord[dim_space*N3];
1192 double y1 = coord[dim_space*N1+1];
1193 double y2 = coord[dim_space*N2+1];
1194 double y3 = coord[dim_space*N3+1];
1198 xarea = - ((x2-x1)*(y3-y1) - (x3-x1)*(y2-y1))/2.0;
1202 double z1 = coord[dim_space*N1+2];
1203 double z2 = coord[dim_space*N2+2];
1204 double z3 = coord[dim_space*N3+2];
1206 xarea = sqrt(((y2-y1)*(z3-z1) - (y3-y1)*(z2-z1))*
1207 ((y2-y1)*(z3-z1) - (y3-y1)*(z2-z1)) +
1208 ((x3-x1)*(z2-z1) - (x2-x1)*(z3-z1))*
1209 ((x3-x1)*(z2-z1) - (x2-x1)*(z3-z1)) +
1210 ((x2-x1)*(y3-y1) - (x3-x1)*(y2-y1))*
1211 ((x2-x1)*(y3-y1) - (x3-x1)*(y2-y1)))/2.0;
1214 area[index] = xarea ;
1219 case MED_QUAD4 : case MED_QUAD8 :
1221 for (int quad=0;quad<nb_entity_type;quad++)
1223 int quad_index = (type%100)*quad;
1225 int N1 = global_connectivity[quad_index]-1;
1226 int N2 = global_connectivity[quad_index+1]-1;
1227 int N3 = global_connectivity[quad_index+2]-1;
1228 int N4 = global_connectivity[quad_index+3]-1;
1230 double x1 = coord[dim_space*N1];
1231 double x2 = coord[dim_space*N2];
1232 double x3 = coord[dim_space*N3];
1233 double x4 = coord[dim_space*N4];
1235 double y1 = coord[dim_space*N1+1];
1236 double y2 = coord[dim_space*N2+1];
1237 double y3 = coord[dim_space*N3+1];
1238 double y4 = coord[dim_space*N4+1];
1242 double a1 = (x2-x1)/4.0, a2 = (y2-y1)/4.0;
1243 double b1 = (x3-x4)/4.0, b2 = (y3-y4)/4.0;
1244 double c1 = (x3-x2)/4.0, c2 = (y3-y2)/4.0;
1245 double d1 = (x4-x1)/4.0, d2 = (y4-y1)/4.0;
1247 xarea = - 4.0*(b1*c2 - c1*b2 + a1*c2 - c1*a2 + b1*d2 -
1248 d1*b2 + a1*d2 - d1*a2);
1252 double z1 = coord[dim_space*N1+2];
1253 double z2 = coord[dim_space*N2+2];
1254 double z3 = coord[dim_space*N3+2];
1255 double z4 = coord[dim_space*N4+2];
1257 xarea = (sqrt(((y2-y1)*(z4-z1) - (y4-y1)*(z2-z1))*
1258 ((y2-y1)*(z4-z1) - (y4-y1)*(z2-z1)) +
1259 ((x4-x1)*(z2-z1) - (x2-x1)*(z4-z1))*
1260 ((x4-x1)*(z2-z1) - (x2-x1)*(z4-z1)) +
1261 ((x2-x1)*(y4-y1) - (x4-x1)*(y2-y1))*
1262 ((x2-x1)*(y4-y1) - (x4-x1)*(y2-y1))) +
1263 sqrt(((y4-y3)*(z2-z3) - (y2-y3)*(z4-z3))*
1264 ((y4-y3)*(z2-z3) - (y2-y3)*(z4-z3)) +
1265 ((x2-x3)*(z4-z3) - (x4-x3)*(z2-z3))*
1266 ((x2-x3)*(z4-z3) - (x4-x3)*(z2-z3)) +
1267 ((x4-x3)*(y2-y3) - (x2-x3)*(y4-y3))*
1268 ((x4-x3)*(y2-y3) - (x2-x3)*(y4-y3))))/2.0;
1271 area[index] = xarea ;
1277 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Areas on it !"));
1281 if (!onAll) delete [] global_connectivity ;
1284 Area->setValue(MED_FULL_INTERLACE,area);
1289 FIELD<double>* MESH::getLength(const SUPPORT * Support) const throw (MEDEXCEPTION)
1291 const char * LOC = "MESH::getLength(SUPPORT*) : ";
1294 // Support must be on 1D elements
1296 // Make sure that the MESH class is the same as the MESH class attribut
1297 // in the class Support
1298 MESH* myMesh = Support->getMesh();
1300 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
1302 int dim_space = getSpaceDimension();
1303 medEntityMesh support_entity = Support->getEntity();
1304 bool onAll = Support->isOnAllElements();
1306 int nb_type, length_values;
1307 const medGeometryElement* types;
1309 // !!!! WARNING : use of nodal global numbering in the mesh !!!!
1310 const int* global_connectivity;
1314 // nb_type = myMesh->getNumberOfTypes(support_entity);
1315 // length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
1316 // types = getTypes(support_entity);
1320 nb_type = Support->getNumberOfTypes();
1321 length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
1322 types = Support->getTypes();
1326 FIELD<double>* Length;
1328 Length = new FIELD<double>::FIELD(Support,1);
1329 // double *length = new double [length_values];
1330 Length->setValueType(MED_REEL64);
1332 //const double *length = Length->getValue(MED_FULL_INTERLACE);
1333 MEDARRAY<double> * length = Length->getvalue();
1335 const double * coord = getCoordinates(MED_FULL_INTERLACE);
1338 for (int i=0;i<nb_type;i++)
1340 medGeometryElement type = types[i] ;
1345 nb_entity_type = getNumberOfElements(support_entity,type);
1346 global_connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,type);
1350 //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all !"));
1352 nb_entity_type = Support->getNumberOfElements(type);
1354 const int * supp_number = Support->getNumber(type);
1355 const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
1356 const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
1357 int* global_connectivity_tmp = new int[(type%100)*nb_entity_type];
1359 for (int k_type = 0; k_type<nb_entity_type; k_type++) {
1360 for (int j_ent = 0; j_ent<(type%100); j_ent++) {
1361 global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
1364 global_connectivity = global_connectivity_tmp ;
1369 case MED_SEG2 : case MED_SEG3 :
1371 for (int edge=0;edge<nb_entity_type;edge++)
1373 int edge_index = (type%100)*edge;
1375 int N1 = global_connectivity[edge_index]-1;
1376 int N2 = global_connectivity[edge_index+1]-1;
1378 double x1 = coord[dim_space*N1];
1379 double x2 = coord[dim_space*N2];
1381 double y1 = coord[dim_space*N1+1];
1382 double y2 = coord[dim_space*N2+1];
1384 double z1, z2 ; z1 = 0.0 ; z2 = 0.0 ; if (dim_space==3) {
1385 z1 = coord[dim_space*N1+2]; z2 = coord[dim_space*N2+2];}
1387 xlength = sqrt((x1 - x2)*(x1 - x2) + (y1 - y2)*(y1 - y2) +
1388 (z1 - z2)*(z1 - z2));
1390 length->setIJ(index,1,xlength) ;
1396 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Lengths on it !"));
1400 if (!onAll) delete [] global_connectivity ;
1406 FIELD<double>* MESH::getNormal(const SUPPORT * Support) const throw (MEDEXCEPTION)
1408 const char * LOC = "MESH::getNormal(SUPPORT*) : ";
1411 // Support must be on 2D or 1D elements
1413 // Make sure that the MESH class is the same as the MESH class attribut
1414 // in the class Support
1415 MESH* myMesh = Support->getMesh();
1417 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh : pointeur problem !"));
1419 int dim_space = getSpaceDimension();
1420 medEntityMesh support_entity = Support->getEntity();
1421 bool onAll = Support->isOnAllElements();
1423 int nb_type, length_values;
1424 const medGeometryElement* types;
1426 // !!!! WARNING : use of nodal global numbering in the mesh !!!!
1427 const int* global_connectivity;
1431 // nb_type = myMesh->getNumberOfTypes(support_entity);
1432 // length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
1433 // types = getTypes(support_entity);
1437 nb_type = Support->getNumberOfTypes();
1438 length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
1439 types = Support->getTypes();
1444 FIELD<double>* Normal = new FIELD<double>::FIELD(Support,dim_space);
1445 Normal->setName("NORMAL");
1446 Normal->setDescription("cells or faces normal");
1447 for (int k=1;k<=dim_space;k++) {
1448 Normal->setComponentName(k,"normal");
1449 Normal->setComponentDescription(k,"desc-comp");
1450 Normal->setMEDComponentUnit(k,"unit");
1453 Normal->setValueType(MED_REEL64);
1455 Normal->setIterationNumber(MED_NOPDT);
1456 Normal->setOrderNumber(MED_NONOR);
1457 Normal->setTime(0.0);
1459 double * normal = new double [dim_space*length_values];
1460 //double *normal = Normal->getValue(MED_FULL_INTERLACE);
1462 const double * coord = getCoordinates(MED_FULL_INTERLACE);
1465 for (int i=0;i<nb_type;i++)
1467 medGeometryElement type = types[i] ;
1469 // Make sure we trying to get Normals on
1470 // TRIA3 or TRIA6 or QUAD4 QUAD8 (2D) cells on a 3D mesh
1471 // or on SEG2 or SEG3 (1D) cells on a 2D mesh
1473 if ( (((type==MED_TRIA3) || (type==MED_TRIA6) ||
1474 (type==MED_QUAD4) || (type==MED_QUAD8)) &&
1475 (dim_space != 3)) || (((type==MED_SEG2) || (type==MED_SEG3)) &&
1477 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh : dimension problem !"));
1479 double xnormal1, xnormal2, xnormal3 ;
1483 nb_entity_type = getNumberOfElements(support_entity,type);
1484 global_connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,type);
1488 // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all for instance !"));
1489 nb_entity_type = Support->getNumberOfElements(type);
1491 const int * supp_number = Support->getNumber(type);
1492 const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
1493 const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
1494 int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
1496 for (int k_type = 0; k_type<nb_entity_type; k_type++) {
1497 for (int j_ent = 0; j_ent<(type%100); j_ent++) {
1498 global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
1502 global_connectivity = global_connectivity_tmp ;
1508 case MED_TRIA3 : case MED_TRIA6 :
1510 for (int tria=0;tria<nb_entity_type;tria++)
1512 int tria_index = (type%100)*tria;
1514 int N1 = global_connectivity[tria_index]-1;
1515 int N2 = global_connectivity[tria_index+1]-1;
1516 int N3 = global_connectivity[tria_index+2]-1;
1518 //double xarea; !! UNUSED VARIABLE !!
1519 double x1 = coord[dim_space*N1];
1520 double x2 = coord[dim_space*N2];
1521 double x3 = coord[dim_space*N3];
1523 double y1 = coord[dim_space*N1+1];
1524 double y2 = coord[dim_space*N2+1];
1525 double y3 = coord[dim_space*N3+1];
1527 double z1 = coord[dim_space*N1+2];
1528 double z2 = coord[dim_space*N2+2];
1529 double z3 = coord[dim_space*N3+2];
1531 xnormal1 = ((y2-y1)*(z3-z1) - (y3-y1)*(z2-z1))/2.0;
1532 xnormal2 = ((x3-x1)*(z2-z1) - (x2-x1)*(z3-z1))/2.0;
1533 xnormal3 = ((x2-x1)*(y3-y1) - (x3-x1)*(y2-y1))/2.0;
1535 normal[3*index] = xnormal1 ;
1536 normal[3*index+1] = xnormal2 ;
1537 normal[3*index+2] = xnormal3 ;
1543 case MED_QUAD4 : case MED_QUAD8 :
1545 for (int quad=0;quad<nb_entity_type;quad++)
1547 int quad_index = (type%100)*quad;
1549 int N1 = global_connectivity[quad_index]-1;
1550 int N2 = global_connectivity[quad_index+1]-1;
1551 int N3 = global_connectivity[quad_index+2]-1;
1552 int N4 = global_connectivity[quad_index+3]-1;
1555 double x1 = coord[dim_space*N1];
1556 double x2 = coord[dim_space*N2];
1557 double x3 = coord[dim_space*N3];
1558 double x4 = coord[dim_space*N4];
1560 double y1 = coord[dim_space*N1+1];
1561 double y2 = coord[dim_space*N2+1];
1562 double y3 = coord[dim_space*N3+1];
1563 double y4 = coord[dim_space*N4+1];
1565 double z1 = coord[dim_space*N1+2];
1566 double z2 = coord[dim_space*N2+2];
1567 double z3 = coord[dim_space*N3+2];
1568 double z4 = coord[dim_space*N4+2];
1570 xnormal1 = (y2-y1)*(z4-z1) - (y4-y1)*(z2-z1);
1571 xnormal2 = (x4-x1)*(z2-z1) - (x2-x1)*(z4-z1);
1572 xnormal3 = (x2-x1)*(y4-y1) - (x4-x1)*(y2-y1);
1573 xarea = sqrt(xnormal1*xnormal1 + xnormal2*xnormal2 +
1576 xnormal1 = xnormal1/xarea;
1577 xnormal2 = xnormal2/xarea;
1578 xnormal3 = xnormal3/xarea;
1580 xarea = (sqrt(((y2-y1)*(z4-z1) - (y4-y1)*(z2-z1))*
1581 ((y2-y1)*(z4-z1) - (y4-y1)*(z2-z1)) +
1582 ((x4-x1)*(z2-z1) - (x2-x1)*(z4-z1))*
1583 ((x4-x1)*(z2-z1) - (x2-x1)*(z4-z1)) +
1584 ((x2-x1)*(y4-y1) - (x4-x1)*(y2-y1))*
1585 ((x2-x1)*(y4-y1) - (x4-x1)*(y2-y1))) +
1586 sqrt(((y3-y2)*(z4-z3) - (y4-y3)*(z3-z2))*
1587 ((y3-y2)*(z4-z3) - (y4-y3)*(z3-z2)) +
1588 ((x4-x3)*(z2-z2) - (x3-x2)*(z4-z3))*
1589 ((x4-x3)*(z4-z2) - (x3-x2)*(z4-z3)) +
1590 ((x3-x2)*(y4-y3) - (x4-x3)*(y3-y2))*
1591 ((x3-x2)*(y4-y3) - (x4-x3)*(y3-y2))))/2.0;
1593 xnormal1 = xnormal1*xarea;
1594 xnormal2 = xnormal2*xarea;
1595 xnormal3 = xnormal3*xarea;
1597 normal[3*index] = xnormal1 ;
1598 normal[3*index+1] = xnormal2 ;
1599 normal[3*index+2] = xnormal3 ;
1605 case MED_SEG2 : case MED_SEG3 :
1607 for (int edge=0;edge<nb_entity_type;edge++)
1609 int edge_index = (type%100)*edge;
1611 int N1 = global_connectivity[edge_index]-1;
1612 int N2 = global_connectivity[edge_index+1]-1;
1614 double x1 = coord[dim_space*N1];
1615 double x2 = coord[dim_space*N2];
1617 double y1 = coord[dim_space*N1+1];
1618 double y2 = coord[dim_space*N2+1];
1620 xnormal1 = -(y2-y1);
1623 normal[2*index] = xnormal1 ;
1624 normal[2*index+1] = xnormal2 ;
1631 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get Normals on it !"));
1635 if (!onAll) delete [] global_connectivity ;
1638 Normal->setValue(MED_FULL_INTERLACE,normal);
1646 FIELD<double>* MESH::getBarycenter(const SUPPORT * Support) const throw (MEDEXCEPTION)
1648 const char * LOC = "MESH::getBarycenter(SUPPORT*) : ";
1651 // Make sure that the MESH class is the same as the MESH class attribut
1652 // in the class Support
1653 MESH* myMesh = Support->getMesh();
1655 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"no compatibility between *this and SUPPORT::_mesh !"));
1657 int dim_space = getSpaceDimension();
1658 medEntityMesh support_entity = Support->getEntity();
1659 bool onAll = Support->isOnAllElements();
1661 int nb_type, length_values;
1662 const medGeometryElement* types;
1664 // !!!! WARNING : use of nodal global numbering in the mesh !!!!
1665 const int* global_connectivity;
1669 // nb_type = myMesh->getNumberOfTypes(support_entity);
1670 // length_values = getNumberOfElements(support_entity,MED_ALL_ELEMENTS);
1671 // types = getTypes(support_entity);
1675 nb_type = Support->getNumberOfTypes();
1676 length_values = Support->getNumberOfElements(MED_ALL_ELEMENTS);
1677 types = Support->getTypes();
1681 FIELD<double>* Barycenter;
1683 Barycenter = new FIELD<double>::FIELD(Support,dim_space);
1684 Barycenter->setName("BARYCENTER");
1685 Barycenter->setDescription("cells or faces barycenter");
1687 for (int k=0;k<dim_space;k++) {
1689 Barycenter->setComponentName(kp1,myMesh->getCoordinatesNames()[k]);
1690 Barycenter->setComponentDescription(kp1,"desc-comp");
1691 Barycenter->setMEDComponentUnit(kp1,myMesh->getCoordinatesUnits()[k]);
1694 Barycenter->setValueType(MED_REEL64);
1696 Barycenter->setIterationNumber(0);
1697 Barycenter->setOrderNumber(0);
1698 Barycenter->setTime(0.0);
1700 double *barycenter = new double [dim_space*length_values];
1701 // double *barycenter = Barycenter->getValue(MED_FULL_INTERLACE);
1703 const double * coord = getCoordinates(MED_FULL_INTERLACE);
1706 for (int i=0;i<nb_type;i++)
1708 medGeometryElement type = types[i] ;
1709 double xbarycenter1, xbarycenter2, xbarycenter3;
1713 nb_entity_type = getNumberOfElements(support_entity,type);
1714 global_connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,type);
1718 // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Support must be on all !"));
1719 nb_entity_type = Support->getNumberOfElements(type);
1721 const int * supp_number = Support->getNumber(type);
1722 const int * connectivity = getConnectivity(MED_FULL_INTERLACE,MED_NODAL,support_entity,MED_ALL_ELEMENTS);
1723 const int * connectivityIndex = getConnectivityIndex(MED_NODAL,support_entity);
1724 int * global_connectivity_tmp = new int[(type%100)*nb_entity_type];
1726 for (int k_type = 0; k_type<nb_entity_type; k_type++) {
1727 for (int j_ent = 0; j_ent<(type%100); j_ent++) {
1728 global_connectivity_tmp[k_type*(type%100)+j_ent] = connectivity[connectivityIndex[supp_number[k_type]-1]+j_ent-1];
1731 global_connectivity = global_connectivity_tmp;
1736 case MED_TETRA4 : case MED_TETRA10 :
1738 for (int tetra =0;tetra<nb_entity_type;tetra++)
1740 int tetra_index = (type%100)*tetra;
1742 int N1 = global_connectivity[tetra_index]-1;
1743 int N2 = global_connectivity[tetra_index+1]-1;
1744 int N3 = global_connectivity[tetra_index+2]-1;
1745 int N4 = global_connectivity[tetra_index+3]-1;
1747 double x1 = coord[dim_space*N1];
1748 double x2 = coord[dim_space*N2];
1749 double x3 = coord[dim_space*N3];
1750 double x4 = coord[dim_space*N4];
1752 double y1 = coord[dim_space*N1+1];
1753 double y2 = coord[dim_space*N2+1];
1754 double y3 = coord[dim_space*N3+1];
1755 double y4 = coord[dim_space*N4+1];
1757 double z1 = coord[dim_space*N1+2];
1758 double z2 = coord[dim_space*N2+2];
1759 double z3 = coord[dim_space*N3+2];
1760 double z4 = coord[dim_space*N4+2];
1762 xbarycenter1 = (x1 + x2 + x3 + x4)/4.0;
1763 xbarycenter2 = (y1 + y2 + y3 + y4)/4.0;
1764 xbarycenter3 = (z1 + z2 + z3 + z4)/4.0;
1765 barycenter[3*index] = xbarycenter1 ;
1766 barycenter[3*index+1] = xbarycenter2 ;
1767 barycenter[3*index+2] = xbarycenter3 ;
1772 case MED_PYRA5 : case MED_PYRA13 :
1774 for (int pyra=0;pyra<nb_entity_type;pyra++)
1776 int pyra_index = (type%100)*pyra;
1778 int N1 = global_connectivity[pyra_index]-1;
1779 int N2 = global_connectivity[pyra_index+1]-1;
1780 int N3 = global_connectivity[pyra_index+2]-1;
1781 int N4 = global_connectivity[pyra_index+3]-1;
1782 int N5 = global_connectivity[pyra_index+4]-1;
1784 double x1 = coord[dim_space*N1];
1785 double x2 = coord[dim_space*N2];
1786 double x3 = coord[dim_space*N3];
1787 double x4 = coord[dim_space*N4];
1788 double x5 = coord[dim_space*N5];
1790 double y1 = coord[dim_space*N1+1];
1791 double y2 = coord[dim_space*N2+1];
1792 double y3 = coord[dim_space*N3+1];
1793 double y4 = coord[dim_space*N4+1];
1794 double y5 = coord[dim_space*N5+1];
1796 double z1 = coord[dim_space*N1+2];
1797 double z2 = coord[dim_space*N2+2];
1798 double z3 = coord[dim_space*N3+2];
1799 double z4 = coord[dim_space*N4+2];
1800 double z5 = coord[dim_space*N5+2];
1802 xbarycenter1 = (x1 + x2 + x3 + x4 + x5)/5.0;
1803 xbarycenter2 = (y1 + y2 + y3 + y4 + y5)/5.0;
1804 xbarycenter3 = (z1 + z2 + z3 + z4 + z5)/5.0;
1805 barycenter[3*index] = xbarycenter1 ;
1806 barycenter[3*index+1] = xbarycenter2 ;
1807 barycenter[3*index+2] = xbarycenter3 ;
1812 case MED_PENTA6 : case MED_PENTA15 :
1814 for (int penta=0;penta<nb_entity_type;penta++)
1816 int penta_index = (type%100)*penta;
1818 int N1 = global_connectivity[penta_index]-1;
1819 int N2 = global_connectivity[penta_index+1]-1;
1820 int N3 = global_connectivity[penta_index+2]-1;
1821 int N4 = global_connectivity[penta_index+3]-1;
1822 int N5 = global_connectivity[penta_index+4]-1;
1823 int N6 = global_connectivity[penta_index+5]-1;
1825 double x1 = coord[dim_space*N1];
1826 double x2 = coord[dim_space*N2];
1827 double x3 = coord[dim_space*N3];
1828 double x4 = coord[dim_space*N4];
1829 double x5 = coord[dim_space*N5];
1830 double x6 = coord[dim_space*N6];
1832 double y1 = coord[dim_space*N1+1];
1833 double y2 = coord[dim_space*N2+1];
1834 double y3 = coord[dim_space*N3+1];
1835 double y4 = coord[dim_space*N4+1];
1836 double y5 = coord[dim_space*N5+1];
1837 double y6 = coord[dim_space*N6+1];
1839 double z1 = coord[dim_space*N1+2];
1840 double z2 = coord[dim_space*N2+2];
1841 double z3 = coord[dim_space*N3+2];
1842 double z4 = coord[dim_space*N4+2];
1843 double z5 = coord[dim_space*N5+2];
1844 double z6 = coord[dim_space*N6+2];
1846 xbarycenter1 = (x1 + x2 + x3 + x4 + x5 + x6)/6.0;
1847 xbarycenter2 = (y1 + y2 + y3 + y4 + y5 + y6)/6.0;
1848 xbarycenter3 = (z1 + z2 + z3 + z4 + z5 + z6)/6.0;
1849 barycenter[3*index] = xbarycenter1 ;
1850 barycenter[3*index+1] = xbarycenter2 ;
1851 barycenter[3*index+2] = xbarycenter3 ;
1856 case MED_HEXA8 : case MED_HEXA20 :
1858 for (int hexa=0;hexa<nb_entity_type;hexa++)
1860 int hexa_index = (type%100)*hexa;
1862 int N1 = global_connectivity[hexa_index]-1;
1863 int N2 = global_connectivity[hexa_index+1]-1;
1864 int N3 = global_connectivity[hexa_index+2]-1;
1865 int N4 = global_connectivity[hexa_index+3]-1;
1866 int N5 = global_connectivity[hexa_index+4]-1;
1867 int N6 = global_connectivity[hexa_index+5]-1;
1868 int N7 = global_connectivity[hexa_index+6]-1;
1869 int N8 = global_connectivity[hexa_index+7]-1;
1871 double x1 = coord[dim_space*N1];
1872 double x2 = coord[dim_space*N2];
1873 double x3 = coord[dim_space*N3];
1874 double x4 = coord[dim_space*N4];
1875 double x5 = coord[dim_space*N5];
1876 double x6 = coord[dim_space*N6];
1877 double x7 = coord[dim_space*N7];
1878 double x8 = coord[dim_space*N8];
1880 double y1 = coord[dim_space*N1+1];
1881 double y2 = coord[dim_space*N2+1];
1882 double y3 = coord[dim_space*N3+1];
1883 double y4 = coord[dim_space*N4+1];
1884 double y5 = coord[dim_space*N5+1];
1885 double y6 = coord[dim_space*N6+1];
1886 double y7 = coord[dim_space*N7+1];
1887 double y8 = coord[dim_space*N8+1];
1889 double z1 = coord[dim_space*N1+2];
1890 double z2 = coord[dim_space*N2+2];
1891 double z3 = coord[dim_space*N3+2];
1892 double z4 = coord[dim_space*N4+2];
1893 double z5 = coord[dim_space*N5+2];
1894 double z6 = coord[dim_space*N6+2];
1895 double z7 = coord[dim_space*N7+2];
1896 double z8 = coord[dim_space*N8+2];
1898 xbarycenter1 = (x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8)/8.0;
1899 xbarycenter2 = (y1 + y2 + y3 + y4 + y5 + y6 + y7 + y8)/8.0;
1900 xbarycenter3 = (z1 + z2 + z3 + z4 + z5 + z6 + z7 + z8)/8.0;
1902 barycenter[3*index] = xbarycenter1 ;
1903 barycenter[3*index+1] = xbarycenter2 ;
1904 barycenter[3*index+2] = xbarycenter3 ;
1910 case MED_TRIA3 : case MED_TRIA6 :
1912 for (int tria=0;tria<nb_entity_type;tria++)
1914 int tria_index = (type%100)*tria;
1916 int N1 = global_connectivity[tria_index]-1;
1917 int N2 = global_connectivity[tria_index+1]-1;
1918 int N3 = global_connectivity[tria_index+2]-1;
1920 double x1 = coord[dim_space*N1];
1921 double x2 = coord[dim_space*N2];
1922 double x3 = coord[dim_space*N3];
1924 double y1 = coord[dim_space*N1+1];
1925 double y2 = coord[dim_space*N2+1];
1926 double y3 = coord[dim_space*N3+1];
1928 xbarycenter1 = (x1 + x2 + x3)/3.0;
1929 xbarycenter2 = (y1 + y2 + y3)/3.0;
1933 barycenter[2*index] = xbarycenter1 ;
1934 barycenter[2*index+1] = xbarycenter2 ;
1939 coord[dim_space*N1+2];
1941 coord[dim_space*N2+2];
1943 coord[dim_space*N3+2];
1945 xbarycenter3 = (z1 + z2 + z3)/3.0;
1947 barycenter[3*index] = xbarycenter1 ;
1948 barycenter[3*index+1] = xbarycenter2 ;
1949 barycenter[3*index+2] = xbarycenter3 ;
1956 case MED_QUAD4 : case MED_QUAD8 :
1958 for (int quad=0;quad<nb_entity_type;quad++)
1960 int quad_index = (type%100)*quad;
1962 int N1 = global_connectivity[quad_index]-1;
1963 int N2 = global_connectivity[quad_index+1]-1;
1964 int N3 = global_connectivity[quad_index+2]-1;
1965 int N4 = global_connectivity[quad_index+3]-1;
1967 double x1 = coord[dim_space*N1];
1968 double x2 = coord[dim_space*N2];
1969 double x3 = coord[dim_space*N3];
1970 double x4 = coord[dim_space*N4];
1972 double y1 = coord[dim_space*N1+1];
1973 double y2 = coord[dim_space*N2+1];
1974 double y3 = coord[dim_space*N3+1];
1975 double y4 = coord[dim_space*N4+1];
1977 xbarycenter1 = (x1 + x2 + x3 + x4)/4.0;
1978 xbarycenter2 = (y1 + y2 + y3 + y4)/4.0;
1982 barycenter[2*index] = xbarycenter1 ;
1983 barycenter[2*index+1] = xbarycenter2 ;
1987 double z1 = coord[dim_space*N1+2];
1988 double z2 = coord[dim_space*N2+2];
1989 double z3 = coord[dim_space*N3+2];
1990 double z4 = coord[dim_space*N4+2];
1992 xbarycenter3 = (z1 + z2 + z3 + z4)/4.0;
1994 barycenter[3*index] = xbarycenter1 ;
1995 barycenter[3*index+1] = xbarycenter2 ;
1996 barycenter[3*index+2] = xbarycenter3 ;
2002 case MED_SEG2 : case MED_SEG3 :
2004 for (int edge=0;edge<nb_entity_type;edge++)
2006 int edge_index = (type%100)*edge;
2008 int N1 = global_connectivity[edge_index]-1;
2009 int N2 = global_connectivity[edge_index+1]-1;
2011 double x1 = coord[dim_space*N1];
2012 double x2 = coord[dim_space*N2];
2014 double y1 = coord[dim_space*N1+1];
2015 double y2 = coord[dim_space*N2+1];
2017 xbarycenter1 = (x1 + x2)/2.0;
2018 xbarycenter2 = (y1 + y2)/2.0;
2022 barycenter[2*index] = xbarycenter1 ;
2023 barycenter[2*index+1] = xbarycenter2 ;
2027 double z1 = coord[dim_space*N1+2];
2028 double z2 = coord[dim_space*N2+2];
2030 xbarycenter3 = (z1 + z2)/2.0;
2032 barycenter[3*index] = xbarycenter1 ;
2033 barycenter[3*index+1] = xbarycenter2 ;
2034 barycenter[3*index+2] = xbarycenter3 ;
2041 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Bad Support to get a barycenter on it (in fact unknown type) !"));
2045 if (!onAll) delete [] global_connectivity ;
2048 Barycenter->setValue(MED_FULL_INTERLACE,barycenter);
2050 delete[] barycenter ;
2057 //=======================================================================
2058 //function : checkGridFillCoords
2059 //purpose : if this->_isAGrid, assure that _coordinate is filled
2060 //=======================================================================
2062 // inline void MESH::checkGridFillCoords() const
2065 // ((GRID *) this)->fillCoordinates();
2068 //=======================================================================
2069 //function : checkGridFillConnectivity
2070 //purpose : if this->_isAGrid, assure that _connectivity is filled
2071 //=======================================================================
2073 // inline void MESH::checkGridFillConnectivity() const
2076 // ((GRID *) this)->fillConnectivity();
2079 bool MESH::isEmpty() const
2081 bool notempty = _name != "" || _coordinate != NULL || _connectivity != NULL ||
2082 _spaceDimension !=MED_INVALID || _meshDimension !=MED_INVALID ||
2083 _numberOfNodes !=MED_INVALID || _numberOfNodesFamilies !=0 ||
2084 _familyNode.size() != 0 || _numberOfCellsFamilies != 0 ||
2085 _familyCell.size() != 0 || _numberOfFacesFamilies != 0 ||
2086 _familyFace.size() != 0 || _numberOfEdgesFamilies !=0 ||
2087 _familyEdge.size() != 0 || _isAGrid != 0 ;
2091 void MESH::read(int index)
2093 const char * LOC = "MESH::read(int index=0) : ";
2096 if (_drivers[index]) {
2097 _drivers[index]->open();
2098 _drivers[index]->read();
2099 _drivers[index]->close();
2102 throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
2103 << "The index given is invalid, index must be between 0 and |"
2108 // ((GRID *) this)->fillMeshAfterRead();
2112 //=======================================================================
2113 //function : getSkin
2115 //=======================================================================
2117 SUPPORT * MESH::getSkin(const SUPPORT * Support3D) throw (MEDEXCEPTION)
2119 const char * LOC = "MESH::getSkin : " ;
2122 if (this != Support3D->getMesh())
2123 throw MEDEXCEPTION(STRING(LOC) << "no compatibility between *this and SUPPORT::_mesh !");
2124 if (_meshDimension != 3 || Support3D->getEntity() != MED_CELL)
2125 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Defined on 3D cells only"));
2127 // well, all rigth !
2128 SUPPORT * mySupport = new SUPPORT(this,"Skin",MED_FACE);
2129 mySupport->setAll(false);
2131 list<int> myElementsList ;
2134 calculateConnectivity(MED_FULL_INTERLACE, MED_DESCENDING, MED_CELL);
2135 if (Support3D->isOnAllElements())
2137 int * myConnectivityValue = const_cast <int*> (getReverseConnectivity(MED_DESCENDING)) ;
2138 int nbFaces = getNumberOfElements(MED_FACE,MED_ALL_ELEMENTS);
2139 for (i=0, j=1 ; j<=nbFaces; ++j, i += 2)
2141 int cellNb1 = myConnectivityValue [i];
2142 int cellNb2 = myConnectivityValue [i+1];
2143 //MESSAGE( " FACE # " << j << " -- Cells: " << cellNb1 << ", " << cellNb2 );
2144 if ((cellNb1 == 0 || cellNb2 == 0) && (cellNb1 + cellNb2 > 0))
2146 myElementsList.push_back( j ) ;
2153 map<int,int> FaceNbEncounterNb;
2155 int * myConnectivityValue = const_cast <int *>
2156 (getConnectivity(MED_FULL_INTERLACE, MED_DESCENDING,
2157 MED_CELL, MED_ALL_ELEMENTS));
2158 int * myConnectivityIndex = const_cast <int *> (getConnectivityIndex(MED_DESCENDING, MED_CELL));
2159 int * myCellNbs = const_cast <int *> (Support3D->getnumber()->getValue());
2160 int nbCells = Support3D->getnumber()->getLength();
2161 for (i=0; i<nbCells; ++i)
2163 int cellNb = myCellNbs [ i ];
2164 int faceFirst = myConnectivityIndex[ cellNb-1 ];
2165 int faceLast = myConnectivityIndex[ cellNb ];
2166 for (j = faceFirst; j < faceLast; ++j)
2168 int faceNb = abs( myConnectivityValue [ j-1 ] );
2169 //MESSAGE( "Cell # " << i << " -- Face: " << faceNb);
2170 if (FaceNbEncounterNb.find( faceNb ) == FaceNbEncounterNb.end())
2171 FaceNbEncounterNb[ faceNb ] = 1;
2173 FaceNbEncounterNb[ faceNb ] ++;
2176 map<int,int>::iterator FaceNbEncounterNbItr;
2177 for (FaceNbEncounterNbItr = FaceNbEncounterNb.begin();
2178 FaceNbEncounterNbItr != FaceNbEncounterNb.end();
2179 FaceNbEncounterNbItr ++)
2180 if ((*FaceNbEncounterNbItr).second == 1)
2182 myElementsList.push_back( (*FaceNbEncounterNbItr).first) ;
2186 // Well, we must know how many geometric type we have found
2187 int * myListArray = new int[size] ;
2189 list<int>::iterator myElementsListIt ;
2190 for (myElementsListIt=myElementsList.begin();myElementsListIt!=myElementsList.end();myElementsListIt++) {
2191 myListArray[id]=(*myElementsListIt) ;
2195 int numberOfGeometricType ;
2196 medGeometryElement* geometricType ;
2197 int * numberOfGaussPoint ;
2198 int * geometricTypeNumber ;
2199 int * numberOfEntities ;
2200 // MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY() ;
2201 int * mySkyLineArrayIndex ;
2203 int numberOfType = getNumberOfTypes(MED_FACE) ;
2204 if (numberOfType == 1) { // wonderfull : it's easy !
2205 numberOfGeometricType = 1 ;
2206 geometricType = new medGeometryElement[1] ;
2207 const medGeometryElement * allType = getTypes(MED_FACE);
2208 geometricType[0] = allType[0] ;
2209 numberOfGaussPoint = new int[1] ;
2210 numberOfGaussPoint[0] = 1 ;
2211 geometricTypeNumber = new int[1] ; // not use, but initialized to nothing
2212 geometricTypeNumber[0] = 0 ;
2213 numberOfEntities = new int[1] ;
2214 numberOfEntities[0] = size ;
2215 mySkyLineArrayIndex = new int[2] ;
2216 mySkyLineArrayIndex[0]=1 ;
2217 mySkyLineArrayIndex[1]=1+size ;
2220 map<medGeometryElement,int> theType ;
2221 for (myElementsListIt=myElementsList.begin();myElementsListIt!=myElementsList.end();myElementsListIt++) {
2222 medGeometryElement myType = getElementType(MED_FACE,*myElementsListIt) ;
2223 if (theType.find(myType) != theType.end() )
2224 theType[myType]+=1 ;
2228 numberOfGeometricType = theType.size() ;
2229 geometricType = new medGeometryElement[numberOfGeometricType] ;
2230 //const medGeometryElement * allType = getTypes(MED_FACE); !! UNUSED VARIABLE !!
2231 numberOfGaussPoint = new int[numberOfGeometricType] ;
2232 geometricTypeNumber = new int[numberOfGeometricType] ; // not use, but initialized to nothing
2233 numberOfEntities = new int[numberOfGeometricType] ;
2234 mySkyLineArrayIndex = new int[numberOfGeometricType+1] ;
2236 mySkyLineArrayIndex[0]=1 ;
2237 map<medGeometryElement,int>::iterator theTypeIt ;
2238 for (theTypeIt=theType.begin();theTypeIt!=theType.end();theTypeIt++) {
2239 geometricType[index] = (*theTypeIt).first ;
2240 numberOfGaussPoint[index] = 1 ;
2241 geometricTypeNumber[index] = 0 ;
2242 numberOfEntities[index] = (*theTypeIt).second ;
2243 mySkyLineArrayIndex[index+1]=mySkyLineArrayIndex[index]+numberOfEntities[index] ;
2247 // mySkyLineArray->setMEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
2248 MEDSKYLINEARRAY * mySkyLineArray = new MEDSKYLINEARRAY(numberOfGeometricType,size,mySkyLineArrayIndex,myListArray) ;
2250 mySupport->setNumberOfGeometricType(numberOfGeometricType) ;
2251 mySupport->setGeometricType(geometricType) ;
2252 mySupport->setNumberOfGaussPoint(numberOfGaussPoint) ;
2253 // mySupport->setGeometricTypeNumber(geometricTypeNumber) ;
2254 mySupport->setNumberOfElements(numberOfEntities) ;
2255 mySupport->setTotalNumberOfElements(size) ;
2256 mySupport->setNumber(mySkyLineArray) ;
2258 delete[] numberOfEntities;
2259 delete[] geometricTypeNumber;
2260 delete[] numberOfGaussPoint;
2261 delete[] geometricType;
2262 delete[] mySkyLineArrayIndex;
2263 delete[] myListArray;
2264 // delete mySkyLineArray;
2272 return a SUPPORT pointer on the union of all SUPPORTs in Supports.
2273 You should delete this pointer after use to avoid memory leaks.
2275 SUPPORT * MESH::mergeSupports(const vector<SUPPORT *> Supports) const throw (MEDEXCEPTION)
2277 const char * LOC = "MESH:::mergeSupports(const vector<SUPPORT *> ) : " ;
2280 SUPPORT * returnedSupport;
2281 string returnedSupportName;
2282 string returnedSupportDescription;
2283 char * returnedSupportNameChar;
2284 char * returnedSupportDescriptionChar;
2285 int size = Supports.size();
2289 MESSAGE(LOC <<" there is only one SUPPORT in the argument list, the method return a copy of this object !");
2290 SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
2292 returnedSupport = new SUPPORT(*obj);
2294 int lenName = strlen((Supports[0]->getName()).c_str()) + 8 + 1;
2295 int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 8 + 1;
2297 returnedSupportNameChar = new char[lenName];
2298 returnedSupportDescriptionChar = new char[lenDescription];
2300 returnedSupportNameChar = strcpy(returnedSupportNameChar,"Copy of ");
2301 returnedSupportNameChar = strcat(returnedSupportNameChar,(Supports[0]->getName()).c_str());
2302 returnedSupportDescriptionChar = strcpy(returnedSupportDescriptionChar,"Copy of ");
2303 returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
2304 (Supports[0]->getDescription()).c_str());
2306 returnedSupportName = string(returnedSupportNameChar);
2307 returnedSupportDescription = string(returnedSupportDescriptionChar);
2309 returnedSupport->setName(returnedSupportName);
2310 returnedSupport->setDescription(returnedSupportDescription);
2312 delete [] returnedSupportNameChar;
2313 delete [] returnedSupportDescriptionChar;
2317 SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
2318 returnedSupport = new SUPPORT(*obj);
2320 int lenName = strlen((Supports[0]->getName()).c_str()) + 9 + 1;
2321 int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 9 + 1;
2323 for(int i = 1;i<size;i++)
2325 obj = const_cast <SUPPORT *> (Supports[i]);
2326 returnedSupport->blending(obj);
2330 lenName = lenName + 5 + strlen((Supports[i]->getName()).c_str());
2331 lenDescription = lenDescription + 5 +
2332 strlen((Supports[i]->getDescription()).c_str());
2336 lenName = lenName + 2 + strlen((Supports[i]->getName()).c_str());
2337 lenDescription = lenDescription + 2 +
2338 strlen((Supports[i]->getDescription()).c_str());
2342 returnedSupportNameChar = new char[lenName];
2343 returnedSupportDescriptionChar = new char[lenDescription];
2345 returnedSupportNameChar = strcpy(returnedSupportNameChar,"Merge of ");
2346 returnedSupportDescriptionChar = strcpy(returnedSupportDescriptionChar,"Merge of ");
2348 returnedSupportNameChar = strcat(returnedSupportNameChar,(Supports[0]->getName()).c_str());
2349 returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
2350 (Supports[0]->getDescription()).c_str());
2352 for(int i = 1;i<size;i++)
2356 returnedSupportNameChar = strcat(returnedSupportNameChar," and ");
2357 returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar," and ");
2359 returnedSupportNameChar = strcat(returnedSupportNameChar,
2360 (Supports[i]->getName()).c_str());
2361 returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
2362 (Supports[i]->getDescription()).c_str());
2366 returnedSupportNameChar = strcat(returnedSupportNameChar,", ");
2367 returnedSupportNameChar = strcat(returnedSupportNameChar,
2368 (Supports[i]->getName()).c_str());
2370 returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,", ");
2371 returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
2372 (Supports[i]->getDescription()).c_str());
2376 returnedSupportName = string(returnedSupportNameChar);
2377 returnedSupport->setName(returnedSupportName);
2379 returnedSupportDescription = string(returnedSupportDescriptionChar);
2380 returnedSupport->setDescription(returnedSupportDescription);
2382 delete [] returnedSupportNameChar;
2383 delete [] returnedSupportDescriptionChar;
2387 return returnedSupport;
2391 return a SUPPORT pointer on the intersection of all SUPPORTs in Supports.
2392 The (SUPPORT *) NULL pointer is returned if the intersection is empty.
2393 You should delete this pointer after use to avois memory leaks.
2395 SUPPORT * MESH::intersectSupports(const vector<SUPPORT *> Supports) const throw (MEDEXCEPTION)
2397 const char * LOC = "MESH:::intersectSupports(const vector<SUPPORT *> ) : " ;
2400 SUPPORT * returnedSupport;
2401 string returnedSupportName;
2402 string returnedSupportDescription;
2403 char * returnedSupportNameChar;
2404 char * returnedSupportDescriptionChar;
2405 int size = Supports.size();
2409 MESSAGE(LOC <<" there is only one SUPPORT in the argument list, the method return a copy of this object !");
2410 SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
2412 returnedSupport = new SUPPORT(*obj);
2414 int lenName = strlen((Supports[0]->getName()).c_str()) + 8 + 1;
2415 int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 8 + 1;
2417 returnedSupportNameChar = new char[lenName];
2418 returnedSupportDescriptionChar = new char[lenDescription];
2420 returnedSupportNameChar = strcpy(returnedSupportNameChar,"Copy of ");
2421 returnedSupportNameChar = strcat(returnedSupportNameChar,(Supports[0]->getName()).c_str());
2422 returnedSupportDescriptionChar = strcpy(returnedSupportDescriptionChar,"Copy of ");
2423 returnedSupportDescriptionChar = strcat(returnedSupportDescriptionChar,
2424 (Supports[0]->getDescription()).c_str());
2426 returnedSupportName = string(returnedSupportNameChar);
2427 returnedSupportDescription = string(returnedSupportDescriptionChar);
2429 returnedSupport->setName(returnedSupportName);
2430 returnedSupport->setDescription(returnedSupportDescription);
2432 delete [] returnedSupportNameChar;
2433 delete [] returnedSupportDescriptionChar;
2437 SUPPORT * obj = const_cast <SUPPORT *> (Supports[0]);
2438 returnedSupport = new SUPPORT(*obj);
2440 int lenName = strlen((Supports[0]->getName()).c_str()) + 16 + 1;
2441 int lenDescription = strlen((Supports[0]->getDescription()).c_str()) + 16 + 1;
2443 for(int i = 1;i<size;i++)
2445 obj = const_cast <SUPPORT *> (Supports[i]);
2446 returnedSupport->intersecting(obj);
2450 lenName = lenName + 5 + strlen((Supports[i]->getName()).c_str());
2451 lenDescription = lenDescription + 5 +
2452 strlen((Supports[i]->getDescription()).c_str());
2456 lenName = lenName + 2 + strlen((Supports[i]->getName()).c_str());
2457 lenDescription = lenDescription + 2 +
2458 strlen((Supports[i]->getDescription()).c_str());
2462 if(returnedSupport != (SUPPORT *) NULL)
2464 returnedSupportNameChar = new char[lenName];
2465 returnedSupportDescriptionChar = new char[lenDescription];
2467 returnedSupportNameChar = strcpy(returnedSupportNameChar,
2468 "Intersection of ");
2469 returnedSupportDescriptionChar =
2470 strcpy(returnedSupportDescriptionChar,"Intersection of ");
2472 returnedSupportNameChar = strcat(returnedSupportNameChar,
2473 (Supports[0]->getName()).c_str());
2474 returnedSupportDescriptionChar =
2475 strcat(returnedSupportDescriptionChar,
2476 (Supports[0]->getDescription()).c_str());
2478 for(int i = 1;i<size;i++)
2482 returnedSupportNameChar = strcat(returnedSupportNameChar,
2484 returnedSupportDescriptionChar =
2485 strcat(returnedSupportDescriptionChar," and ");
2487 returnedSupportNameChar =
2488 strcat(returnedSupportNameChar,
2489 (Supports[i]->getName()).c_str());
2490 returnedSupportDescriptionChar =
2491 strcat(returnedSupportDescriptionChar,
2492 (Supports[i]->getDescription()).c_str());
2496 returnedSupportNameChar = strcat(returnedSupportNameChar,
2498 returnedSupportNameChar =
2499 strcat(returnedSupportNameChar,
2500 (Supports[i]->getName()).c_str());
2502 returnedSupportDescriptionChar =
2503 strcat(returnedSupportDescriptionChar,", ");
2504 returnedSupportDescriptionChar =
2505 strcat(returnedSupportDescriptionChar,
2506 (Supports[i]->getDescription()).c_str());
2510 returnedSupportName = string(returnedSupportNameChar);
2511 returnedSupport->setName(returnedSupportName);
2513 returnedSupportDescription = string(returnedSupportDescriptionChar);
2514 returnedSupport->setDescription(returnedSupportDescription);
2516 delete [] returnedSupportNameChar;
2517 delete [] returnedSupportDescriptionChar;
2522 return returnedSupport;