1 #include "MEDMEM_Connectivity.hxx"
2 #include "MEDMEM_Family.hxx"
3 #include "MEDMEM_Group.hxx"
4 #include "MEDMEM_CellModel.hxx"
6 #include "MEDMEM_SkyLineArray.hxx"
7 #include "MEDMEM_ModulusArray.hxx"
9 #include "MEDMEM_STRING.hxx"
13 using namespace MEDMEM;
14 using namespace MED_EN;
16 // Enlarge the vector if necessary to insert the element
17 static inline void insert_vector(vector<int> &Vect, int Indice, int Element)
19 if (Indice >= Vect.capacity())
20 Vect.reserve(Indice + 1000);
22 if (Indice >= Vect.size())
23 Vect.resize(Indice+1);
25 Vect[Indice] = Element;
29 Default Constructor. /n
30 Default for Entity is MED_CELL and type of Connectivity is MED_NODAL */
31 //--------------------------------------------------------------//
32 CONNECTIVITY::CONNECTIVITY(medEntityMesh Entity /* =MED_CELL */) :
33 //--------------------------------------------------------------//
35 _typeConnectivity(MED_NODAL),
37 _geometricTypes((medGeometryElement*)NULL),
38 _type((CELLMODEL*)NULL),
42 _nodal((MEDSKYLINEARRAY*)NULL),
43 _polygonsNodal((MEDSKYLINEARRAY*)NULL),
44 _polyhedronNodal((POLYHEDRONARRAY*)NULL),
45 _descending((MEDSKYLINEARRAY*)NULL),
46 _polygonsDescending((MEDSKYLINEARRAY*)NULL),
47 _polyhedronDescending((MEDSKYLINEARRAY*)NULL),
48 _reverseNodalConnectivity((MEDSKYLINEARRAY*)NULL),
49 _reverseDescendingConnectivity((MEDSKYLINEARRAY*)NULL),
50 _neighbourhood((MEDSKYLINEARRAY*)NULL),
51 _constituent((CONNECTIVITY*)NULL)
53 MESSAGE("CONNECTIVITY(medEntityMesh Entity=MED_CELL)");
58 Default for Entity is MED_CELL */
59 //------------------------------------------------------------------------------//
60 CONNECTIVITY::CONNECTIVITY(int numberOfTypes,medEntityMesh Entity /* =MED_CELL */):
61 //------------------------------------------------------------------------------//
63 _typeConnectivity(MED_NODAL),
64 _numberOfTypes(numberOfTypes),
67 _nodal((MEDSKYLINEARRAY*)NULL),
68 _polygonsNodal((MEDSKYLINEARRAY*)NULL),
69 _polyhedronNodal((POLYHEDRONARRAY*)NULL),
70 _descending((MEDSKYLINEARRAY*)NULL),
71 _polygonsDescending((MEDSKYLINEARRAY*)NULL),
72 _polyhedronDescending((MEDSKYLINEARRAY*)NULL),
73 _reverseNodalConnectivity((MEDSKYLINEARRAY*)NULL),
74 _reverseDescendingConnectivity((MEDSKYLINEARRAY*)NULL),
75 _neighbourhood((MEDSKYLINEARRAY*)NULL),
76 _constituent((CONNECTIVITY*)NULL)
78 MESSAGE("CONNECTIVITY(int numberOfTypes,medEntityMesh Entity=MED_CELL)");
79 _geometricTypes = new medGeometryElement[numberOfTypes];
80 _type = new CELLMODEL[numberOfTypes];
81 _count = new int[numberOfTypes+1];
87 //------------------------------------------------------//
88 CONNECTIVITY::CONNECTIVITY (const CONNECTIVITY & m):
89 //----------------------------------------------------//
91 _typeConnectivity (m._typeConnectivity),
92 _numberOfTypes (m._numberOfTypes),
93 _entityDimension (m._entityDimension),
94 _numberOfNodes (m._numberOfNodes)
96 if (m._geometricTypes != NULL)
98 _geometricTypes = new medGeometryElement[_numberOfTypes];
99 memcpy(_geometricTypes,m._geometricTypes,_numberOfTypes*sizeof(medGeometryElement));
102 _geometricTypes = (medGeometryElement *) NULL;
106 _type = new CELLMODEL[_numberOfTypes];
107 for (int i=0;i<_numberOfTypes;i++)
108 _type[i] = CELLMODEL(m._type[i]);
111 _type = (CELLMODEL *) NULL;
113 if (m._count != NULL)
115 _count = new int[_numberOfTypes+1];
116 memcpy(_count,m._count,(_numberOfTypes+1)*sizeof(int));
119 _count = (int *) NULL;
121 if (m._nodal != NULL)
122 _nodal = new MEDSKYLINEARRAY(* m._nodal);
124 _nodal = (MEDSKYLINEARRAY *) NULL;
126 if (m._polygonsNodal != NULL)
127 _polygonsNodal = new MEDSKYLINEARRAY(* m._polygonsNodal);
129 _polygonsNodal = (MEDSKYLINEARRAY *) NULL;
131 if (m._polyhedronNodal != NULL)
132 _polyhedronNodal = new POLYHEDRONARRAY(* m._polyhedronNodal);
134 _polyhedronNodal = (POLYHEDRONARRAY *) NULL;
136 if (m._descending != NULL)
137 _descending = new MEDSKYLINEARRAY(* m._descending);
139 _descending = (MEDSKYLINEARRAY *) NULL;
141 if (m._polygonsDescending != NULL)
142 _polygonsDescending = new MEDSKYLINEARRAY(* m._polygonsDescending);
144 _polygonsDescending = (MEDSKYLINEARRAY *) NULL;
146 if (m._polyhedronDescending != NULL)
147 _polyhedronDescending = new MEDSKYLINEARRAY(* m._polyhedronDescending);
149 _polyhedronDescending = (MEDSKYLINEARRAY *) NULL;
151 if (m._reverseNodalConnectivity != NULL)
152 _reverseNodalConnectivity = new MEDSKYLINEARRAY(* m._reverseNodalConnectivity);
154 _reverseNodalConnectivity = (MEDSKYLINEARRAY *) NULL;
156 if (m._reverseDescendingConnectivity != NULL)
157 _reverseDescendingConnectivity = new MEDSKYLINEARRAY(* m._reverseDescendingConnectivity);
159 _reverseDescendingConnectivity = (MEDSKYLINEARRAY *) NULL;
161 if (m._neighbourhood != NULL)
162 _neighbourhood = new MEDSKYLINEARRAY(* m._neighbourhood);
164 _neighbourhood = (MEDSKYLINEARRAY *) NULL;
166 if (m._constituent != NULL)
167 _constituent = new CONNECTIVITY(* m._constituent);
169 _constituent = (CONNECTIVITY *) NULL;
174 desallocates existing pointers */
175 //----------------------------//
176 CONNECTIVITY::~CONNECTIVITY()
177 //----------------------------//
179 MESSAGE("Destructeur de CONNECTIVITY()");
181 if (_geometricTypes != NULL)
182 delete [] _geometricTypes;
189 if (_polygonsNodal != NULL)
190 delete _polygonsNodal;
191 if (_polyhedronNodal != NULL)
192 delete _polyhedronNodal;
193 if (_descending != NULL)
195 if (_polygonsDescending != NULL)
196 delete _polygonsDescending;
197 if (_polyhedronDescending != NULL)
198 delete _polyhedronDescending;
199 if (_reverseNodalConnectivity != NULL)
200 delete _reverseNodalConnectivity;
201 if (_reverseDescendingConnectivity != NULL)
202 delete _reverseDescendingConnectivity;
203 if (_constituent != NULL)
208 set _constituent to Constituent
209 be aware desallocation of _constituent is done by CONNECTIVITY:~CONNECTIVITY
210 throws an exception if Constituent = MED_CELL
213 //----------------------------------------------------------//
214 void CONNECTIVITY::setConstituent(CONNECTIVITY * Constituent)
216 //----------------------------------------------------------//
218 medEntityMesh Entity = Constituent->getEntity();
219 if (Entity == MED_CELL)
220 throw MEDEXCEPTION(LOCALIZED("CONNECTIVITY::setConstituent : could not set constituent on MED_CELL !"));
222 if ((Entity == MED_EDGE)&(_entityDimension == 3))
224 if (_constituent == NULL)
225 throw MEDEXCEPTION(LOCALIZED("CONNECTIVITY::setConstituent : Entity not found !"));
226 _constituent->setConstituent(Constituent);
229 _constituent = Constituent;
232 /*! Duplicated Types array in CONNECTIVITY object. */
233 //--------------------------------------------------------------------//
234 void CONNECTIVITY::setGeometricTypes(const medGeometryElement * Types,
235 const medEntityMesh Entity)
237 //--------------------------------------------------------------------//
239 if (Entity == _entity)
240 for (int i=0; i<_numberOfTypes; i++)
242 _geometricTypes[i] = Types[i];
243 _type[i] = CELLMODEL(_geometricTypes[i]);
247 if (_constituent == NULL)
248 throw MEDEXCEPTION(LOCALIZED("CONNECTIVITY::setGeometricTypes : Entity not found !"));
249 _constituent->setGeometricTypes(Types,Entity);
254 //--------------------------------------------------------------------//
255 void CONNECTIVITY::setCount(const int * Count, const medEntityMesh Entity)
257 //--------------------------------------------------------------------//
259 if (Entity == _entity)
261 int * index = new int[Count[_numberOfTypes]];
264 for (int i=0; i<_numberOfTypes; i++) {
265 _count[i+1] = Count[i+1];
266 int NumberOfNodesPerElement = _type[i].getNumberOfNodes() ;
267 for (int j=_count[i]; j<_count[i+1]; j++)
268 index[j] = index[j-1]+NumberOfNodesPerElement;
271 if (_nodal != NULL) delete _nodal;
272 _nodal = new MEDSKYLINEARRAY(_count[_numberOfTypes]-1,index[_count[_numberOfTypes]-1]-1);
273 _nodal->setIndex(index);
278 if (_constituent == NULL)
279 throw MEDEXCEPTION(LOCALIZED("CONNECTIVITY::setCount : Entity not found !"));
280 _constituent->setCount(Count,Entity);
284 //--------------------------------------------------------------------//
285 void CONNECTIVITY::setNodal(const int * Connectivity,
286 const medEntityMesh Entity,
287 const medGeometryElement Type)
289 //--------------------------------------------------------------------//
291 if (_entity == Entity)
293 // find geometric type
295 for (int i=0; i<_numberOfTypes; i++)
296 if (_geometricTypes[i] == Type) {
298 int NumberOfNodePerElements = _type[i].getNumberOfNodes() ;
299 //_nodal->setI(i+1,Connectivity);
300 for( int j=_count[i];j<_count[i+1]; j++)
301 _nodal->setI(j,Connectivity+(j-_count[i])*NumberOfNodePerElements);
304 throw MEDEXCEPTION(LOCALIZED("CONNECTIVITY::setNodal : geometric type not found"));
307 if (_constituent == NULL)
308 throw MEDEXCEPTION(LOCALIZED("CONNECTIVITY::setNodal : Entity not found !"));
309 _constituent->setNodal(Connectivity,Entity,Type);
314 //--------------------------------------------------------------------//
315 void CONNECTIVITY::setPolygonsConnectivity(medConnectivity ConnectivityType, medEntityMesh Entity, const int* PolygonsConnectivity, const int* PolygonsConnectivityIndex, int ConnectivitySize, int NumberOfPolygons)
316 //--------------------------------------------------------------------//
318 const char* LOC = "CONNECTIVITY::setPolygonsConnectivity";
321 if (_entity == Entity)
323 MEDSKYLINEARRAY* Connectivity = new MEDSKYLINEARRAY(NumberOfPolygons,ConnectivitySize,PolygonsConnectivityIndex,PolygonsConnectivity);
325 if (ConnectivityType == MED_NODAL)
327 if (_polygonsNodal != (MEDSKYLINEARRAY*) NULL)
328 delete _polygonsNodal;
329 _polygonsNodal = Connectivity;
333 if (_typeConnectivity != MED_DESCENDING)
334 _typeConnectivity = MED_DESCENDING; //by default it is set to MED_NODAL
335 if (_polygonsDescending != (MEDSKYLINEARRAY*) NULL)
336 delete _polygonsDescending;
337 _polygonsDescending = Connectivity;
342 if (_constituent == (CONNECTIVITY*) NULL)
343 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not found !"));
344 _constituent->setPolygonsConnectivity(ConnectivityType, Entity, PolygonsConnectivity, PolygonsConnectivityIndex, ConnectivitySize, NumberOfPolygons);
349 //--------------------------------------------------------------------//
350 void CONNECTIVITY::setPolyhedronConnectivity(medConnectivity ConnectivityType, const int* PolyhedronConnectivity, const int* PolyhedronIndex, int ConnectivitySize, int NumberOfPolyhedron, const int* PolyhedronFacesIndex /* =(int*)NULL */, int NumberOfFaces /* =0 */)
351 //--------------------------------------------------------------------//
353 const char* LOC = "CONNECTIVITY::setPolyhedronConnectivity";
356 if (_entity == MED_CELL)
358 if (ConnectivityType == MED_NODAL)
360 if (_polyhedronNodal != (POLYHEDRONARRAY*) NULL)
361 delete _polyhedronNodal;
362 _polyhedronNodal = new POLYHEDRONARRAY(NumberOfPolyhedron,NumberOfFaces,ConnectivitySize);
363 _polyhedronNodal->setPolyhedronIndex(PolyhedronIndex);
364 _polyhedronNodal->setFacesIndex(PolyhedronFacesIndex);
365 _polyhedronNodal->setNodes(PolyhedronConnectivity);
369 if (_typeConnectivity != MED_DESCENDING)
370 _typeConnectivity = MED_DESCENDING; //by default it is set to MED_NODAL
371 if (_polyhedronDescending != (MEDSKYLINEARRAY*) NULL)
372 delete _polyhedronDescending;
373 _polyhedronDescending = new MEDSKYLINEARRAY(NumberOfPolyhedron,ConnectivitySize,PolyhedronIndex,PolyhedronConnectivity);
377 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : _entity must be MED_CELL to set polyhedron !"));
382 //------------------------------------------------------------------------------------------//
383 void CONNECTIVITY::calculateConnectivity(medConnectivity ConnectivityType, medEntityMesh Entity)
384 //------------------------------------------------------------------------------------------//
386 MESSAGE("CONNECTIVITY::calculateConnectivity");
388 // a temporary limitation due to calculteDescendingConnectivity function !!!
389 if ((_entityDimension==3) & (Entity==MED_EDGE))
390 throw MEDEXCEPTION("CONNECTIVITY::calculateConnectivity : Could not build EDGE !");
393 if (ConnectivityType==MED_NODAL)
394 calculateNodalConnectivity();
396 if (Entity==MED_CELL)
397 calculateDescendingConnectivity();
399 throw MEDEXCEPTION("CONNECTIVITY::calculateConnectivity : Could not build DESCENDING !");
400 if (Entity!=_entity) {
401 calculateDescendingConnectivity();
402 if (_entityDimension == 2 || _entityDimension == 3)
403 _constituent->calculateConnectivity(ConnectivityType,Entity);
407 /*! Give, in full or no interlace mode (for nodal connectivity),
408 descending or nodal connectivity.
410 You must give a %medEntityMesh (ie:MED_EDGE)
411 and a %medGeometryElement (ie:MED_SEG3). */
413 //------------------------------------------------------------//
414 void CONNECTIVITY::updateFamily(const vector<FAMILY*>& myFamilies)
415 //------------------------------------------------------------//
417 const char * LOC = "CONNECTIVITY::updateFamily(vector<FAMILY*>) ";
418 int numberOfFamilies = myFamilies.size();
419 if (numberOfFamilies == 0 || _constituent == NULL)
421 // does we do an update ?
422 if ((_constituent != NULL) && (_descending != NULL))
424 if (myFamilies[0]->getEntity() != _constituent->getEntity())
426 CONNECTIVITY * oldConstituent = _constituent;
427 _constituent = (CONNECTIVITY *)NULL;
428 if (oldConstituent->_nodal==NULL)
429 throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"We have no nodal connectivity of sub cell"));
431 //Loc vars defined to treat polygons exactly the same as classic types. Not nice but necessary.
432 int oldNumberOfFaceTab[2];
433 const int * oldConstituentValueTab[2];
434 const int * oldConstituentIndexTab[2];
435 int * renumberingFromOldToNewTab[2];//Final mapping array between old numbers and new numbers.;
437 int oldNumberOfFace = oldConstituent->_nodal->getNumberOf(); oldNumberOfFaceTab[0]=oldNumberOfFace;
438 const int * oldConstituentValue = oldConstituent->_nodal->getValue(); oldConstituentValueTab[0]=oldConstituentValue;
439 const int * oldConstituentIndex = oldConstituent->_nodal->getIndex(); oldConstituentIndexTab[0]=oldConstituentIndex;
440 int * renumberingFromOldToNew= new int [oldNumberOfFace]; renumberingFromOldToNewTab[0]=renumberingFromOldToNew;
442 int oldNumberOfFacePoly = oldConstituent->getNumberOfPolygons();
443 const int * oldConstituentValuePoly=0;
444 const int * oldConstituentIndexPoly=0;
445 int * renumberingFromOldToNewPoly=0;
447 int nbOfTurnInGlobalLoop=1;//Defined to know if a second search on polygons is needed.
448 if(oldNumberOfFacePoly>0)
450 oldNumberOfFaceTab[1]=oldNumberOfFacePoly;
451 oldConstituentValuePoly = oldConstituent->_polygonsNodal->getValue(); oldConstituentValueTab[1]=oldConstituentValuePoly;
452 oldConstituentIndexPoly = oldConstituent->_polygonsNodal->getIndex(); oldConstituentIndexTab[1]=oldConstituentIndexPoly;
453 renumberingFromOldToNewPoly=new int[oldNumberOfFacePoly]; renumberingFromOldToNewTab[1]=renumberingFromOldToNewPoly;
454 nbOfTurnInGlobalLoop++;
457 calculateDescendingConnectivity();//perform calculation of descending connectivity to update all connectivities after having taken into account of direction of faces...
458 _constituent->calculateReverseNodalConnectivity();//getting d-1 nodal connectivity to find get new face numbers from nodes numbers...
460 const int * reverseFaceNodal = _constituent->getReverseNodalConnectivity(); //Common to polygons and classic geometric types
461 const int * reverseFaceNodalIndex = _constituent->getReverseNodalConnectivityIndex(); //Common to polygons and classic geometric types
463 for(int loop=0;loop<nbOfTurnInGlobalLoop;loop++)
465 int oldNumberOfFaceLoop=oldNumberOfFaceTab[loop];
466 const int * oldConstituentValueLoop=oldConstituentValueTab[loop];
467 const int * oldConstituentIndexLoop= oldConstituentIndexTab[loop];
468 int * renumberingFromOldToNewLoop=renumberingFromOldToNewTab[loop];
469 for(int iOldFace=0;iOldFace<oldNumberOfFaceLoop;iOldFace++)
471 const int *nodesOfCurrentFaceOld=oldConstituentValueLoop+oldConstituentIndexLoop[iOldFace]-1;
472 int nbOfNodesOfCurrentFaceOld=oldConstituentIndexLoop[iOldFace+1]-oldConstituentIndexLoop[iOldFace];
474 //retrieving new number of polygon face...
475 int sizeOfNewFaceNb=reverseFaceNodalIndex[nodesOfCurrentFaceOld[0]]-reverseFaceNodalIndex[nodesOfCurrentFaceOld[0]-1];
476 int *newFaceNb=new int[sizeOfNewFaceNb];
477 memcpy(newFaceNb,reverseFaceNodal+reverseFaceNodalIndex[nodesOfCurrentFaceOld[0]-1]-1,sizeOfNewFaceNb*sizeof(int));//all faces candidates are in newFaceNb !!!
478 for(int curNode=1;curNode<nbOfNodesOfCurrentFaceOld && sizeOfNewFaceNb>1;curNode++)
480 const int *newFacesNbForCurNode=reverseFaceNodal+reverseFaceNodalIndex[nodesOfCurrentFaceOld[curNode]-1]-1;
481 int sizeOfNewFacesNbForCurNode=reverseFaceNodalIndex[nodesOfCurrentFaceOld[curNode]]-reverseFaceNodalIndex[nodesOfCurrentFaceOld[curNode]-1];
482 for(int i=0;i<sizeOfNewFaceNb;)
485 for(int j=0;j<sizeOfNewFacesNbForCurNode && !found;j++)
487 if(newFacesNbForCurNode[j]==newFaceNb[i])
493 newFaceNb[i]=newFaceNb[--sizeOfNewFaceNb];
496 if(sizeOfNewFaceNb!=1)
497 throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"unexisting face specified!!!"));
498 renumberingFromOldToNewLoop[iOldFace]=newFaceNb[0];
500 //end of retrieving new number of polygon face...
501 int nbOfNodesOfCurrentFaceNew;
502 const int *nodesOfCurrentFaceNew=_constituent->getConnectivityOfAnElementWithPoly(MED_NODAL,_constituent->getEntity(),
503 renumberingFromOldToNewLoop[iOldFace],nbOfNodesOfCurrentFaceNew);
504 MEDMODULUSARRAY modulusArrayOld(nbOfNodesOfCurrentFaceOld,nodesOfCurrentFaceOld);
505 MEDMODULUSARRAY modulusArrayNew(nbOfNodesOfCurrentFaceNew,nodesOfCurrentFaceNew);
506 int retCompareNewOld=modulusArrayNew.compare(modulusArrayOld);
507 if(retCompareNewOld==0)
508 throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"Uncompatible given user face with calculated existing faces"));
509 if(retCompareNewOld==-1)
510 invertConnectivityForAFace(renumberingFromOldToNewLoop[iOldFace],nodesOfCurrentFaceOld,loop==1);
513 // Updating the Family
514 for(vector<FAMILY*>::const_iterator iter=myFamilies.begin();iter!=myFamilies.end();iter++)
515 (*iter)->changeElementsNbs(_constituent->getEntity(),renumberingFromOldToNew,oldNumberOfFace,renumberingFromOldToNewPoly);
516 delete oldConstituent ;
517 delete [] renumberingFromOldToNew;
518 if(oldNumberOfFacePoly>0)
519 delete [] renumberingFromOldToNewPoly;
523 //------------------------------------------------------------------------------------------------------------------//
524 const int * MEDMEM::CONNECTIVITY::getConnectivity(medConnectivity ConnectivityType, medEntityMesh Entity, medGeometryElement Type)
525 //------------------------------------------------------------------------------------------------------------------//
527 const char * LOC = "CONNECTIVITY::getConnectivity";
530 MEDSKYLINEARRAY * Connectivity;
531 if (Entity==_entity) {
533 if (ConnectivityType==MED_NODAL)
535 calculateNodalConnectivity();
540 calculateDescendingConnectivity();
541 Connectivity=_descending;
544 if (Connectivity!=NULL)
545 if (Type==MED_ALL_ELEMENTS)
546 return Connectivity->getValue();
548 for (int i=0; i<_numberOfTypes; i++)
549 if (_geometricTypes[i]==Type)
550 //return Connectivity->getI(i+1);
551 return Connectivity->getI(_count[i]);
552 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Type not found !"));
555 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Connectivity not defined !"));
557 if (_constituent != NULL)
558 return _constituent->getConnectivity(ConnectivityType,Entity,Type);
560 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !"));
563 //------------------------------------------------------------------------------------------------------------------//
564 int CONNECTIVITY::getConnectivityLength(medConnectivity ConnectivityType, medEntityMesh Entity, medGeometryElement Type)
565 //------------------------------------------------------------------------------------------------------------------//
567 const char * LOC = "CONNECTIVITY::getConnectivity";
570 MEDSKYLINEARRAY * Connectivity;
571 if (Entity==_entity) {
573 if (ConnectivityType==MED_NODAL)
575 calculateNodalConnectivity();
580 calculateDescendingConnectivity();
581 Connectivity=_descending;
584 if (Connectivity!=NULL)
585 if (Type==MED_ALL_ELEMENTS)
586 return Connectivity->getLength();
588 for (int i=0; i<_numberOfTypes; i++)
589 if (_geometricTypes[i]==Type)
590 return Connectivity->getNumberOfI(_count[i]);
591 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Type not found !"));
594 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Connectivity not defined !"));
597 if (_constituent != NULL)
598 return _constituent->getConnectivityLength(ConnectivityType,Entity,Type);
600 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !"));
603 /*! Give morse index array to use with
604 getConnectivity(MED_FULL_INTERLACE,mode,entity,MED_ALL_ELEMENTS).
606 Each value give start index for corresponding entity in connectivity array.
608 Example : i-th element, j-th node of it :
609 - In C mode : Connectivity[ConnectivityIndex[i]-1+j-1]
610 - In fortran mode : Connectivity[ConnectivityIndex[i]+j] */
611 //-----------------------------------------------------------------------------------------------//
612 const int * CONNECTIVITY::getConnectivityIndex(medConnectivity ConnectivityType, medEntityMesh Entity)
613 //-----------------------------------------------------------------------------------------------//
615 const char * LOC = "CONNECTIVITY::getConnectivityIndex";
618 MEDSKYLINEARRAY * Connectivity;
619 if (Entity==_entity) {
621 if (ConnectivityType==MED_NODAL)
624 Connectivity=_descending;
626 if (Connectivity!=NULL)
627 return Connectivity->getIndex();
629 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Connectivity not defined !"));
632 if (_constituent != NULL)
633 return _constituent->getConnectivityIndex(ConnectivityType,Entity);
635 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !"));
639 //-------------------------------------------------------------//
640 const int* CONNECTIVITY::getPolygonsConnectivity(medConnectivity ConnectivityType, medEntityMesh Entity)
641 //-------------------------------------------------------------//
643 const char* LOC = "CONNECTIVITY::getPolygonsConnectivity";
646 MEDSKYLINEARRAY* Connectivity;
647 if (Entity == _entity)
649 if (ConnectivityType == MED_NODAL)
651 calculateNodalConnectivity();
652 Connectivity = _polygonsNodal;
656 calculateDescendingConnectivity();
657 Connectivity = _polygonsDescending;
659 if (Connectivity != NULL)
660 return Connectivity->getValue();
662 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Polygons Connectivity not defined !"));
665 if (_constituent != NULL)
666 return _constituent->getPolygonsConnectivity(ConnectivityType, Entity);
667 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !"));
671 //-------------------------------------------------------------//
672 const int* CONNECTIVITY::getPolygonsConnectivityIndex(medConnectivity ConnectivityType, medEntityMesh Entity)
673 //-------------------------------------------------------------//
675 const char* LOC = "CONNECTIVITY::getPolygonsConnectivityIndex";
678 MEDSKYLINEARRAY* Connectivity;
679 if (Entity == _entity)
681 if (ConnectivityType == MED_NODAL)
683 // calculateNodalConnectivity();
684 Connectivity = _polygonsNodal;
688 // calculateDescendingConnectivity();
689 Connectivity = _polygonsDescending;
691 if (Connectivity != NULL)
692 return Connectivity->getIndex();
694 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Polygons Connectivity not defined !"));
697 if (_constituent != NULL)
698 return _constituent->getPolygonsConnectivityIndex(ConnectivityType, Entity);
699 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !"));
703 /*! We suppose in this method that nodal and descending connectivities
705 //-------------------------------------------------------------//
706 int CONNECTIVITY::getNumberOfPolygons() const
707 //-------------------------------------------------------------//
709 if (_polygonsNodal != (MEDSKYLINEARRAY*) NULL)
710 return _polygonsNodal->getNumberOf();
711 else if (_polygonsDescending != (MEDSKYLINEARRAY*) NULL)
712 return _polygonsDescending->getNumberOf();
718 //--------------------------------------------------------------//
719 const int* CONNECTIVITY::getPolyhedronConnectivity(medConnectivity ConnectivityType) const
720 //--------------------------------------------------------------//
722 const char* LOC = "CONNECTIVITY::getPolyhedronConnectivity";
725 if (_entity == MED_CELL) //polyhedron can only be MED_CELL
727 if (ConnectivityType == MED_NODAL)
729 ((CONNECTIVITY *)(this))->calculateNodalConnectivity();
730 if (_polyhedronNodal != (POLYHEDRONARRAY*) NULL)
731 return _polyhedronNodal->getNodes();
733 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Polyhedron Nodal Connectivity not defined !"));
737 ((CONNECTIVITY *)(this))->calculateDescendingConnectivity();
738 if (_polyhedronDescending != (MEDSKYLINEARRAY*) NULL)
739 return _polyhedronDescending->getValue();
741 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Polyhedron Descending Connectivity not defined !"));
744 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : No polyhedron for _entity different from MED_CELL !"));
748 //---------------------------------------------------------------//
749 const int* CONNECTIVITY::getPolyhedronFacesIndex() const
750 //---------------------------------------------------------------//
752 const char* LOC = "CONNECTIVITY::getPolyhedronFacesIndex";
755 if (_entity == MED_CELL) //polyhedron can only be MED_CELL
757 // calculateNodalConnectivity();
758 if (_polyhedronNodal != (POLYHEDRONARRAY*) NULL)
759 return _polyhedronNodal->getFacesIndex();
761 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : No Polyhedron in that Connectivity !"));
763 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : No polyhedron for _entity different from MED_CELL !"));
767 //--------------------------------------------------------------//
768 const int* CONNECTIVITY::getPolyhedronIndex(medConnectivity ConnectivityType) const
769 //--------------------------------------------------------------//
771 const char* LOC = "CONNECTIVITY::getPolyhedronIndex";
774 if (_entity == MED_CELL) //polyhedron can only be MED_CELL
776 if (ConnectivityType == MED_NODAL)
778 // calculateNodalConnectivity();
779 if (_polyhedronNodal != (POLYHEDRONARRAY*) NULL)
780 return _polyhedronNodal->getPolyhedronIndex();
782 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Polyhedron Nodal Connectivity not defined !"));
786 // calculateDescendingConnectivity();
787 if (_polyhedronDescending != (MEDSKYLINEARRAY*) NULL)
788 return _polyhedronDescending->getIndex();
790 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Polyhedron Descending Connectivity not defined !"));
793 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : No polyhedron for _entity different from MED_CELL !"));
797 /*! We suppose in this method that nodal and descending connectivities
799 //-------------------------------------------------------------//
800 int CONNECTIVITY::getNumberOfPolyhedronFaces() const
801 //-------------------------------------------------------------//
803 // if (_polyhedronNodal == (POLYHEDRONARRAY*) NULL)
804 // calculateNodalConnectivity();
805 if (_polyhedronNodal != (POLYHEDRONARRAY*) NULL)
806 return _polyhedronNodal->getNumberOfFaces();
812 /*! We suppose in this method that nodal and descending connectivities
814 //--------------------------------------------------------------//
815 int CONNECTIVITY::getNumberOfPolyhedron() const
816 //--------------------------------------------------------------//
818 if (_polyhedronNodal != (POLYHEDRONARRAY*) NULL)
819 return _polyhedronNodal->getNumberOfPolyhedron();
820 else if (_polyhedronDescending != (MEDSKYLINEARRAY*) NULL)
821 return _polyhedronDescending->getNumberOf();
828 //--------------------------------------------------------------//
829 const CELLMODEL & CONNECTIVITY::getType(medGeometryElement Type) const
830 //--------------------------------------------------------------//
832 const char * LOC = "CONNECTIVITY::getType";
835 if ((Type==MED_ALL_ELEMENTS) || (Type==MED_NONE))
836 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement must be different of MED_ALL_ELEMENTS and MED_NONE !"));
837 for (int i=0; i<_numberOfTypes; i++)
838 if (_geometricTypes[i]==Type)
840 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement not found !"));
843 /*! Returns the number of elements of type %medGeometryElement.
844 Note : not implemented for MED_ALL_ELEMENTS nor for MED_NONE */
845 //------------------------------------------------------------------------//
846 int CONNECTIVITY::getNumberOfNodesInType(medGeometryElement Type) const
847 //------------------------------------------------------------------------//
849 const char * LOC = "CONNECTIVITY::getNumberOfNodesInType";
852 if ((Type==MED_ALL_ELEMENTS) || (Type==MED_NONE))
853 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement must be different of MED_ALL_ELEMENTS and MED_NONE!"));
854 for (int i=0; i<_numberOfTypes; i++)
855 if (_geometricTypes[i]==Type)
856 return _type[i].getNumberOfNodes();
857 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement not found !"));
860 /*! Returns the number of geometric sub cells of %medGeometryElement type.
861 not implemented for MED_ALL_ELEMENTS nor for MED_NONE */
862 //------------------------------------------------------------------------//
863 int CONNECTIVITY::getNumberOfSubCellInType(medGeometryElement Type) const
864 //------------------------------------------------------------------------//
866 if ((Type==MED_ALL_ELEMENTS) || (Type==MED_NONE))
867 throw MEDEXCEPTION("CONNECTIVITY::getNumberOfSubCell : medGeometryElement must be different of MED_ALL_ELEMENTS and MED_NONE!");
868 for (int i=0; i<_numberOfTypes; i++)
869 if (_geometricTypes[i]==Type)
870 return _type[i].getNumberOfConstituents(1);
871 throw MEDEXCEPTION("CONNECTIVITY::getNumberOfSubCell : medGeometryElement not found !");
874 /*! Returns the number of elements of type %medGeometryElement.
877 - Implemented for MED_ALL_ELEMENTS
878 - Not implemented for MED_ALL_ENTITIES (A VERIFIER)
879 - Not implemented for MED_NONE */
880 //-----------------------------------------------------------------------------------//
881 int CONNECTIVITY::getNumberOf(medEntityMesh Entity, medGeometryElement Type) const
882 //-----------------------------------------------------------------------------------//
884 //const char * LOC = "CONNECTIVITY::getNumberOf";
885 if (Entity==_entity) {
887 return 0; // not defined !
888 //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement must be different of MED_NONE"));
889 if (!existConnectivity(MED_NODAL,Entity) && !existConnectivity(MED_DESCENDING,Entity)) return 0; //case with only polygons for example
890 if (Type==MED_ALL_ELEMENTS)
891 return _count[_numberOfTypes]-1;
892 for (int i=0; i<_numberOfTypes; i++)
893 if (_geometricTypes[i]==Type)
894 return (_count[i+1] - _count[i]);
895 //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement not found !"));
897 if (_constituent != NULL)
898 return _constituent->getNumberOf(Entity,Type);
900 return 0; // valid if they are nothing else !
901 //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !"));
905 //--------------------------------------------------------------//
906 const int* CONNECTIVITY::getValue(medConnectivity TypeConnectivity,
907 medGeometryElement Type)
908 //--------------------------------------------------------------//
910 if (TypeConnectivity == MED_NODAL)
912 calculateNodalConnectivity();
913 if (Type==MED_ALL_ELEMENTS)
914 return _nodal->getValue();
915 for (int i=0; i<_numberOfTypes; i++)
916 if (_geometricTypes[i]==Type)
917 return _nodal->getI(_count[i]);
921 calculateDescendingConnectivity();
922 if (Type==MED_ALL_ELEMENTS)
923 return _descending->getValue();
924 for (int i=0; i<_numberOfTypes; i++)
925 if (_geometricTypes[i]==Type)
926 return _descending->getI(Type);
928 throw MEDEXCEPTION("Not found");
932 //---------------------------------------------------------------------//
933 const int* CONNECTIVITY:: getValueIndex(medConnectivity TypeConnectivity)
934 //---------------------------------------------------------------------//
936 if (TypeConnectivity == MED_NODAL)
938 calculateNodalConnectivity();
939 return _nodal->getIndex();
943 calculateDescendingConnectivity();
944 return _descending->getIndex();
948 /*! Not yet implemented */
949 //----------------------------------------------//
950 const int* CONNECTIVITY:: getNeighbourhood() const
951 //----------------------------------------------//
953 throw MEDEXCEPTION("CONNECTIVITY::getNeighbourhood : Not implemented");
956 /*! Returns an array which contains, for each node, all cells
958 //-------------------------------------------------//
959 const int* CONNECTIVITY::getReverseNodalConnectivity()
960 //-------------------------------------------------//
962 calculateReverseNodalConnectivity();
963 return _reverseNodalConnectivity->getValue();
966 /*! Give index array to use with getReverseConnectivity(MED_NODAL).
967 It is unusefull with MED_DESCENDING mode, because we have allways two cells. */
968 //-------------------------------------------------------//
969 const int* CONNECTIVITY::getReverseNodalConnectivityIndex()
970 //-------------------------------------------------------//
972 calculateReverseNodalConnectivity();
973 return _reverseNodalConnectivity->getIndex();
976 /*! Returns an array which contains, for each face (or edge),
977 the 2 cells of each side. First is cell which face normal is outgoing.
979 //------------------------------------------------------//
980 const int* CONNECTIVITY::getReverseDescendingConnectivity()
981 //------------------------------------------------------//
983 // it is in _constituent connectivity only if we are in MED_CELL
984 // (we could not for instance calculate face-edge connectivity !)
985 if (_entity!=MED_CELL)
986 throw MEDEXCEPTION("CONNECTIVITY::getReverseDescendingConnectivity : Error Only in MED_CELL connectivity");
988 // we need descending connectivity
989 calculateDescendingConnectivity();
990 return _reverseDescendingConnectivity->getValue();
993 /*! calculate the reverse descending Connectivity
994 and returns the index ( A DOCUMENTER MIEUX)*/
995 //-----------------------------------------------------------//
996 const int* CONNECTIVITY::getReverseDescendingConnectivityIndex()
997 //-----------------------------------------------------------//
999 // it is in _constituent connectivity only if we are in MED_CELL
1000 if (_entity!=MED_CELL)
1001 throw MEDEXCEPTION("CONNECTIVITY::getReverseDescendingConnectivityIndex : Error Only in MED_CELL connectivity");
1003 // we need descending connectivity
1004 calculateDescendingConnectivity();
1005 return _reverseDescendingConnectivity->getIndex();
1008 /*! A DOCUMENTER (et a finir ???) */
1009 //--------------------------------------------//
1010 void CONNECTIVITY::calculateNodalConnectivity()
1011 //--------------------------------------------//
1013 if (_nodal==NULL && _polygonsNodal==NULL && _polyhedronNodal==NULL)
1015 if (_descending==NULL && _polygonsDescending==NULL && _polyhedronDescending==NULL)
1016 throw MEDEXCEPTION("CONNECTIVITY::calculateNodalConnectivity : No connectivity found !");
1017 // calculate _nodal, _polygonsNodal and _polyhedronNodal from _descending, _polygonsDescending and _polyhedronDescending
1021 /*! If not yet done, calculate the nodal Connectivity
1022 and the reverse nodal Connectivity*/
1023 //---------------------------------------------------//
1024 void CONNECTIVITY::calculateReverseNodalConnectivity()
1025 //---------------------------------------------------//
1027 const char * LOC = "CONNECTIVITY::calculateReverseNodalConnectivity : ";
1031 SCRUTE(_reverseNodalConnectivity);
1034 calculateNodalConnectivity();
1036 if(_reverseNodalConnectivity==NULL) {
1038 int node_number = 0;
1039 vector <vector <int> > reverse_connectivity;
1040 reverse_connectivity.resize(_numberOfNodes+1);
1042 // Treat all cells types
1044 for (int j = 0; j < _numberOfTypes; j++)
1046 // node number of the cell type
1047 node_number = _type[j].getNumberOfNodes();
1048 // treat all cells of a particular type
1049 for (int k = _count[j]; k < _count[j+1]; k++)
1050 // treat all nodes of the cell type
1051 for (int local_node_number = 1; local_node_number < node_number+1; local_node_number++)
1053 int global_node = _nodal->getIJ(k,local_node_number);
1054 reverse_connectivity[global_node].push_back(k);
1058 if(_polygonsNodal != NULL && _entityDimension==2)
1060 int nbOfPolygons=_polygonsNodal->getNumberOf();
1061 int offset=getNumberOf(_entity,MED_ALL_ELEMENTS);
1062 const int *index=_polygonsNodal->getIndex();
1063 const int *value=_polygonsNodal->getValue();
1064 for (int local_polyg_number = 0; local_polyg_number < nbOfPolygons; local_polyg_number++)
1066 for(int i=index[local_polyg_number];i<index[local_polyg_number+1];i++)
1067 reverse_connectivity[value[i-1]].push_back(offset+local_polyg_number+1);
1071 if(_polyhedronNodal != NULL && _entityDimension==3)
1073 int nbOfPolyhedra=_polyhedronNodal->getNumberOfPolyhedron();
1074 int offset=getNumberOf(_entity,MED_ALL_ELEMENTS);
1075 for (int local_polyh_number = 0; local_polyh_number < nbOfPolyhedra; local_polyh_number++)
1078 int global_polyh_number=offset+local_polyh_number+1;
1079 int *nodes=getNodesOfPolyhedron(global_polyh_number,nbOfNodes);
1080 for(int i=0;i<nbOfNodes;i++)
1081 reverse_connectivity[nodes[i]].push_back(global_polyh_number);
1086 // Full reverse_nodal_connectivity and reverse_nodal_connectivity_index
1088 //calculate size of reverse_nodal_connectivity array
1089 int size_reverse_nodal_connectivity = 0;
1090 for (int i = 1; i < _numberOfNodes+1; i++)
1091 size_reverse_nodal_connectivity += reverse_connectivity[i].size();
1093 int * reverse_nodal_connectivity_index = new int[_numberOfNodes+1];
1094 int * reverse_nodal_connectivity = new int[size_reverse_nodal_connectivity];
1096 reverse_nodal_connectivity_index[0] = 1;
1097 for (int i = 1; i < _numberOfNodes+1; i++)
1099 int size = reverse_connectivity[i].size();
1100 reverse_nodal_connectivity_index[i] = reverse_nodal_connectivity_index[i-1] + size;
1101 for (int j = 0; j < size; j++)
1102 reverse_nodal_connectivity[reverse_nodal_connectivity_index[i-1]-1+j]= reverse_connectivity[i][j];
1105 _reverseNodalConnectivity = new MEDSKYLINEARRAY(_numberOfNodes,size_reverse_nodal_connectivity,
1106 reverse_nodal_connectivity_index,
1107 reverse_nodal_connectivity,true);
1112 /*! If not yet done, calculate the Descending Connectivity */
1113 //-------------------------------------------------//
1114 void CONNECTIVITY::calculateDescendingConnectivity()
1115 //-------------------------------------------------//
1117 const char * LOC = "CONNECTIVITY::calculateDescendingConnectivity() : ";
1120 if (_descending==NULL && _polygonsDescending==NULL && _polyhedronDescending==NULL)
1124 MESSAGE(LOC<<"No connectivity found !");
1125 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"No connectivity found !"));
1127 // calcul _descending from _nodal
1128 // we need CONNECTIVITY for constituent
1130 if (_constituent != NULL)
1131 // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"ERROR : No descending connectivity defined, but constituent exist !"));
1132 MESSAGE(LOC<<": No descending connectivity defined, but constituent exist !");
1134 if (_entityDimension == 3)
1135 _constituent = new CONNECTIVITY(MED_FACE);
1136 else if (_entityDimension == 2)
1137 _constituent = new CONNECTIVITY(MED_EDGE);
1140 MESSAGE(LOC<<"We are in 1D");
1144 _constituent->_typeConnectivity = MED_NODAL;
1145 _constituent->_numberOfNodes = _numberOfNodes;
1146 // foreach cells, we built array of constituent
1147 int DescendingSize = 0;
1148 for(int i=0; i<_numberOfTypes; i++)
1149 DescendingSize+=(_count[i+1]-_count[i])*_type[i].getNumberOfConstituents(1);
1150 //_descending = new MEDSKYLINEARRAY(_count[_numberOfTypes]-1,DescendingSize);
1151 //const int * descend_connectivity = _descending->getValue();
1152 int * descend_connectivity = new int[DescendingSize];
1153 for (int i=0; i<DescendingSize; i++)
1154 descend_connectivity[i]=0;
1155 //const int * descend_connectivity_index = _descending->getIndex();
1156 int * descend_connectivity_index = new int[_count[_numberOfTypes]];
1157 descend_connectivity_index[0]=1;
1158 medGeometryElement* ConstituentsTypes = new medGeometryElement[2];
1159 ConstituentsTypes[0]=MED_NONE;
1160 ConstituentsTypes[1]=MED_NONE;
1161 int * NumberOfConstituentsForeachType = new int[2];
1162 NumberOfConstituentsForeachType[0]=0;
1163 NumberOfConstituentsForeachType[1]=0;
1164 for(int i=0; i<_numberOfTypes; i++)
1166 // initialize descend_connectivity_index array :
1167 int NumberOfConstituents = _type[i].getNumberOfConstituents(1);
1168 for (int j=_count[i];j<_count[i+1];j++)
1170 descend_connectivity_index[j]=descend_connectivity_index[j-1]+NumberOfConstituents;
1171 // compute number of constituent of all cell for each type
1172 for(int k=1;k<NumberOfConstituents+1;k++)
1174 medGeometryElement MEDType = _type[i].getConstituentType(1,k);
1175 if (ConstituentsTypes[0]==MED_NONE)
1177 ConstituentsTypes[0]=MEDType;
1178 NumberOfConstituentsForeachType[0]++;
1179 } else if (ConstituentsTypes[0]==MEDType)
1180 NumberOfConstituentsForeachType[0]++;
1181 else if (ConstituentsTypes[1]==MED_NONE)
1183 ConstituentsTypes[1]=MEDType;
1184 NumberOfConstituentsForeachType[1]++;
1185 } else if (ConstituentsTypes[1]==MEDType)
1186 NumberOfConstituentsForeachType[1]++;
1188 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<MEDType<<" is different of "<<ConstituentsTypes[0]<<" and "<<ConstituentsTypes[1]));
1193 // we could built _constituent !
1194 int TotalNumberOfConstituents = NumberOfConstituentsForeachType[0]+NumberOfConstituentsForeachType[1];
1195 int TotalNumberOfNodes = NumberOfConstituentsForeachType[0]*(ConstituentsTypes[0]%100)+NumberOfConstituentsForeachType[1]*(ConstituentsTypes[1]%100);
1197 //_constituent->_nodal = new MEDSKYLINEARRAY(TotalNumberOfConstituents,TotalNumberOfNodes);
1199 // we use _constituent->_nodal
1200 //const int * ConstituentNodalConnectivity = _constituent->_nodal->getValue();
1201 //const int * ConstituentNodalConnectivityIndex = _constituent->_nodal->getIndex();
1202 int * ConstituentNodalConnectivity = new int[TotalNumberOfNodes];
1203 int * ConstituentNodalConnectivityIndex = new int[TotalNumberOfConstituents+1];
1204 ConstituentNodalConnectivityIndex[0]=1;
1206 _constituent->_entityDimension = _entityDimension-1;
1207 if (ConstituentsTypes[1]==MED_NONE)
1208 _constituent->_numberOfTypes = 1;
1210 _constituent->_numberOfTypes = 2;
1211 _constituent->_geometricTypes = new medGeometryElement[_constituent->_numberOfTypes];
1212 _constituent->_type = new CELLMODEL[_constituent->_numberOfTypes];
1213 _constituent->_count = new med_int[_constituent->_numberOfTypes+1];
1214 _constituent->_count[0]=1;
1215 int* tmp_NumberOfConstituentsForeachType = new med_int[_constituent->_numberOfTypes+1];
1216 tmp_NumberOfConstituentsForeachType[0]=0; // to count constituent of each type
1217 for (int i=0; i<_constituent->_numberOfTypes;i++) {
1218 _constituent->_geometricTypes[i]=ConstituentsTypes[i];
1219 _constituent->_count[i+1]=_constituent->_count[i]+NumberOfConstituentsForeachType[i];
1220 CELLMODEL Type(ConstituentsTypes[i]);
1221 _constituent->_type[i]=Type;
1222 tmp_NumberOfConstituentsForeachType[i+1]=NumberOfConstituentsForeachType[i];
1223 for (int j=tmp_NumberOfConstituentsForeachType[i]; j<tmp_NumberOfConstituentsForeachType[i+1]+tmp_NumberOfConstituentsForeachType[i]; j++)
1224 ConstituentNodalConnectivityIndex[j+1]=ConstituentNodalConnectivityIndex[j]+(ConstituentsTypes[i]%100);
1226 delete [] ConstituentsTypes;
1227 delete [] NumberOfConstituentsForeachType;
1229 // we need reverse nodal connectivity
1230 if (! _reverseNodalConnectivity)
1231 calculateReverseNodalConnectivity();
1232 const int * ReverseNodalConnectivityIndex = _reverseNodalConnectivity->getIndex();
1233 const int * ReverseNodalConnectivityValue = _reverseNodalConnectivity->getValue();
1235 // array to keep reverse descending connectivity
1236 int * ReverseDescendingConnectivityValue = new int[TotalNumberOfConstituents*2];
1238 int TotalNumberOfSubCell = 0;
1239 for (int i=0; i<_numberOfTypes; i++)
1240 { // loop on all cell type
1241 CELLMODEL Type = _type[i];
1242 // int NumberOfNodesPerCell = Type.getNumberOfNodes();
1243 int NumberOfConstituentPerCell = Type.getNumberOfConstituents(1);
1244 for (int j=_count[i]; j<_count[i+1]; j++) // we loop on all cell of this type
1245 for (int k=1; k<=NumberOfConstituentPerCell; k++)
1246 { // we loop on all sub cell of it
1247 if (descend_connectivity[descend_connectivity_index[j-1]+k-2]==0) {
1248 // it is a new sub cell !
1249 // TotalNumberOfSubCell++;
1251 if (Type.getConstituentType(1,k)==_constituent->_geometricTypes[0]){
1252 tmp_NumberOfConstituentsForeachType[0]++;
1253 TotalNumberOfSubCell=tmp_NumberOfConstituentsForeachType[0];
1255 tmp_NumberOfConstituentsForeachType[1]++;
1256 TotalNumberOfSubCell=tmp_NumberOfConstituentsForeachType[1];
1258 //we have maximum two types
1260 descend_connectivity[descend_connectivity_index[j-1]+k-2]=TotalNumberOfSubCell;
1261 ReverseDescendingConnectivityValue[(TotalNumberOfSubCell-1)*2]=j;
1262 int NumberOfNodesPerConstituent = Type.getConstituentType(1,k)%100;
1264 int * NodesLists = new int[NumberOfNodesPerConstituent];
1265 for (int l=0; l<NumberOfNodesPerConstituent; l++) {
1266 NodesLists[l]=_nodal->getIJ(j,Type.getNodeConstituent(1,k,l+1));
1267 ConstituentNodalConnectivity[ConstituentNodalConnectivityIndex[TotalNumberOfSubCell-1]-1+l]=NodesLists[l];
1269 // we use reverse_nodal_connectivity to find the other element which contain this sub cell
1271 // all elements which contains first node of sub cell :
1272 int ReverseNodalConnectivityIndex_0 = ReverseNodalConnectivityIndex[NodesLists[0]-1];
1273 int ReverseNodalConnectivityIndex_1 = getIndexOfEndClassicElementInReverseNodal(ReverseNodalConnectivityValue,ReverseNodalConnectivityIndex,NodesLists[0]-1);
1274 //ReverseNodalConnectivityIndex[NodesLists[0]];
1275 int NumberOfCellsInList = ReverseNodalConnectivityIndex_1-ReverseNodalConnectivityIndex_0;
1277 if (NumberOfCellsInList > 0)
1278 { // we could have no element !
1279 int * CellsList = new int[NumberOfCellsInList];
1280 for (int l=ReverseNodalConnectivityIndex_0; l<ReverseNodalConnectivityIndex_1; l++)
1281 CellsList[l-ReverseNodalConnectivityIndex_0]=ReverseNodalConnectivityValue[l-1];
1283 // foreach node in sub cell, we search elements which are in common
1284 // at the end, we must have only one !
1286 for (int l=1; l<NumberOfNodesPerConstituent; l++) {
1287 int NewNumberOfCellsInList = 0;
1288 int * NewCellsList = new int[NumberOfCellsInList];
1289 for (int l1=0; l1<NumberOfCellsInList; l1++)
1290 for (int l2=ReverseNodalConnectivityIndex[NodesLists[l]-1]; l2<getIndexOfEndClassicElementInReverseNodal(ReverseNodalConnectivityValue,ReverseNodalConnectivityIndex,NodesLists[l]-1); l2++)
1291 //for (int l2=ReverseNodalConnectivityIndex[NodesLists[l]-1]; l2<ReverseNodalConnectivityIndex[NodesLists[l]]; l2++)
1293 if (CellsList[l1]<ReverseNodalConnectivityValue[l2-1])
1294 // increasing order : CellsList[l1] are not in elements list of node l
1296 if ((CellsList[l1]==ReverseNodalConnectivityValue[l2-1])&(CellsList[l1]!=j)) {
1297 // we have found one
1298 NewCellsList[NewNumberOfCellsInList]=CellsList[l1];
1299 NewNumberOfCellsInList++;
1303 NumberOfCellsInList = NewNumberOfCellsInList;
1305 delete [] CellsList;
1306 CellsList = NewCellsList;
1309 if (NumberOfCellsInList > 0) { // We have found some elements !
1310 int CellNumber = CellsList[0];
1311 //delete [] CellsList;
1312 if (NumberOfCellsInList>1)
1313 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"More than one other Cell ("<<NumberOfCellsInList<<") !"));
1315 if (NumberOfCellsInList==1)
1317 ReverseDescendingConnectivityValue[(TotalNumberOfSubCell-1)*2+1]=CellNumber;
1319 // we search sub cell number in this cell to not calculate it another time
1322 for (int l=0; l<_numberOfTypes; l++)
1323 if (CellNumber < _count[l+1]) {
1327 //int sub_cell_count2 = Type2.get_entities_count(1);
1328 //int nodes_cell_count2 = Type2.get_nodes_count();
1330 for (int l=1; l<=Type2.getNumberOfConstituents(1);l++) { // on all sub cell
1332 for (int m=1; m<=Type2.getConstituentType(1,l)%100; m++)
1333 for (int n=1; n<=Type.getConstituentType(1,k)%100; n++)
1335 if (_nodal->getIJ(CellNumber,Type2.getNodeConstituent(1,l,m)) == NodesLists[n-1])
1338 if (counter==Type.getConstituentType(1,k)%100)
1340 descend_connectivity[descend_connectivity_index[CellNumber-1]+l-2]=-1*TotalNumberOfSubCell; // because, see it in other side !
1347 MESSAGE(LOC<<"ERROR ERROR ERROR ERROR ERROR : we find any subcell !!!"); // exception ?
1350 ReverseDescendingConnectivityValue[(TotalNumberOfSubCell-1)*2+1]=0;
1352 delete [] CellsList;
1355 delete [] NodesLists;
1359 // we adjust _constituent
1360 int NumberOfConstituent=0;
1361 int SizeOfConstituentNodal=0;
1362 for (int i=0;i<_constituent->_numberOfTypes; i++) {
1363 NumberOfConstituent += tmp_NumberOfConstituentsForeachType[i]-_constituent->_count[i]+1;
1364 SizeOfConstituentNodal += (tmp_NumberOfConstituentsForeachType[i]-_constituent->_count[i]+1)*_constituent->_type[i].getNumberOfNodes();
1366 // we built new _nodal attribute in _constituent
1367 //MEDSKYLINEARRAY *ConstituentNodal = new MEDSKYLINEARRAY(NumberOfConstituent,SizeOfConstituentNodal);
1368 //const int *ConstituentNodalValue = ConstituentNodal->getValue();
1369 //const int *ConstituentNodalIndex = ConstituentNodal->getIndex();
1370 int *ConstituentNodalValue = new int[SizeOfConstituentNodal];
1371 int *ConstituentNodalIndex = new int[NumberOfConstituent+1];
1372 ConstituentNodalIndex[0]=1;
1373 // we build _reverseDescendingConnectivity
1374 //_reverseDescendingConnectivity = new MEDSKYLINEARRAY(NumberOfConstituent,2*NumberOfConstituent);
1375 //const int *reverseDescendingConnectivityValue = _reverseDescendingConnectivity->getValue();
1376 //const int *reverseDescendingConnectivityIndex = _reverseDescendingConnectivity->getIndex();
1377 int *reverseDescendingConnectivityValue = new int[2*NumberOfConstituent];
1378 int *reverseDescendingConnectivityIndex = new int[NumberOfConstituent+1];
1379 reverseDescendingConnectivityIndex[0]=1;
1381 // first constituent type
1382 for(int j=0; j<tmp_NumberOfConstituentsForeachType[0]; j++) {
1383 ConstituentNodalIndex[j+1]=ConstituentNodalConnectivityIndex[j+1];
1384 for(int k=ConstituentNodalIndex[j]; k<ConstituentNodalIndex[j+1]; k++){
1385 ConstituentNodalValue[k-1]=ConstituentNodalConnectivity[k-1];
1387 reverseDescendingConnectivityIndex[j+1]=reverseDescendingConnectivityIndex[j]+2;
1388 for(int k=reverseDescendingConnectivityIndex[j]; k<reverseDescendingConnectivityIndex[j+1]; k++){
1389 reverseDescendingConnectivityValue[k-1]=ReverseDescendingConnectivityValue[k-1];
1392 // second type if any
1393 if (_constituent->_numberOfTypes==2) {
1394 int offset = _constituent->_count[1]-tmp_NumberOfConstituentsForeachType[0]-1;
1395 int offset1=offset*_constituent->_type[0].getNumberOfNodes();
1396 int offset2=offset*2;
1397 int NumberOfNodesPerConstituent = _constituent->_type[1].getNumberOfNodes();
1398 for(int j=tmp_NumberOfConstituentsForeachType[0]; j<NumberOfConstituent; j++) {
1399 ConstituentNodalIndex[j+1]=ConstituentNodalIndex[j]+NumberOfNodesPerConstituent;
1400 for(int k=ConstituentNodalIndex[j]; k<ConstituentNodalIndex[j+1]; k++){
1401 ConstituentNodalValue[k-1]=ConstituentNodalConnectivity[offset1+k-1];
1403 reverseDescendingConnectivityIndex[j+1]=reverseDescendingConnectivityIndex[j]+2;
1404 for(int k=reverseDescendingConnectivityIndex[j]; k<reverseDescendingConnectivityIndex[j+1]; k++){
1405 reverseDescendingConnectivityValue[k-1]=ReverseDescendingConnectivityValue[offset2+k-1];
1408 _constituent->_count[2]=NumberOfConstituent+1;
1409 // we correct _descending to adjust face number
1410 for(int j=0;j<DescendingSize;j++)
1411 if (descend_connectivity[j]>tmp_NumberOfConstituentsForeachType[0])
1412 descend_connectivity[j]-=offset;
1415 delete [] ConstituentNodalConnectivityIndex;
1416 delete [] ConstituentNodalConnectivity;
1417 delete [] ReverseDescendingConnectivityValue;
1418 _constituent->_count[1]=tmp_NumberOfConstituentsForeachType[0]+1;
1419 delete [] tmp_NumberOfConstituentsForeachType;
1421 _descending = new MEDSKYLINEARRAY(_count[_numberOfTypes]-1,
1423 descend_connectivity_index,
1424 descend_connectivity);
1425 delete [] descend_connectivity_index;
1426 delete [] descend_connectivity;
1429 vector<int> PolyDescending;
1430 vector<int> Reversedescendingconnectivityvalue(reverseDescendingConnectivityValue,reverseDescendingConnectivityValue + 2*NumberOfConstituent);
1431 vector<int> Reversedescendingconnectivityindex(reverseDescendingConnectivityIndex,reverseDescendingConnectivityIndex + NumberOfConstituent);
1432 delete [] reverseDescendingConnectivityValue;
1433 delete [] reverseDescendingConnectivityIndex;
1436 // polygons (2D mesh)
1438 vector<int> Constituentnodalvalue(ConstituentNodalValue,ConstituentNodalValue + SizeOfConstituentNodal);
1439 vector<int> Constituentnodalindex(ConstituentNodalIndex,ConstituentNodalIndex + NumberOfConstituent+1);
1440 delete [] ConstituentNodalValue;
1441 delete [] ConstituentNodalIndex;
1442 int NumberOfNewSeg = 0;
1444 for (int i=0; i <getNumberOfPolygons(); i++) // for each polygon
1446 const int * vector_begin = &_polygonsNodal->getValue()[_polygonsNodal->getIndex()[i]-1];
1447 int vector_size = _polygonsNodal->getIndex()[i+1]-_polygonsNodal->getIndex()[i]+1;
1448 vector<int> myPolygon(vector_begin, vector_begin+vector_size);
1449 myPolygon[myPolygon.size()-1] = myPolygon[0]; // because first and last point make a segment
1451 for (int j=0; j<myPolygon.size()-1; j++) // for each segment of polygon
1453 MEDMODULUSARRAY segment_poly(2,&myPolygon[j]);
1454 int ret_compare = 0;
1456 // we search it in existing segments
1458 for (int k=0; k<Constituentnodalindex.size()-1; k++)
1460 MEDMODULUSARRAY segment(2,&Constituentnodalvalue[0] + Constituentnodalindex[k]-1);
1461 ret_compare = segment_poly.compare(segment);
1462 if (ret_compare) // segment_poly already exists
1464 PolyDescending.push_back(ret_compare*(k+1)); // we had it to the connectivity
1465 insert_vector(Reversedescendingconnectivityvalue, 2*k+1, i+1 + getNumberOf(MED_CELL,MED_ALL_ELEMENTS)); // add polygon i to reverse descending connectivity for segment_poly (in 2sd place)
1470 // segment_poly must be created
1475 PolyDescending.push_back(NumberOfConstituent+NumberOfNewSeg); // we had it to the connectivity
1476 Constituentnodalvalue.push_back(segment_poly[0]);
1477 Constituentnodalvalue.push_back(segment_poly[1]);
1478 Constituentnodalindex.push_back(Constituentnodalindex[Constituentnodalindex.size()-1] + 2); // we have only segments
1479 insert_vector(Reversedescendingconnectivityvalue, 2*(NumberOfConstituent+NumberOfNewSeg-1), i+1 + getNumberOf(MED_CELL,MED_ALL_ELEMENTS)); // add polygon i to reverse descending connectivity for segment_poly (in 1st place)
1480 insert_vector(Reversedescendingconnectivityindex, NumberOfConstituent+NumberOfNewSeg-1, 2*(NumberOfConstituent+NumberOfNewSeg-1)+1); // idem with index
1485 if (getNumberOfPolygons() > 0)
1487 _polygonsDescending = new MEDSKYLINEARRAY(getNumberOfPolygons(),_polygonsNodal->getLength(),_polygonsNodal->getIndex(),&PolyDescending[0]); // index are the same for polygons nodal and descending connectivities
1489 NumberOfConstituent += NumberOfNewSeg;
1490 SizeOfConstituentNodal += 2*NumberOfNewSeg;
1491 _constituent->_count[1] = NumberOfConstituent+1;
1495 // polyhedron (3D mesh)
1497 vector<int> Constituentpolygonsnodalvalue;
1498 vector<int> Constituentpolygonsnodalindex(1,1);
1499 int NumberOfNewFaces = 0; // by convention new faces are polygons
1501 for (int i=0; i<getNumberOfPolyhedron(); i++) // for each polyhedron
1503 // we create a POLYHEDRONARRAY containing only this polyhedra
1504 int myNumberOfFaces = _polyhedronNodal->getPolyhedronIndex()[i+1]-_polyhedronNodal->getPolyhedronIndex()[i];
1505 int myNumberOfNodes = _polyhedronNodal->getFacesIndex()[_polyhedronNodal->getPolyhedronIndex()[i+1]-1]-_polyhedronNodal->getFacesIndex()[_polyhedronNodal->getPolyhedronIndex()[i]-1];
1506 POLYHEDRONARRAY myPolyhedra(1,myNumberOfFaces,myNumberOfNodes);
1507 vector<int> myFacesIndex(_polyhedronNodal->getFacesIndex() + _polyhedronNodal->getPolyhedronIndex()[i]-1, _polyhedronNodal->getFacesIndex() + _polyhedronNodal->getPolyhedronIndex()[i]-1 + myNumberOfFaces+1);
1508 for (int j=myNumberOfFaces; j>=0; j--)
1509 myFacesIndex[j] -= myFacesIndex[0]-1;
1510 myPolyhedra.setFacesIndex(&myFacesIndex[0]);
1511 myPolyhedra.setNodes(_polyhedronNodal->getNodes() + _polyhedronNodal->getFacesIndex()[_polyhedronNodal->getPolyhedronIndex()[i]-1]-1);
1513 for (int j=0; j<myPolyhedra.getNumberOfFaces(); j++) // for each face of polyhedra
1515 int myFaceNumberOfNodes = myPolyhedra.getFacesIndex()[j+1]-myPolyhedra.getFacesIndex()[j];
1516 MEDMODULUSARRAY face_poly(myFaceNumberOfNodes,myPolyhedra.getNodes() + myPolyhedra.getFacesIndex()[j]-1);
1517 int ret_compare = 0;
1519 // we search it in existing faces
1521 // we search first in TRIA3 and QUAD4
1522 if (myFaceNumberOfNodes == 3 || myFaceNumberOfNodes == 4)
1524 int Begin = -1; // first TRIA3 or QUAD4
1525 int Number = 0; // number of TRIA3 or QUAD4
1526 for (int k=0; k<Constituentnodalindex.size()-1; k++)
1527 if (Constituentnodalindex[k+1]-Constituentnodalindex[k] == myFaceNumberOfNodes)
1534 for (int k=0; k<Number; k++)
1536 MEDMODULUSARRAY face(myFaceNumberOfNodes,&Constituentnodalvalue[0] + Constituentnodalindex[Begin+k]-1);
1537 ret_compare = face_poly.compare(face);
1540 PolyDescending.push_back(ret_compare*(Begin+k+1)); // we had it to the connectivity
1541 insert_vector(Reversedescendingconnectivityvalue, 2*(Begin+k)+1, i+1 + getNumberOf(MED_CELL,MED_ALL_ELEMENTS)); // add polyhedra i to reverse descending connectivity for face_poly (in 2sd place)
1547 // we search last in POLYGONS
1550 for (int k=0; k<static_cast<int>(Constituentpolygonsnodalindex.size())-1; k++) // we must cast the unsigned int into int before doing -1
1552 if (Constituentpolygonsnodalindex[k+1]-Constituentpolygonsnodalindex[k] == myFaceNumberOfNodes)
1554 MEDMODULUSARRAY face(myFaceNumberOfNodes,&Constituentpolygonsnodalvalue[0] + Constituentpolygonsnodalindex[k]-1);
1555 ret_compare = face_poly.compare(face);
1558 PolyDescending.push_back(ret_compare*(NumberOfConstituent+k+1)); // we had it to the connectivity
1559 insert_vector(Reversedescendingconnectivityvalue, 2*(NumberOfConstituent+k)+1, i+1 + getNumberOf(MED_CELL,MED_ALL_ELEMENTS)); // add polyhedra i to reverse descending connectivity for face_poly (in 2sd place)
1566 // if not found, face_poly must be created
1571 PolyDescending.push_back(NumberOfConstituent+NumberOfNewFaces); // we had it to the connectivity
1572 for (int k=0; k<myFaceNumberOfNodes; k++)
1573 Constituentpolygonsnodalvalue.push_back(face_poly[k]);
1574 Constituentpolygonsnodalindex.push_back(Constituentpolygonsnodalindex[Constituentpolygonsnodalindex.size()-1]+myFaceNumberOfNodes);
1575 insert_vector(Reversedescendingconnectivityvalue, 2*(NumberOfConstituent+NumberOfNewFaces-1), i+1 + getNumberOf(MED_CELL,MED_ALL_ELEMENTS)); // add polyhedra i to reverse descending connectivity for face_poly (in 1st place)
1576 insert_vector(Reversedescendingconnectivityindex, NumberOfConstituent+NumberOfNewFaces-1, 2*(NumberOfConstituent+NumberOfNewFaces-1)+1); // idem with index
1581 if (getNumberOfPolyhedron() > 0)
1583 _polyhedronDescending = new MEDSKYLINEARRAY(getNumberOfPolyhedron(),_polyhedronNodal->getNumberOfFaces(),_polyhedronNodal->getPolyhedronIndex(),&PolyDescending[0]); // polyhedron index are the same for nodal and descending connectivities
1585 if (_constituent->_polygonsNodal != NULL)
1586 delete [] _constituent->_polygonsNodal;
1587 _constituent->_polygonsNodal = new MEDSKYLINEARRAY(Constituentpolygonsnodalindex.size()-1,Constituentpolygonsnodalvalue.size(),&Constituentpolygonsnodalindex[0],&Constituentpolygonsnodalvalue[0]);
1590 // delete _constituent->_nodal
1591 _constituent->_nodal = new MEDSKYLINEARRAY(NumberOfConstituent,
1592 SizeOfConstituentNodal,
1593 &Constituentnodalindex[0],
1594 &Constituentnodalvalue[0]);
1596 int NumberOfConstituentWithPolygons = NumberOfConstituent + NumberOfNewFaces;
1597 Reversedescendingconnectivityindex.push_back(Reversedescendingconnectivityindex[Reversedescendingconnectivityindex.size()-1]+2); // we complete the index
1598 _reverseDescendingConnectivity = new MEDSKYLINEARRAY(NumberOfConstituentWithPolygons+1,
1599 2*NumberOfConstituentWithPolygons,
1600 &Reversedescendingconnectivityindex[0],
1601 &Reversedescendingconnectivityvalue[0]);
1607 /*! Not implemented yet */
1608 //--------------------------------------------------------------------//
1609 void CONNECTIVITY::calculateNeighbourhood(CONNECTIVITY &myConnectivity)
1610 //--------------------------------------------------------------------//
1612 const char * LOC = "CONNECTIVITY::calculateNeighbourhood(CONNECTIVITY &myConnectivity) : ";
1615 MESSAGE(LOC<<"method not yet implemented " << myConnectivity._entity);
1624 Returns the geometry of an element given by its entity type & its global number.
1626 Example : medGeometryElement myType = myConnectivity.getElementType(MED_CELL,35);
1628 //--------------------------------------------------------------------//
1629 medGeometryElement CONNECTIVITY::getElementType(medEntityMesh Entity,int globalNumber) const
1630 //--------------------------------------------------------------------//
1632 const char * LOC = "medGeometryElement CONNECTIVITY::getElementType(medEntityMesh Entity, int globalNumber) const : ";
1634 int globalNumberMin = 1;
1635 int globalNumberMax ;
1637 if (_entity==Entity) globalNumberMax = _count[_numberOfTypes];
1638 else if (_constituent!=NULL) globalNumberMax = _constituent->_count[_constituent->_numberOfTypes];
1640 throw MEDEXCEPTION(LOCALIZED("getElementType : Entity not defined !"));
1642 // The globalNumber must verify : 1 <= globalNumber < _count[_numberOfTypes] (== totalNumberOfElement+1)
1644 if ( (globalNumber < globalNumberMin) || (globalNumber > globalNumberMax-1 ) )
1645 throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "globalNumber must be between >= |"
1646 << globalNumberMin <<"| and <= |" << globalNumberMax-1 ));
1648 if (_entity==Entity) {
1649 for(int i=1; i<=_numberOfTypes;i++)
1650 if (globalNumber<_count[i])
1651 return _geometricTypes[i-1];
1653 else if (_constituent!=NULL)
1654 return _constituent->getElementType(Entity,globalNumber);
1656 throw MEDEXCEPTION(LOCALIZED("getElementType : Entity not defined !"));
1657 throw MEDEXCEPTION(LOCALIZED("getElementType : Wrong Number !"));
1663 Idem as getElementType method except that it manages polygon and polyhedron.
1665 //--------------------------------------------------------------------//
1666 medGeometryElement CONNECTIVITY::getElementTypeWithPoly(medEntityMesh Entity,int globalNumber) const
1667 //--------------------------------------------------------------------//
1669 int globalNumberMaxOfClassicType;
1672 globalNumberMaxOfClassicType=_count[_numberOfTypes];
1675 if(globalNumber<globalNumberMaxOfClassicType)
1677 for(int i=1; i<=_numberOfTypes;i++)
1678 if (globalNumber<_count[i])
1679 return _geometricTypes[i-1];
1680 throw MEDEXCEPTION("never happens just for compilo");
1684 int localNumberInPolyArray=globalNumber-globalNumberMaxOfClassicType;
1685 int nbOfPol=getNumberOfElementOfPolyType(_entity);
1686 if(localNumberInPolyArray<nbOfPol)
1687 return getPolyTypeRelativeTo();
1689 throw MEDEXCEPTION("getElementTypeWithPoly : unexisting element type");
1693 throw MEDEXCEPTION("getElementTypeWithPoly : globalNumber < 1");
1697 if(_constituent!=NULL)
1698 return _constituent->getElementTypeWithPoly(Entity,globalNumber);
1700 throw MEDEXCEPTION("getElementTypeWithPoly : Entity not defined !");
1704 ostream & MEDMEM::operator<<(ostream &os, CONNECTIVITY &co)
1706 os << endl << "------------- Entity = ";
1721 case MED_ALL_ENTITIES:
1722 os << "MED_ALL_ENTITIES";
1728 os << " -------------\n\nMedConnectivity : ";
1729 switch (co._typeConnectivity)
1732 os << "MED_NODAL\n";
1734 case MED_DESCENDING:
1735 os << "MED_DESCENDING\n";
1740 os << "Entity dimension : " << co._entityDimension << endl;
1741 os << "Number of nodes : " << co._numberOfNodes << endl;
1742 os << "Number of types : " << co._numberOfTypes << endl;
1743 for (int i=0; i!=co._numberOfTypes ; ++i)
1744 os << " -> Type " << co._geometricTypes[i] << " (" << co._type[i].getName() << ") : "
1745 << co._count[i+1]-co._count[i] << " elements" << endl;
1747 if (co._typeConnectivity == MED_NODAL )
1749 for (int i=0; i<co._numberOfTypes; i++)
1751 os << endl << co._type[i].getName() << " : " << endl;
1752 int numberofelements = co._count[i+1]-co._count[i];
1753 const int * connectivity = co.getConnectivity(co._typeConnectivity, co._entity, co._geometricTypes[i]);
1754 int numberofnodespercell = co._geometricTypes[i]%100;
1755 for (int j=0;j<numberofelements;j++)
1757 os << setw(6) << j+1 << " : ";
1758 for (int k=0;k<numberofnodespercell;k++)
1759 os << connectivity[j*numberofnodespercell+k]<<" ";
1764 os << endl << "MED_POLYGON : " << co.getNumberOfPolygons() << " polygons" << endl;
1765 if (co.getNumberOfPolygons() > 0)
1767 const int* polygonsconnectivity = co.getPolygonsConnectivity(MED_NODAL,co.getEntity());
1768 const int* polygonsconnectivityindex = co.getPolygonsConnectivityIndex(MED_NODAL,co.getEntity());
1770 for (int i=0; i<co.getNumberOfPolygons(); i++)
1772 int numberofnodesperpolygon = polygonsconnectivityindex[i+1]-polygonsconnectivityindex[i];
1773 for (int j=0; j<numberofnodesperpolygon; j++)
1774 os << polygonsconnectivity[polygonsconnectivityindex[i]-1+j] << " ";
1779 os << endl << "MED_POLYHEDRA : " << co.getNumberOfPolyhedron() << " polyhedron" << endl;
1780 if (co.getNumberOfPolyhedron() > 0)
1782 const int* polyhedronconnectivity = co.getPolyhedronConnectivity(MED_NODAL);
1783 const int* polyhedronfacesindex = co.getPolyhedronFacesIndex();
1784 const int* polyhedronindex = co.getPolyhedronIndex(MED_NODAL);
1786 for (int i=0; i<co.getNumberOfPolyhedron(); i++)
1788 int numberoffacesperpolyhedra = polyhedronindex[i+1]-polyhedronindex[i];
1789 for (int j=0; j<numberoffacesperpolyhedra; j++)
1791 int numberofnodesperface = polyhedronfacesindex[polyhedronindex[i]-1+j+1]-polyhedronfacesindex[polyhedronindex[i]-1+j];
1792 for (int k=0; k<numberofnodesperface; k++)
1793 os << polyhedronconnectivity[polyhedronfacesindex[polyhedronindex[i]-1+j]-1+k] << " ";
1794 if (j != numberoffacesperpolyhedra-1) os << "| ";
1801 else if (co._typeConnectivity == MED_DESCENDING)
1803 int numberofelements = co.getNumberOf( co._entity , MED_ALL_ELEMENTS);
1804 if (numberofelements > 0)
1806 const int *connectivity = co.getConnectivity( co._typeConnectivity, co._entity, MED_ALL_ELEMENTS);
1807 const int *connectivity_index = co.getConnectivityIndex( co._typeConnectivity, co._entity );
1809 for ( int j=0; j!=numberofelements; j++ )
1811 os << "Element " << j+1 << " : ";
1812 for ( int k=connectivity_index[j]; k!=connectivity_index[j+1]; k++ )
1813 os << connectivity[k-1] << " ";
1818 os << endl << "MED_POLYGON : " << co.getNumberOfPolygons() << " polygons" << endl;
1819 if (co.getNumberOfPolygons() > 0)
1821 const int* polygonsconnectivity = co.getPolygonsConnectivity(MED_DESCENDING,co.getEntity());
1822 const int* polygonsconnectivityindex = co.getPolygonsConnectivityIndex(MED_DESCENDING,co.getEntity());
1824 for (int i=0; i<co.getNumberOfPolygons(); i++)
1826 int numberofedgesperpolygon = polygonsconnectivityindex[i+1]-polygonsconnectivityindex[i];
1827 for (int j=0; j<numberofedgesperpolygon; j++)
1828 os << polygonsconnectivity[polygonsconnectivityindex[i]-1+j] << " ";
1833 os << endl << "MED_POLYHEDRA : " << co.getNumberOfPolyhedron() << " polyhedron" << endl;
1834 if (co.getNumberOfPolyhedron() > 0)
1836 const int* polyhedronconnectivity = co.getPolyhedronConnectivity(MED_DESCENDING);
1837 const int* polyhedronindex = co.getPolyhedronIndex(MED_DESCENDING);
1839 for (int i=0; i<co.getNumberOfPolyhedron(); i++)
1841 int numberoffacesperpolyhedra = polyhedronindex[i+1]-polyhedronindex[i];
1842 for (int j=0; j<numberoffacesperpolyhedra; j++)
1843 os << polyhedronconnectivity[polyhedronindex[i]-1+j] << " ";
1850 if (co._constituent)
1851 os << endl << *co._constituent << endl;
1857 method that adds to vector 'nodes' all the nodes of polyhedron with id 'polyhedronId'.
1858 WARNING the returned pointer should be deallocated. Returned nodes and polyhedronId are in form [1,...]
1860 int *CONNECTIVITY::getNodesOfPolyhedron(int polyhedronId, int& lgthOfTab) const
1862 const int *nodes=getPolyhedronConnectivity(MED_EN::MED_NODAL);
1863 const int *faces=getPolyhedronFacesIndex();
1864 const int *glob=getPolyhedronIndex(MED_EN::MED_NODAL);
1865 int offsetWithClassicType=getNumberOf(_entity,MED_ALL_ELEMENTS);
1867 int startFace=glob[polyhedronId-offsetWithClassicType-1]-1;
1868 int endFace=glob[polyhedronId-offsetWithClassicType]-1;
1870 for(i=startFace;i!=endFace;i++)
1872 for(int j=faces[i];j<faces[i+1];j++)
1873 retInSet.insert(nodes[j-1]);
1875 lgthOfTab=retInSet.size();
1876 int *ret=new int[lgthOfTab];
1877 set<int>::iterator iter=retInSet.begin();
1879 for(iter=retInSet.begin();iter!=retInSet.end();iter++)
1885 Idem as MESH::getNodesOfPolyhedron except that returned nodes are sorted by face. 'nbOfNodesPerFaces' is an array of size 'nbOfFaces'.
1886 Returned int** has a size of 'nbOfNodesPerFaces' too, and for each element j in int** the size is nbOfNodesPerFaces[j].
1887 Warning both returned 'nbOfNodesPerFaces' and returned value should be deallocated. Returned nodes and 'polyhedronId' are in form [1,...]
1889 int **CONNECTIVITY::getNodesPerFaceOfPolyhedron(int polyhedronId, int& nbOfFaces, int* & nbOfNodesPerFaces) const
1892 const int *nodes=getPolyhedronConnectivity(MED_EN::MED_NODAL);
1893 const int *faces=getPolyhedronFacesIndex();
1894 const int *glob=getPolyhedronIndex(MED_EN::MED_NODAL);
1895 int offsetWithClassicType=getNumberOf(_entity,MED_ALL_ELEMENTS);
1897 int startFace=glob[polyhedronId-offsetWithClassicType-1]-1;
1898 nbOfFaces=glob[polyhedronId-offsetWithClassicType]-startFace-1;
1899 nbOfNodesPerFaces=new int[nbOfFaces];
1900 int **ret=new int* [nbOfFaces];
1901 for(i=0;i<nbOfFaces;i++)
1903 int startNode=faces[startFace+i]-1;
1904 nbOfNodesPerFaces[i]=faces[startFace+i+1]-startNode-1;
1905 ret[i]=(int *)(nodes)+startNode;
1910 int CONNECTIVITY::getNumberOfTypesWithPoly(MED_EN::medEntityMesh Entity) const
1912 if (_entity==Entity)
1913 return _numberOfTypes+getNumberOfPolyType();
1914 else if (_constituent!=NULL)
1915 return _constituent->getNumberOfTypesWithPoly(Entity);
1920 int CONNECTIVITY::getNumberOfPolyType() const
1922 if (_entity==MED_CELL && _entityDimension==3)
1924 if(getNumberOfPolyhedron()>0)
1927 else if ((_entity==MED_CELL && _entityDimension==2) || (_entity==MED_FACE && _entityDimension==2))
1928 if(getNumberOfPolygons()>0)
1933 int CONNECTIVITY::getNumberOfElementOfPolyType(MED_EN::medEntityMesh Entity) const
1937 if (_entity==MED_CELL && _entityDimension==3)
1938 return getNumberOfPolyhedron();
1939 else if ((_entity==MED_CELL && _entityDimension==2) || (_entity==MED_FACE && _entityDimension==2))
1940 return getNumberOfPolygons();
1945 if(_constituent!=NULL)
1946 return _constituent->getNumberOfElementOfPolyType(Entity);
1948 throw MEDEXCEPTION("_constituent required to evaluate getNumberOfElementOfPolyType");
1953 Method equivalent to getGeometricTypes except that it includes not only classical Types but polygons/polyhedra also.
1954 WARNING the returned array MUST be deallocated.
1956 MED_EN::medGeometryElement * CONNECTIVITY::getGeometricTypesWithPoly(MED_EN::medEntityMesh Entity) const throw (MEDEXCEPTION)
1960 int nbOfTypesTotal=getNumberOfTypesWithPoly(Entity);
1961 int nbOfTypesWithoutPoly=getNumberOfTypes(Entity);
1962 medGeometryElement *ret=new medGeometryElement[nbOfTypesTotal];
1963 memcpy(ret,getGeometricTypes(Entity),nbOfTypesWithoutPoly*sizeof(medGeometryElement));
1964 if(nbOfTypesTotal>nbOfTypesWithoutPoly)
1966 if (Entity==MED_CELL && _entityDimension==3)
1967 ret[nbOfTypesWithoutPoly]=MED_EN::MED_POLYHEDRA;
1968 else if((Entity==MED_CELL && _entityDimension==2) || (Entity==MED_FACE && _entityDimension==2))
1969 ret[nbOfTypesWithoutPoly]=MED_EN::MED_POLYGON;
1976 return _constituent->getGeometricTypesWithPoly(Entity);
1977 throw MEDEXCEPTION("_constituent required to evaluate getGeometricTypesWithPoly");
1982 Method used in CalculateDescendingConnectivity. So it's typically a private method.
1983 The aim of this method is to hide to CalculateDescendingConnectivity algorithm the fact that in reverse connectivity polygons and polyhedrons
1984 are not in separate data structure, contrary to not reverse connectivities.
1986 int CONNECTIVITY::getIndexOfEndClassicElementInReverseNodal(const int *reverseNodalValue, const int *reverseNodalIndex, int rk) const
1988 int nbOfLastElt=getNumberOf(_entity,MED_ALL_ELEMENTS);
1989 int ret=reverseNodalIndex[rk];
1990 for(int i=reverseNodalIndex[rk];i<reverseNodalIndex[rk+1];i++)
1992 if(reverseNodalValue[i-1]<=nbOfLastElt)
1999 Method that inverts the face with faceId 'faceId' in the data structure. If it's a polygon face 'faceId' is a value between 1 and nbOfPolygons.
2000 This method has to be applied on a instance of MEDMEM::CONNECTIVITY with _entityDimension==3.
2001 WARNING calculateDescendingConnectivity must have been called before.
2003 void CONNECTIVITY::invertConnectivityForAFace(int faceId, const int *nodalConnForFace, bool polygonFace)
2005 // we have always 2 neighbourings
2006 int cell1 = _reverseDescendingConnectivity->getIJ(faceId,1);
2007 int cell2 = _reverseDescendingConnectivity->getIJ(faceId,2);
2009 if (cell2 != 0) { // we are not on border, update compulsary. If we aren't on border no update necessary so WARNING because user specified a bad oriented face
2010 // Updating _reverseDescendingConnectivity
2011 _reverseDescendingConnectivity->setIJ(faceId,1,cell2);
2012 _reverseDescendingConnectivity->setIJ(faceId,2,cell1);
2013 // Updating _constituent->_nodal because of reversity
2014 MEDSKYLINEARRAY *currentNodal=(!polygonFace)?_constituent->_nodal:_constituent->_polygonsNodal;
2015 MEDSKYLINEARRAY *currentDescending=(!polygonFace)?_descending:_polygonsDescending;
2016 const int *descendingNodalIndex=(!polygonFace)?_constituent->_nodal->getIndex():_constituent->_polygonsNodal->getIndex();
2017 const int *newDescendingIndex=(!polygonFace)?_descending->getIndex():_polygonsDescending->getIndex();
2018 for(int iarray=1;iarray<=(descendingNodalIndex[faceId]-descendingNodalIndex[faceId-1]);iarray++)
2019 currentNodal->setIJ(faceId,iarray,nodalConnForFace[iarray-1]);
2021 // Updating _descending for cell1 and cell2
2022 for(int iface=newDescendingIndex[cell1-1];iface<=newDescendingIndex[cell1];iface++)
2023 if (currentDescending->getIndexValue(iface)==faceId)
2024 currentDescending->setIndexValue(iface,-faceId);
2025 else if (currentDescending->getIndexValue(iface)==-faceId)
2026 currentDescending->setIndexValue(iface,faceId);
2028 for(int iface=newDescendingIndex[cell2-1];iface<newDescendingIndex[cell2];iface++)
2029 if (currentDescending->getIndexValue(iface)==faceId)
2030 currentDescending->setIndexValue(iface,-faceId);
2031 else if (_descending->getIndexValue(iface)==-faceId)
2032 currentDescending->setIndexValue(iface,faceId);
2037 Method with 2 output : the connectivity required and its length 'lgth'. This method gives the connectivity independently it is a polygons/polyhedrons or normal
2040 const int * CONNECTIVITY::getConnectivityOfAnElementWithPoly(MED_EN::medConnectivity ConnectivityType, MED_EN::medEntityMesh Entity, int Number, int &lgth)
2042 if(Entity==MED_EN::MED_NODE)
2043 throw MEDEXCEPTION("No connectivity attached to a node entity");
2046 if(_entity==MED_EDGE && _entityDimension==1)
2048 const int * newConstituentValue = 0;
2049 const int * newConstituentIndex = 0;
2050 newConstituentValue = getConnectivity(ConnectivityType,Entity,MED_ALL_ELEMENTS);
2051 newConstituentIndex = getConnectivityIndex(ConnectivityType,Entity);
2052 lgth=newConstituentIndex[Number]-newConstituentIndex[Number-1];
2053 return newConstituentValue+newConstituentIndex[Number-1]-1;
2055 int nbOfClassicalElements=getNumberOf(_entity,MED_ALL_ELEMENTS);
2056 if((_entity==MED_FACE && _entityDimension==2) || (_entity==MED_CELL && _entityDimension==2))
2058 const int * newConstituentValue = 0;
2059 const int * newConstituentIndex = 0;
2060 if(Number<=nbOfClassicalElements)
2062 newConstituentValue = getConnectivity(ConnectivityType,Entity,MED_ALL_ELEMENTS);
2063 newConstituentIndex = getConnectivityIndex(ConnectivityType,Entity);
2064 lgth=newConstituentIndex[Number]-newConstituentIndex[Number-1];
2065 return newConstituentValue+newConstituentIndex[Number-1]-1;
2069 int localNumber=Number-nbOfClassicalElements-1;
2070 if(localNumber<getNumberOfPolygons())
2072 newConstituentValue = getPolygonsConnectivity(ConnectivityType,Entity);
2073 newConstituentIndex = getPolygonsConnectivityIndex(ConnectivityType,Entity);
2074 lgth=newConstituentIndex[localNumber+1]-newConstituentIndex[localNumber];
2075 return newConstituentValue+newConstituentIndex[localNumber]-1;
2078 throw MEDEXCEPTION("Unknown number");
2081 else if(_entity==MED_CELL && _entityDimension==3)
2083 const int * newConstituentValue = 0;
2084 const int * newConstituentIndex = 0;
2085 if(Number<=nbOfClassicalElements)
2087 newConstituentValue = getConnectivity(ConnectivityType,Entity,MED_ALL_ELEMENTS);
2088 newConstituentIndex = getConnectivityIndex(ConnectivityType,Entity);
2089 lgth=newConstituentIndex[Number]-newConstituentIndex[Number-1];
2090 return newConstituentValue+newConstituentIndex[Number-1]-1;
2094 int localNumber=Number-nbOfClassicalElements-1;
2095 if(localNumber<getNumberOfPolyhedron())
2097 if(ConnectivityType==MED_NODAL)
2098 throw MEDEXCEPTION("NODAL Connectivity required for a polyhedron");
2099 newConstituentValue = _polyhedronDescending->getValue();
2100 newConstituentIndex = _polyhedronDescending->getIndex();
2101 lgth=newConstituentIndex[localNumber+1]-newConstituentIndex[localNumber];
2102 return newConstituentValue+newConstituentIndex[localNumber]-1;
2105 throw MEDEXCEPTION("Unknown number");
2109 throw MEDEXCEPTION("No connectivity available");
2113 if(_constituent==NULL)
2114 calculateDescendingConnectivity();
2115 return _constituent->getConnectivityOfAnElementWithPoly(ConnectivityType,Entity,Number,lgth);
2119 int CONNECTIVITY::getNumberOfElementsWithPoly(MED_EN::medEntityMesh Entity, MED_EN::medGeometryElement Type) const
2121 if(Entity==MED_EN::MED_NODE)
2122 return _numberOfNodes;
2125 if(Type==MED_POLYGON || Type==MED_POLYHEDRA)
2126 return getNumberOfElementOfPolyType(_entity);
2128 return getNumberOf(_entity,Type);
2133 return _constituent->getNumberOfElementsWithPoly(Entity,Type);
2135 throw MEDEXCEPTION("CONNECTIVITY::getNumberOfElementsWithPoly : _constituent needed");
2140 Perform a deep comparison of the 2 connectivities in NODAL mode and on all elements.
2142 bool CONNECTIVITY::deepCompare(const CONNECTIVITY& other) const
2144 CONNECTIVITY* temp=(CONNECTIVITY* )this;
2145 const int *conn1=temp->getConnectivity(MED_NODAL,_entity,MED_ALL_ELEMENTS);
2146 int size1=temp->getConnectivityLength(MED_NODAL,_entity,MED_ALL_ELEMENTS);
2147 temp=(CONNECTIVITY* )(&other);
2148 const int *conn2=temp->getConnectivity(MED_NODAL,_entity,MED_ALL_ELEMENTS);
2149 int size2=temp->getConnectivityLength(MED_NODAL,_entity,MED_ALL_ELEMENTS);
2153 for(int i=0;i<size1 && ret;i++)
2155 ret=(conn1[i]==conn2[i]);