Salome HOME
Addition of MEDCoupling1SGTUMesh::structurizeMe
[tools/medcoupling.git] / src / MEDCoupling / MEDCoupling1GTUMesh.cxx
1 // Copyright (C) 2007-2014  CEA/DEN, EDF R&D
2 //
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
7 //
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
11 // Lesser General Public License for more details.
12 //
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
16 //
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
18 //
19 // Author : Anthony Geay (CEA/DEN)
20
21 #include "MEDCoupling1GTUMesh.hxx"
22 #include "MEDCouplingUMesh.hxx"
23 #include "MEDCouplingFieldDouble.hxx"
24 #include "MEDCouplingCMesh.hxx"
25
26 #include "SplitterTetra.hxx"
27
28 using namespace ParaMEDMEM;
29
30 const int MEDCoupling1SGTUMesh::HEXA8_FACE_PAIRS[6]={0,1,2,4,3,5};
31
32 MEDCoupling1GTUMesh::MEDCoupling1GTUMesh():_cm(0)
33 {
34 }
35
36 MEDCoupling1GTUMesh::MEDCoupling1GTUMesh(const std::string& name, const INTERP_KERNEL::CellModel& cm):_cm(&cm)
37 {
38   setName(name);
39 }
40
41 MEDCoupling1GTUMesh::MEDCoupling1GTUMesh(const MEDCoupling1GTUMesh& other, bool recDeepCpy):MEDCouplingPointSet(other,recDeepCpy),_cm(other._cm)
42 {
43 }
44
45 MEDCoupling1GTUMesh *MEDCoupling1GTUMesh::New(const std::string& name, INTERP_KERNEL::NormalizedCellType type)
46 {
47   if(type==INTERP_KERNEL::NORM_ERROR)
48     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::New : NORM_ERROR is not a valid type to be used as base geometric type for a mesh !");
49   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
50   if(!cm.isDynamic())
51     return MEDCoupling1SGTUMesh::New(name,type);
52   else
53     return MEDCoupling1DGTUMesh::New(name,type);
54 }
55
56 MEDCoupling1GTUMesh *MEDCoupling1GTUMesh::New(const MEDCouplingUMesh *m)
57 {
58   if(!m)
59     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::New : input mesh is null !");
60   std::set<INTERP_KERNEL::NormalizedCellType> gts(m->getAllGeoTypes());
61   if(gts.size()!=1)
62     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::New : input mesh must have exactly one geometric type !");
63   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(*gts.begin());
64   if(!cm.isDynamic())
65     return MEDCoupling1SGTUMesh::New(m);
66   else
67     return MEDCoupling1DGTUMesh::New(m);
68 }
69
70 const INTERP_KERNEL::CellModel& MEDCoupling1GTUMesh::getCellModel() const
71 {
72   return *_cm;
73 }
74
75 INTERP_KERNEL::NormalizedCellType MEDCoupling1GTUMesh::getCellModelEnum() const
76 {
77   return _cm->getEnum();
78 }
79
80 int MEDCoupling1GTUMesh::getMeshDimension() const
81 {
82   return (int)_cm->getDimension();
83 }
84
85 /*!
86  * This method returns a newly allocated array containing cell ids (ascendingly sorted) whose geometric type are equal to type.
87  * This method does not throw exception if geometric type \a type is not in \a this.
88  * This method throws an INTERP_KERNEL::Exception if meshdimension of \b this is not equal to those of \b type.
89  * The coordinates array is not considered here.
90  *
91  * \param [in] type the geometric type
92  * \return cell ids in this having geometric type \a type.
93  */
94 DataArrayInt *MEDCoupling1GTUMesh::giveCellsWithType(INTERP_KERNEL::NormalizedCellType type) const
95 {
96   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
97   if(type==getCellModelEnum())
98     ret->alloc(getNumberOfCells(),1);
99   else
100     ret->alloc(0,1);
101   ret->iota();
102   return ret.retn();
103 }
104
105 /*!
106  * Returns nb of cells having the geometric type \a type. No throw if no cells in \a this has the geometric type \a type.
107  */
108 int MEDCoupling1GTUMesh::getNumberOfCellsWithType(INTERP_KERNEL::NormalizedCellType type) const
109 {
110   return type==getCellModelEnum()?getNumberOfCells():0;
111 }
112
113 /*!
114  * Returns a type of a cell by its id.
115  *  \param [in] cellId - the id of the cell of interest.
116  *  \return INTERP_KERNEL::NormalizedCellType - enumeration item describing the cell type.
117  *  \throw If \a cellId is invalid. Valid range is [0, \a this->getNumberOfCells() ).
118  */
119 INTERP_KERNEL::NormalizedCellType MEDCoupling1GTUMesh::getTypeOfCell(int cellId) const
120 {
121   if(cellId>=0 && cellId<getNumberOfCells())
122     return getCellModelEnum();
123   std::ostringstream oss; oss << "MEDCoupling1GTUMesh::getTypeOfCell : Requesting type of cell #" << cellId << " but it should be in [0," << getNumberOfCells() << ") !";
124   throw INTERP_KERNEL::Exception(oss.str().c_str());
125 }
126
127 /*!
128  * Returns a set of all cell types available in \a this mesh.
129  * \return std::set<INTERP_KERNEL::NormalizedCellType> - the set of cell types.
130  * \warning this method does not throw any exception even if \a this is not defined.
131  */
132 std::set<INTERP_KERNEL::NormalizedCellType> MEDCoupling1GTUMesh::getAllGeoTypes() const
133 {
134   std::set<INTERP_KERNEL::NormalizedCellType> ret;
135   ret.insert(getCellModelEnum());
136   return ret;
137 }
138
139 /*!
140  * This method expects that \a this is sorted by types. If not an exception will be thrown.
141  * This method returns in the same format as code (see MEDCouplingUMesh::checkTypeConsistencyAndContig or MEDCouplingUMesh::splitProfilePerType) how
142  * \a this is composed in cell types.
143  * The returned array is of size 3*n where n is the number of different types present in \a this. 
144  * For every k in [0,n] ret[3*k+2]==-1 because it has no sense here. 
145  * This parameter is kept only for compatibility with other methode listed above.
146  */
147 std::vector<int> MEDCoupling1GTUMesh::getDistributionOfTypes() const
148 {
149   std::vector<int> ret(3);
150   ret[0]=(int)getCellModelEnum(); ret[1]=getNumberOfCells(); ret[2]=-1;
151   return ret;
152 }
153
154 /*!
155  * This method is the opposite of MEDCouplingUMesh::checkTypeConsistencyAndContig method. Given a list of cells in \a profile it returns a list of sub-profiles sorted by geo type.
156  * The result is put in the array \a idsPerType. In the returned parameter \a code, foreach i \a code[3*i+2] refers (if different from -1) to a location into the \a idsPerType.
157  * This method has 1 input \a profile and 3 outputs \a code \a idsInPflPerType and \a idsPerType.
158  * 
159  * \param [out] code is a vector of size 3*n where n is the number of different geometric type in \a this \b reduced to the profile \a profile. \a code has exactly the same semantic than in MEDCouplingUMesh::checkTypeConsistencyAndContig method.
160  * \param [out] idsInPflPerType is a vector of size of different geometric type in the subpart defined by \a profile of \a this ( equal to \a code.size()/3). For each i,
161  *              \a idsInPflPerType[i] stores the tuple ids in \a profile that correspond to the geometric type code[3*i+0]
162  * \param [out] idsPerType is a vector of size of different sub profiles needed to be defined to represent the profile \a profile for a given geometric type.
163  *              This vector can be empty in case of all geometric type cells are fully covered in ascending in the given input \a profile.
164  * 
165  * \warning for performance reasons no deep copy will be performed, if \a profile can been used as this in output parameters \a idsInPflPerType and \a idsPerType.
166  *
167  * \throw if \a profile has not exactly one component. It throws too, if \a profile contains some values not in [0,getNumberOfCells()) or if \a this is not fully defined
168  *
169  *  \b Example1: <br>
170  *          - Before \a this has 3 cells \a profile contains [0,1,2]
171  *          - After \a code contains [NORM_...,nbCells,-1], \a idsInPflPerType [[0,1,2]] and \a idsPerType is empty <br>
172  * 
173  *  \b Example2: <br>
174  *          - Before \a this has 3 cells \a profile contains [1,2]
175  *          - After \a code contains [NORM_...,nbCells,0], \a idsInPflPerType [[0,1]] and \a idsPerType is [[1,2]] <br>
176
177  */
178 void MEDCoupling1GTUMesh::splitProfilePerType(const DataArrayInt *profile, std::vector<int>& code, std::vector<DataArrayInt *>& idsInPflPerType, std::vector<DataArrayInt *>& idsPerType) const
179 {
180   if(!profile)
181     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::splitProfilePerType : input profile is NULL !");
182   if(profile->getNumberOfComponents()!=1)
183     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::splitProfilePerType : input profile should have exactly one component !");
184   int nbTuples=profile->getNumberOfTuples();
185   int nbOfCells=getNumberOfCells();
186   code.resize(3); idsInPflPerType.resize(1);
187   code[0]=(int)getCellModelEnum(); code[1]=nbTuples;
188   idsInPflPerType.resize(1);
189   if(profile->isIdentity() && nbTuples==nbOfCells)
190     {
191       code[2]=-1;
192       idsInPflPerType[0]=const_cast<DataArrayInt *>(profile); idsInPflPerType[0]->incrRef();
193       idsPerType.clear();
194       return ;
195     }
196   code[2]=0;
197   profile->checkAllIdsInRange(0,nbOfCells);
198   idsPerType.resize(1);
199   idsPerType[0]=const_cast<DataArrayInt *>(profile); idsPerType[0]->incrRef();
200   idsInPflPerType[0]=DataArrayInt::Range(0,nbTuples,1);
201 }
202
203 /*!
204  * This method tries to minimize at most the number of deep copy.
205  * So if \a idsPerType is not empty it can be returned directly (without copy, but with ref count incremented) in return.
206  * 
207  * \sa MEDCouplingUMesh::checkTypeConsistencyAndContig
208  */
209 DataArrayInt *MEDCoupling1GTUMesh::checkTypeConsistencyAndContig(const std::vector<int>& code, const std::vector<const DataArrayInt *>& idsPerType) const
210 {
211   int nbOfCells=getNumberOfCells();
212   if(code.size()!=3)
213     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : invalid input code should be exactly of size 3 !");
214   if(code[0]!=(int)getCellModelEnum())
215     {
216       std::ostringstream oss; oss << "MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : Mismatch of geometric type ! Asking for " << code[0] << " whereas the geometric type is \a this is " << getCellModelEnum() << " (" << _cm->getRepr() << ") !";
217       throw INTERP_KERNEL::Exception(oss.str().c_str());
218     }
219   if(code[2]==-1)
220     {
221       if(code[1]==nbOfCells)
222         return 0;
223       else
224         {
225           std::ostringstream oss; oss << "MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : mismatch between the number of cells in this (" << nbOfCells << ") and the number of non profile (" << code[1] << ") !";
226           throw INTERP_KERNEL::Exception(oss.str().c_str());
227         }
228     }
229   if(code[2]!=0)
230     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : single geo type mesh ! 0 or -1 is expected at pos #2 of input code !");
231   if(idsPerType.size()!=1)
232     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : input code points to DataArrayInt #0 whereas the size of idsPerType is not equal to 1 !");
233   const DataArrayInt *pfl=idsPerType[0];
234   if(!pfl)
235     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : the input code points to a NULL DataArrayInt at rank 0 !");
236   if(pfl->getNumberOfComponents()!=1)
237     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : input profile should have exactly one component !");
238   pfl->checkAllIdsInRange(0,nbOfCells);
239   pfl->incrRef();
240   return const_cast<DataArrayInt *>(pfl);
241 }
242
243 void MEDCoupling1GTUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const
244 {
245   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
246   m->writeVTKLL(ofs,cellData,pointData,byteData);
247 }
248
249 std::string MEDCoupling1GTUMesh::getVTKDataSetType() const
250 {
251   return std::string("UnstructuredGrid");
252 }
253
254 std::size_t MEDCoupling1GTUMesh::getHeapMemorySizeWithoutChildren() const
255 {
256   return MEDCouplingPointSet::getHeapMemorySizeWithoutChildren();
257 }
258
259 bool MEDCoupling1GTUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
260 {
261   if(!MEDCouplingPointSet::isEqualIfNotWhy(other,prec,reason))
262     return false;
263   if(!other)
264     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::isEqualIfNotWhy : input other pointer is null !");
265   const MEDCoupling1GTUMesh *otherC=dynamic_cast<const MEDCoupling1GTUMesh *>(other);
266   if(!otherC)
267     {
268       reason="mesh given in input is not castable in MEDCouplingSGTUMesh !";
269       return false;
270     }
271   if(_cm!=otherC->_cm)
272     {
273       reason="mismatch in geometric type !";
274       return false;
275     }
276   return true;
277 }
278
279 bool MEDCoupling1GTUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
280 {
281   if(!MEDCouplingPointSet::isEqualWithoutConsideringStr(other,prec))
282     return false;
283   if(!other)
284     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::isEqualWithoutConsideringStr : input other pointer is null !");
285   const MEDCoupling1GTUMesh *otherC=dynamic_cast<const MEDCoupling1GTUMesh *>(other);
286   if(!otherC)
287     return false;
288   if(_cm!=otherC->_cm)
289     return false;
290   return true;
291 }
292
293 void MEDCoupling1GTUMesh::checkCoherency() const
294 {
295   MEDCouplingPointSet::checkCoherency();
296 }
297
298 DataArrayDouble *MEDCoupling1GTUMesh::getBarycenterAndOwner() const
299 {
300   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
301   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=m->getBarycenterAndOwner();
302   return ret.retn();
303 }
304
305 MEDCouplingFieldDouble *MEDCoupling1GTUMesh::getMeasureField(bool isAbs) const
306 {
307   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
308   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=m->getMeasureField(isAbs);
309   ret->setMesh(this);
310   return ret.retn();
311 }
312
313 MEDCouplingFieldDouble *MEDCoupling1GTUMesh::getMeasureFieldOnNode(bool isAbs) const
314 {
315   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
316   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=m->getMeasureFieldOnNode(isAbs);
317   ret->setMesh(this);
318   return ret.retn();
319 }
320
321 /*!
322  * to improve perf !
323  */
324 int MEDCoupling1GTUMesh::getCellContainingPoint(const double *pos, double eps) const
325 {
326   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
327   return m->getCellContainingPoint(pos,eps);
328 }
329
330 MEDCouplingFieldDouble *MEDCoupling1GTUMesh::buildOrthogonalField() const
331 {
332   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
333   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=m->buildOrthogonalField();
334   ret->setMesh(this);
335   return ret.retn();
336 }
337
338 DataArrayInt *MEDCoupling1GTUMesh::getCellsInBoundingBox(const double *bbox, double eps) const
339 {
340   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
341   return m->getCellsInBoundingBox(bbox,eps);
342 }
343
344 DataArrayInt *MEDCoupling1GTUMesh::getCellsInBoundingBox(const INTERP_KERNEL::DirectedBoundingBox& bbox, double eps)
345 {
346   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
347   return m->getCellsInBoundingBox(bbox,eps);
348 }
349
350 MEDCouplingPointSet *MEDCoupling1GTUMesh::buildFacePartOfMySelfNode(const int *start, const int *end, bool fullyIn) const
351 {
352   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
353   return m->buildFacePartOfMySelfNode(start,end,fullyIn);
354 }
355
356 DataArrayInt *MEDCoupling1GTUMesh::findBoundaryNodes() const
357 {
358   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
359   return m->findBoundaryNodes();
360 }
361
362 MEDCouplingPointSet *MEDCoupling1GTUMesh::buildBoundaryMesh(bool keepCoords) const
363 {
364   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
365   return m->buildBoundaryMesh(keepCoords);
366 }
367
368 void MEDCoupling1GTUMesh::findCommonCells(int compType, int startCellId, DataArrayInt *& commonCellsArr, DataArrayInt *& commonCellsIArr) const
369 {
370   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
371   m->findCommonCells(compType,startCellId,commonCellsArr,commonCellsIArr);
372 }
373
374 int MEDCoupling1GTUMesh::getNodalConnectivityLength() const
375 {
376   const DataArrayInt *c1(getNodalConnectivity());
377   if(!c1)
378     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::getNodalConnectivityLength : no connectivity set !");
379   if(c1->getNumberOfComponents()!=1)
380     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::getNodalConnectivityLength : Nodal connectivity array set must have exactly one component !");
381   if(!c1->isAllocated())
382     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::getNodalConnectivityLength : Nodal connectivity array must be allocated !");
383   return c1->getNumberOfTuples();
384 }
385
386 /*!
387  * This method aggregates all the meshes in \a parts to put them in a single unstructured mesh (those returned).
388  * The order of cells is the returned instance is those in the order of instances in \a parts.
389  *
390  * \param [in] parts - all not null parts of single geo type meshes to be aggreagated having the same mesh dimension and same coordinates.
391  * \return MEDCouplingUMesh * - new object to be dealt by the caller.
392  *
393  * \throw If one element is null in \a parts.
394  * \throw If not all the parts do not have the same mesh dimension.
395  * \throw If not all the parts do not share the same coordinates.
396  * \throw If not all the parts have their connectivity set properly.
397  * \throw If \a parts is empty.
398  */
399 MEDCouplingUMesh *MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh(const std::vector< const MEDCoupling1GTUMesh *>& parts)
400 {
401   if(parts.empty())
402     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : input parts vector is empty !");
403   const MEDCoupling1GTUMesh *firstPart(parts[0]);
404   if(!firstPart)
405     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : the first instance in input parts is null !");
406   const DataArrayDouble *coords(firstPart->getCoords());
407   int meshDim(firstPart->getMeshDimension());
408   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret(MEDCouplingUMesh::New(firstPart->getName(),meshDim)); ret->setDescription(firstPart->getDescription());
409   ret->setCoords(coords);
410   int nbOfCells(0),connSize(0);
411   for(std::vector< const MEDCoupling1GTUMesh *>::const_iterator it=parts.begin();it!=parts.end();it++)
412     {
413       if(!(*it))
414         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : presence of null pointer in input vector !");
415       if((*it)->getMeshDimension()!=meshDim)
416         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : all the instances in input vector must have same mesh dimension !");
417       if((*it)->getCoords()!=coords)
418         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : all the instances must share the same coordinates pointer !");
419       nbOfCells+=(*it)->getNumberOfCells();
420       connSize+=(*it)->getNodalConnectivityLength();
421     }
422   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()),connI(DataArrayInt::New());
423   connI->alloc(nbOfCells+1,1); conn->alloc(connSize+nbOfCells,1);
424   int *c(conn->getPointer()),*ci(connI->getPointer()); *ci=0;
425   for(std::vector< const MEDCoupling1GTUMesh *>::const_iterator it=parts.begin();it!=parts.end();it++)
426     {
427       int curNbCells((*it)->getNumberOfCells());
428       int geoType((int)(*it)->getCellModelEnum());
429       const int *cinPtr((*it)->getNodalConnectivity()->begin());
430       const MEDCoupling1SGTUMesh *ps(dynamic_cast<const MEDCoupling1SGTUMesh *>(*it));
431       const MEDCoupling1DGTUMesh *pd(dynamic_cast<const MEDCoupling1DGTUMesh *>(*it));
432       if(ps && !pd)
433         {
434           int nNodesPerCell(ps->getNumberOfNodesPerCell());
435           for(int i=0;i<curNbCells;i++,ci++,cinPtr+=nNodesPerCell)
436             {
437               *c++=geoType;
438               c=std::copy(cinPtr,cinPtr+nNodesPerCell,c);
439               ci[1]=ci[0]+nNodesPerCell+1;
440             }
441         }
442       else if(!ps && pd)
443         {
444           const int *ciinPtr(pd->getNodalConnectivityIndex()->begin());
445           for(int i=0;i<curNbCells;i++,ci++,ciinPtr++)
446             {
447               *c++=geoType;
448               c=std::copy(cinPtr+ciinPtr[0],cinPtr+ciinPtr[1],c);
449               ci[1]=ci[0]+ciinPtr[1]-ciinPtr[0]+1;
450             }
451         }
452       else
453         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : presence of instance which type is not in [MEDCoupling1SGTUMesh,MEDCoupling1DGTUMesh] !");
454     }
455   ret->setConnectivity(conn,connI,true);
456   return ret.retn();
457 }
458
459 //==
460
461 MEDCoupling1SGTUMesh::MEDCoupling1SGTUMesh(const MEDCoupling1SGTUMesh& other, bool recDeepCpy):MEDCoupling1GTUMesh(other,recDeepCpy),_conn(other._conn)
462 {
463   if(recDeepCpy)
464     {
465       const DataArrayInt *c(other._conn);
466       if(c)
467         _conn=c->deepCpy();
468     }
469 }
470
471 MEDCoupling1SGTUMesh::MEDCoupling1SGTUMesh(const std::string& name, const INTERP_KERNEL::CellModel& cm):MEDCoupling1GTUMesh(name,cm)
472 {
473 }
474
475 MEDCoupling1SGTUMesh::MEDCoupling1SGTUMesh()
476 {
477 }
478
479 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::New()
480 {
481   return new MEDCoupling1SGTUMesh;
482 }
483
484 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::New(const std::string& name, INTERP_KERNEL::NormalizedCellType type)
485 {
486   if(type==INTERP_KERNEL::NORM_ERROR)
487     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::New : NORM_ERROR is not a valid type to be used as base geometric type for a mesh !");
488   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
489   if(cm.isDynamic())
490     {
491       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::New : the input geometric type " << cm.getRepr() << " is dynamic ! Only static types are allowed here !";
492       throw INTERP_KERNEL::Exception(oss.str().c_str());
493     }
494   return new MEDCoupling1SGTUMesh(name,cm);
495 }
496
497 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::New(const MEDCouplingUMesh *m)
498 {
499   if(!m)
500     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::New : input mesh is null !");
501   std::set<INTERP_KERNEL::NormalizedCellType> gts(m->getAllGeoTypes());
502   if(gts.size()!=1)
503     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::New : input mesh must have exactly one geometric type !");
504   int geoType((int)*gts.begin());
505   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(MEDCoupling1SGTUMesh::New(m->getName(),*gts.begin()));
506   ret->setCoords(m->getCoords()); ret->setDescription(m->getDescription());
507   int nbCells(m->getNumberOfCells());
508   int nbOfNodesPerCell(ret->getNumberOfNodesPerCell());
509   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()); conn->alloc(nbCells*nbOfNodesPerCell,1);
510   int *c(conn->getPointer());
511   const int *cin(m->getNodalConnectivity()->begin()),*ciin(m->getNodalConnectivityIndex()->begin());
512   for(int i=0;i<nbCells;i++,ciin++)
513     {
514       if(cin[ciin[0]]==geoType)
515         {
516           if(ciin[1]-ciin[0]==nbOfNodesPerCell+1)
517             c=std::copy(cin+ciin[0]+1,cin+ciin[1],c);
518           else
519             {
520               std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The size of cell is not those expected (" << nbOfNodesPerCell << ") !";
521               throw INTERP_KERNEL::Exception(oss.str().c_str());
522             }
523         }
524       else
525         {
526           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The geometric type is not those expected !";
527           throw INTERP_KERNEL::Exception(oss.str().c_str());
528         }
529     }
530   ret->setNodalConnectivity(conn);
531   try
532   { ret->copyTinyInfoFrom(m); }
533   catch(INTERP_KERNEL::Exception&) { }
534   return ret.retn();
535 }
536
537 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::clone(bool recDeepCpy) const
538 {
539   return new MEDCoupling1SGTUMesh(*this,recDeepCpy);
540 }
541
542 /*!
543  * This method behaves mostly like MEDCoupling1SGTUMesh::deepCpy method, except that only nodal connectivity arrays are deeply copied.
544  * The coordinates are shared between \a this and the returned instance.
545  * 
546  * \return MEDCouplingUMesh * - A new object instance holding the copy of \a this (deep for connectivity, shallow for coordiantes)
547  * \sa MEDCoupling1SGTUMesh::deepCpy
548  */
549 MEDCouplingPointSet *MEDCoupling1SGTUMesh::deepCpyConnectivityOnly() const
550 {
551   checkCoherency();
552   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(clone(false));
553   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(_conn->deepCpy());
554   ret->setNodalConnectivity(c);
555   return ret.retn();
556 }
557
558 void MEDCoupling1SGTUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other)
559 {
560   if(!other)
561     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::shallowCopyConnectivityFrom : input pointer is null !");
562   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
563   if(!otherC)
564     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::shallowCopyConnectivityFrom : input pointer is not an MEDCoupling1SGTUMesh instance !");
565   setNodalConnectivity(otherC->getNodalConnectivity());
566 }
567
568 void MEDCoupling1SGTUMesh::updateTime() const
569 {
570   MEDCoupling1GTUMesh::updateTime();
571   const DataArrayInt *c(_conn);
572   if(c)
573     updateTimeWith(*c);
574 }
575
576 std::size_t MEDCoupling1SGTUMesh::getHeapMemorySizeWithoutChildren() const
577 {
578   return MEDCoupling1GTUMesh::getHeapMemorySizeWithoutChildren();
579 }
580
581 std::vector<const BigMemoryObject *> MEDCoupling1SGTUMesh::getDirectChildren() const
582 {
583   std::vector<const BigMemoryObject *> ret(MEDCoupling1GTUMesh::getDirectChildren());
584   const DataArrayInt *c(_conn);
585   if(c)
586     ret.push_back(c);
587   return ret;
588 }
589
590 MEDCouplingMesh *MEDCoupling1SGTUMesh::deepCpy() const
591 {
592   return clone(true);
593 }
594
595 bool MEDCoupling1SGTUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
596 {
597   if(!other)
598     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::isEqualIfNotWhy : input other pointer is null !");
599   std::ostringstream oss; oss.precision(15);
600   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
601   if(!otherC)
602     {
603       reason="mesh given in input is not castable in MEDCoupling1SGTUMesh !";
604       return false;
605     }
606   if(!MEDCoupling1GTUMesh::isEqualIfNotWhy(other,prec,reason))
607     return false;
608   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
609   if(c1==c2)
610     return true;
611   if(!c1 || !c2)
612     {
613       reason="in connectivity of single static geometric type exactly one among this and other is null !";
614       return false;
615     }
616   if(!c1->isEqualIfNotWhy(*c2,reason))
617     {
618       reason.insert(0,"Nodal connectivity DataArrayInt differ : ");
619       return false;
620     }
621   return true;
622 }
623
624 bool MEDCoupling1SGTUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
625 {
626   if(!other)
627     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::isEqualWithoutConsideringStr : input other pointer is null !");
628   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
629   if(!otherC)
630     return false;
631   if(!MEDCoupling1GTUMesh::isEqualWithoutConsideringStr(other,prec))
632     return false;
633   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
634   if(c1==c2)
635     return true;
636   if(!c1 || !c2)
637     return false;
638   if(!c1->isEqualWithoutConsideringStr(*c2))
639     return false;
640   return true;
641 }
642
643 void MEDCoupling1SGTUMesh::checkCoherencyOfConnectivity() const
644 {
645   const DataArrayInt *c1(_conn);
646   if(c1)
647     {
648       if(c1->getNumberOfComponents()!=1)
649         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to be with number of components set to one !");
650       if(c1->getInfoOnComponent(0)!="")
651         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to have no info on its single component !");
652       c1->checkAllocated();
653     }
654   else
655     throw INTERP_KERNEL::Exception("Nodal connectivity array not defined !");
656 }
657
658 void MEDCoupling1SGTUMesh::checkCoherency() const
659 {
660   MEDCouplingPointSet::checkCoherency();
661   checkCoherencyOfConnectivity();
662 }
663
664 void MEDCoupling1SGTUMesh::checkCoherency1(double eps) const
665 {
666   checkCoherency();
667   const DataArrayInt *c1(_conn);
668   int nbOfTuples=c1->getNumberOfTuples();
669   int nbOfNodesPerCell=(int)_cm->getNumberOfNodes();
670   if(nbOfTuples%nbOfNodesPerCell!=0)
671     {
672       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::checkCoherency1 : the nb of tuples in conn is " << nbOfTuples << " and number of nodes per cell is " << nbOfNodesPerCell << ". But " << nbOfTuples << "%" << nbOfNodesPerCell << " !=0 !";
673       throw INTERP_KERNEL::Exception(oss.str().c_str());
674     }
675   int nbOfNodes=getNumberOfNodes();
676   int nbOfCells=nbOfTuples/nbOfNodesPerCell;
677   const int *w(c1->begin());
678   for(int i=0;i<nbOfCells;i++)
679     for(int j=0;j<nbOfNodesPerCell;j++,w++)
680       {
681         if(*w<0 || *w>=nbOfNodes)
682           {
683             std::ostringstream oss; oss << "At node #" << j << " of  cell #" << i << ", is equal to " << *w << " must be in [0," << nbOfNodes << ") !";
684             throw INTERP_KERNEL::Exception(oss.str().c_str());
685           }
686       }
687 }
688
689 void MEDCoupling1SGTUMesh::checkCoherency2(double eps) const
690 {
691   checkCoherency1(eps);
692 }
693
694 int MEDCoupling1SGTUMesh::getNumberOfCells() const
695 {
696   int nbOfTuples=getNodalConnectivityLength();
697   int nbOfNodesPerCell=getNumberOfNodesPerCell();
698   if(nbOfTuples%nbOfNodesPerCell!=0)
699     {
700       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh:getNumberOfCells: : the nb of tuples in conn is " << nbOfTuples << " and number of nodes per cell is " << nbOfNodesPerCell << ". But " << nbOfTuples << "%" << nbOfNodesPerCell << " !=0 !";
701       throw INTERP_KERNEL::Exception(oss.str().c_str());
702     }
703   return nbOfTuples/nbOfNodesPerCell;
704 }
705
706 int MEDCoupling1SGTUMesh::getNumberOfNodesInCell(int cellId) const
707 {
708   return getNumberOfNodesPerCell();
709 }
710
711 int MEDCoupling1SGTUMesh::getNumberOfNodesPerCell() const
712 {
713   checkNonDynamicGeoType();
714   return (int)_cm->getNumberOfNodes();
715 }
716
717 DataArrayInt *MEDCoupling1SGTUMesh::computeNbOfNodesPerCell() const
718 {
719   checkNonDynamicGeoType();
720   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
721   ret->alloc(getNumberOfCells(),1);
722   ret->fillWithValue((int)_cm->getNumberOfNodes());
723   return ret.retn();
724 }
725
726 DataArrayInt *MEDCoupling1SGTUMesh::computeNbOfFacesPerCell() const
727 {
728   checkNonDynamicGeoType();
729   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
730   ret->alloc(getNumberOfCells(),1);
731   ret->fillWithValue((int)_cm->getNumberOfSons());
732   return ret.retn();
733 }
734
735 DataArrayInt *MEDCoupling1SGTUMesh::computeEffectiveNbOfNodesPerCell() const
736 {
737   checkNonDynamicGeoType();
738   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
739   int nbCells(getNumberOfCells());
740   ret->alloc(nbCells,1);
741   int *retPtr(ret->getPointer());
742   int nbNodesPerCell(getNumberOfNodesPerCell());
743   const int *conn(_conn->begin());
744   for(int i=0;i<nbCells;i++,conn+=nbNodesPerCell,retPtr++)
745     {
746       std::set<int> s(conn,conn+nbNodesPerCell);
747       *retPtr=(int)s.size();
748     }
749   return ret.retn();
750 }
751
752 void MEDCoupling1SGTUMesh::getNodeIdsOfCell(int cellId, std::vector<int>& conn) const
753 {
754   int sz=getNumberOfNodesPerCell();
755   conn.resize(sz);
756   if(cellId>=0 && cellId<getNumberOfCells())
757     std::copy(_conn->begin()+cellId*sz,_conn->begin()+(cellId+1)*sz,conn.begin());
758   else
759     {
760       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getNodeIdsOfCell : request for cellId #" << cellId << " must be in [0," << getNumberOfCells() << ") !";
761       throw INTERP_KERNEL::Exception(oss.str().c_str());
762     }
763 }
764
765 void MEDCoupling1SGTUMesh::checkNonDynamicGeoType() const
766 {
767   if(_cm->isDynamic())
768     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkNonDynamicGeoType : internal error ! the internal geo type is dynamic ! should be static !");
769 }
770
771 std::string MEDCoupling1SGTUMesh::simpleRepr() const
772 {
773   static const char msg0[]="No coordinates specified !";
774   std::ostringstream ret;
775   ret << "Single static geometic type (" << _cm->getRepr() << ") unstructured mesh with name : \"" << getName() << "\"\n";
776   ret << "Description of mesh : \"" << getDescription() << "\"\n";
777   int tmpp1,tmpp2;
778   double tt=getTime(tmpp1,tmpp2);
779   ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
780   ret << "Iteration : " << tmpp1  << " Order : " << tmpp2 << "\n";
781   ret << "Mesh dimension : " << getMeshDimension() << "\nSpace dimension : ";
782   if(_coords!=0)
783     {
784       const int spaceDim=getSpaceDimension();
785       ret << spaceDim << "\nInfo attached on space dimension : ";
786       for(int i=0;i<spaceDim;i++)
787         ret << "\"" << _coords->getInfoOnComponent(i) << "\" ";
788       ret << "\n";
789     }
790   else
791     ret << msg0 << "\n";
792   ret << "Number of nodes : ";
793   if(_coords!=0)
794     ret << getNumberOfNodes() << "\n";
795   else
796     ret << msg0 << "\n";
797   ret << "Number of cells : ";
798   if((const DataArrayInt *)_conn)
799     {
800       if(_conn->isAllocated())
801         {
802           if(_conn->getNumberOfComponents()==1)
803             ret << getNumberOfCells() << "\n";
804           else
805             ret << "Nodal connectivity array specified and allocated but with not exactly one component !" << "\n";
806         }
807       else
808         ret << "Nodal connectivity array specified but not allocated !" << "\n";
809     }
810   else
811     ret << "No connectivity specified !" << "\n";
812   ret << "Cell type : " << _cm->getRepr() << "\n";
813   return ret.str();
814 }
815
816 std::string MEDCoupling1SGTUMesh::advancedRepr() const
817 {
818   std::ostringstream ret;
819   ret << simpleRepr();
820   ret << "\nCoordinates array : \n___________________\n\n";
821   if(_coords)
822     _coords->reprWithoutNameStream(ret);
823   else
824     ret << "No array set !\n";
825   ret << "\n\nConnectivity array : \n____________________\n\n";
826   //
827   if((const DataArrayInt *)_conn)
828     {
829       if(_conn->isAllocated())
830         {
831           if(_conn->getNumberOfComponents()==1)
832             {
833               int nbOfCells=getNumberOfCells();
834               int sz=getNumberOfNodesPerCell();
835               const int *connPtr=_conn->begin();
836               for(int i=0;i<nbOfCells;i++,connPtr+=sz)
837                 {
838                   ret << "Cell #" << i << " : ";
839                   std::copy(connPtr,connPtr+sz,std::ostream_iterator<int>(ret," "));
840                   ret << "\n";
841                 }
842             }
843           else
844             ret << "Nodal connectivity array specified and allocated but with not exactly one component !" << "\n";
845         }
846       else
847         ret << "Nodal connectivity array specified but not allocated !" << "\n";
848     }
849   else
850     ret << "No connectivity specified !" << "\n";
851   return ret.str();
852 }
853
854 DataArrayDouble *MEDCoupling1SGTUMesh::computeIsoBarycenterOfNodesPerCell() const
855 {
856   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
857   int spaceDim=getSpaceDimension();
858   int nbOfCells=getNumberOfCells();//checkCoherency()
859   int nbOfNodes=getNumberOfNodes();
860   ret->alloc(nbOfCells,spaceDim);
861   double *ptToFill=ret->getPointer();
862   const double *coor=_coords->begin();
863   const int *nodal=_conn->begin();
864   int sz=getNumberOfNodesPerCell();
865   double coeff=1./(double)sz;
866   for(int i=0;i<nbOfCells;i++,ptToFill+=spaceDim)
867     {
868       std::fill(ptToFill,ptToFill+spaceDim,0.);
869       for(int j=0;j<sz;j++,nodal++)
870         if(*nodal>=0 && *nodal<nbOfNodes)
871           std::transform(coor+spaceDim*nodal[0],coor+spaceDim*(nodal[0]+1),ptToFill,ptToFill,std::plus<double>());
872         else
873           {
874             std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeIsoBarycenterOfNodesPerCell : on cell #" << i << " presence of nodeId #" << *nodal << " should be in [0," <<   nbOfNodes << ") !";
875             throw INTERP_KERNEL::Exception(oss.str().c_str());
876           }
877       std::transform(ptToFill,ptToFill+spaceDim,ptToFill,std::bind2nd(std::multiplies<double>(),coeff));
878     }
879   return ret.retn();
880 }
881
882 void MEDCoupling1SGTUMesh::renumberCells(const int *old2NewBg, bool check)
883 {
884   int nbCells=getNumberOfCells();
885   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> o2n=DataArrayInt::New();
886   o2n->useArray(old2NewBg,false,C_DEALLOC,nbCells,1);
887   if(check)
888     o2n=o2n->checkAndPreparePermutation();
889   //
890   const int *conn=_conn->begin();
891   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> n2o=o2n->invertArrayO2N2N2O(nbCells);
892   const int *n2oPtr=n2o->begin();
893   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New();
894   newConn->alloc(_conn->getNumberOfTuples(),1);
895   newConn->copyStringInfoFrom(*_conn);
896   int sz=getNumberOfNodesPerCell();
897   //
898   int *newC=newConn->getPointer();
899   for(int i=0;i<nbCells;i++,newC+=sz)
900     {
901       int pos=n2oPtr[i];
902       std::copy(conn+pos*sz,conn+(pos+1)*sz,newC);
903     }
904   _conn=newConn;
905 }
906
907 /*!
908  * Keeps from \a this only cells which constituing point id are in the ids specified by [\a begin,\a end).
909  * The resulting cell ids are stored at the end of the 'cellIdsKept' parameter.
910  * Parameter \a fullyIn specifies if a cell that has part of its nodes in ids array is kept or not.
911  * If \a fullyIn is true only cells whose ids are \b fully contained in [\a begin,\a end) tab will be kept.
912  *
913  * \param [in] begin input start of array of node ids.
914  * \param [in] end input end of array of node ids.
915  * \param [in] fullyIn input that specifies if all node ids must be in [\a begin,\a end) array to consider cell to be in.
916  * \param [in,out] cellIdsKeptArr array where all candidate cell ids are put at the end.
917  */
918 void MEDCoupling1SGTUMesh::fillCellIdsToKeepFromNodeIds(const int *begin, const int *end, bool fullyIn, DataArrayInt *&cellIdsKeptArr) const
919 {
920   int nbOfCells=getNumberOfCells();
921   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellIdsKept=DataArrayInt::New(); cellIdsKept->alloc(0,1);
922   int tmp=-1;
923   int sz=_conn->getMaxValue(tmp); sz=std::max(sz,0)+1;
924   std::vector<bool> fastFinder(sz,false);
925   for(const int *work=begin;work!=end;work++)
926     if(*work>=0 && *work<sz)
927       fastFinder[*work]=true;
928   const int *conn=_conn->begin();
929   int nbNodesPerCell=getNumberOfNodesPerCell();
930   for(int i=0;i<nbOfCells;i++,conn+=nbNodesPerCell)
931     {
932       int ref=0,nbOfHit=0;
933       for(int j=0;j<nbNodesPerCell;j++)
934         if(conn[j]>=0)
935           {
936             ref++;
937             if(fastFinder[conn[j]])
938               nbOfHit++;
939           }
940       if((ref==nbOfHit && fullyIn) || (nbOfHit!=0 && !fullyIn))
941         cellIdsKept->pushBackSilent(i);
942     }
943   cellIdsKeptArr=cellIdsKept.retn();
944 }
945
946 MEDCouplingMesh *MEDCoupling1SGTUMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
947 {
948   if(other->getType()!=SINGLE_STATIC_GEO_TYPE_UNSTRUCTURED)
949     throw INTERP_KERNEL::Exception("Merge of umesh only available with umesh single static geo type each other !");
950   const MEDCoupling1SGTUMesh *otherC=static_cast<const MEDCoupling1SGTUMesh *>(other);
951   return Merge1SGTUMeshes(this,otherC);
952 }
953
954 MEDCouplingUMesh *MEDCoupling1SGTUMesh::buildUnstructured() const
955 {
956   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName(),getMeshDimension());
957   ret->setCoords(getCoords());
958   const int *nodalConn=_conn->begin();
959   int nbCells=getNumberOfCells();
960   int nbNodesPerCell=getNumberOfNodesPerCell();
961   int geoType=(int)getCellModelEnum();
962   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=DataArrayInt::New(); c->alloc(nbCells*(nbNodesPerCell+1),1);
963   int *cPtr=c->getPointer();
964   for(int i=0;i<nbCells;i++,nodalConn+=nbNodesPerCell)
965     {
966       *cPtr++=geoType;
967       cPtr=std::copy(nodalConn,nodalConn+nbNodesPerCell,cPtr);
968     }
969   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cI=DataArrayInt::Range(0,(nbCells+1)*(nbNodesPerCell+1),nbNodesPerCell+1);
970   ret->setConnectivity(c,cI,true);
971   try
972   { ret->copyTinyInfoFrom(this); }
973   catch(INTERP_KERNEL::Exception&) { }
974   return ret.retn();
975 }
976
977 DataArrayInt *MEDCoupling1SGTUMesh::simplexize(int policy)
978 {
979   switch(policy)
980   {
981     case 0:
982       return simplexizePol0();
983     case 1:
984       return simplexizePol1();
985     case (int) INTERP_KERNEL::PLANAR_FACE_5:
986         return simplexizePlanarFace5();
987     case (int) INTERP_KERNEL::PLANAR_FACE_6:
988         return simplexizePlanarFace6();
989     default:
990       throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::simplexize : unrecognized policy ! Must be :\n  - 0 or 1 (only available for meshdim=2) \n  - PLANAR_FACE_5, PLANAR_FACE_6  (only for meshdim=3)");
991   }
992 }
993
994 /// @cond INTERNAL
995
996 struct MEDCouplingAccVisit
997 {
998   MEDCouplingAccVisit():_new_nb_of_nodes(0) { }
999   int operator()(int val) { if(val!=-1) return _new_nb_of_nodes++; else return -1; }
1000   int _new_nb_of_nodes;
1001 };
1002
1003 /// @endcond
1004
1005 /*!
1006  * Finds nodes not used in any cell and returns an array giving a new id to every node
1007  * by excluding the unused nodes, for which the array holds -1. The result array is
1008  * a mapping in "Old to New" mode. 
1009  *  \param [out] nbrOfNodesInUse - number of node ids present in the nodal connectivity.
1010  *  \return DataArrayInt * - a new instance of DataArrayInt. Its length is \a
1011  *          this->getNumberOfNodes(). It holds for each node of \a this mesh either -1
1012  *          if the node is unused or a new id else. The caller is to delete this
1013  *          array using decrRef() as it is no more needed.  
1014  *  \throw If the coordinates array is not set.
1015  *  \throw If the nodal connectivity of cells is not defined.
1016  *  \throw If the nodal connectivity includes an invalid id.
1017  */
1018 DataArrayInt *MEDCoupling1SGTUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const
1019 {
1020   nbrOfNodesInUse=-1;
1021   int nbOfNodes=getNumberOfNodes();
1022   int nbOfCells=getNumberOfCells();
1023   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
1024   ret->alloc(nbOfNodes,1);
1025   int *traducer=ret->getPointer();
1026   std::fill(traducer,traducer+nbOfNodes,-1);
1027   const int *conn=_conn->begin();
1028   int nbNodesPerCell=getNumberOfNodesPerCell();
1029   for(int i=0;i<nbOfCells;i++)
1030     for(int j=0;j<nbNodesPerCell;j++,conn++)
1031       if(*conn>=0 && *conn<nbOfNodes)
1032         traducer[*conn]=1;
1033       else
1034         {
1035           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getNodeIdsInUse : In cell #" << i  << " presence of node id " <<  conn[j] << " not in [0," << nbOfNodes << ") !";
1036           throw INTERP_KERNEL::Exception(oss.str().c_str());
1037         }
1038   nbrOfNodesInUse=(int)std::count(traducer,traducer+nbOfNodes,1);
1039   std::transform(traducer,traducer+nbOfNodes,traducer,MEDCouplingAccVisit());
1040   return ret.retn();
1041 }
1042
1043 /*!
1044  * Changes ids of nodes within the nodal connectivity arrays according to a permutation
1045  * array in "Old to New" mode. The node coordinates array is \b not changed by this method.
1046  * This method is a generalization of shiftNodeNumbersInConn().
1047  *  \warning This method performs no check of validity of new ids. **Use it with care !**
1048  *  \param [in] newNodeNumbersO2N - a permutation array, of length \a
1049  *         this->getNumberOfNodes(), in "Old to New" mode. 
1050  *         See \ref MEDCouplingArrayRenumbering for more info on renumbering modes.
1051  *  \throw If the nodal connectivity of cells is not defined.
1052  */
1053 void MEDCoupling1SGTUMesh::renumberNodesInConn(const int *newNodeNumbersO2N)
1054 {
1055   getNumberOfCells();//only to check that all is well defined.
1056   _conn->transformWithIndArr(newNodeNumbersO2N,newNodeNumbersO2N+getNumberOfNodes());
1057   updateTime();
1058 }
1059
1060 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshes(const MEDCoupling1SGTUMesh *mesh1, const MEDCoupling1SGTUMesh *mesh2)
1061 {
1062   std::vector<const MEDCoupling1SGTUMesh *> tmp(2);
1063   tmp[0]=const_cast<MEDCoupling1SGTUMesh *>(mesh1); tmp[1]=const_cast<MEDCoupling1SGTUMesh *>(mesh2);
1064   return Merge1SGTUMeshes(tmp);
1065 }
1066
1067 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshes(std::vector<const MEDCoupling1SGTUMesh *>& a)
1068 {
1069   std::size_t sz=a.size();
1070   if(sz==0)
1071     return Merge1SGTUMeshesLL(a);
1072   for(std::size_t ii=0;ii<sz;ii++)
1073     if(!a[ii])
1074       {
1075         std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::Merge1SGTUMeshes : item #" << ii << " in input array of size "<< sz << " is empty !";
1076         throw INTERP_KERNEL::Exception(oss.str().c_str());
1077       }
1078   const INTERP_KERNEL::CellModel *cm=&(a[0]->getCellModel());
1079   for(std::size_t ii=0;ii<sz;ii++)
1080     if(&(a[ii]->getCellModel())!=cm)
1081       throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshes : all items must have the same geo type !");
1082   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > bb(sz);
1083   std::vector< const MEDCoupling1SGTUMesh * > aa(sz);
1084   int spaceDim=-3;
1085   for(std::size_t i=0;i<sz && spaceDim==-3;i++)
1086     {
1087       const MEDCoupling1SGTUMesh *cur=a[i];
1088       const DataArrayDouble *coo=cur->getCoords();
1089       if(coo)
1090         spaceDim=coo->getNumberOfComponents();
1091     }
1092   if(spaceDim==-3)
1093     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshes : no spaceDim specified ! unable to perform merge !");
1094   for(std::size_t i=0;i<sz;i++)
1095     {
1096       bb[i]=a[i]->buildSetInstanceFromThis(spaceDim);
1097       aa[i]=bb[i];
1098     }
1099   return Merge1SGTUMeshesLL(aa);
1100 }
1101
1102 /*!
1103  * \throw If presence of a null instance in the input vector \a a.
1104  * \throw If a is empty
1105  */
1106 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords(std::vector<const MEDCoupling1SGTUMesh *>& a)
1107 {
1108   if(a.empty())
1109     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : input array must be NON EMPTY !");
1110   std::vector<const MEDCoupling1SGTUMesh *>::const_iterator it=a.begin();
1111   if(!(*it))
1112     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : null instance in the first element of input vector !");
1113   std::vector<const DataArrayInt *> ncs(a.size());
1114   (*it)->getNumberOfCells();//to check that all is OK
1115   const DataArrayDouble *coords=(*it)->getCoords();
1116   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
1117   ncs[0]=(*it)->getNodalConnectivity();
1118   it++;
1119   for(int i=1;it!=a.end();i++,it++)
1120     {
1121       if(!(*it))
1122         throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : presence of a null instance in the input vector !");
1123       if(cm!=&((*it)->getCellModel()))
1124         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1SGTUMeshes impossible !");
1125       (*it)->getNumberOfCells();//to check that all is OK
1126       ncs[i]=(*it)->getNodalConnectivity();
1127       if(coords!=(*it)->getCoords())
1128         throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : not lying on same coords !");
1129     }
1130   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh("merge",*cm));
1131   ret->setCoords(coords);
1132   ret->_conn=DataArrayInt::Aggregate(ncs);
1133   return ret.retn();
1134 }
1135
1136 /*!
1137  * Assume that all instances in \a a are non null. If null it leads to a crash. That's why this method is assigned to be low level (LL)
1138  */
1139 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshesLL(std::vector<const MEDCoupling1SGTUMesh *>& a)
1140 {
1141   if(a.empty())
1142     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshes : input array must be NON EMPTY !");
1143   std::vector<const MEDCoupling1SGTUMesh *>::const_iterator it=a.begin();
1144   int nbOfCells=(*it)->getNumberOfCells();
1145   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
1146   int nbNodesPerCell=(*it)->getNumberOfNodesPerCell();
1147   it++;
1148   for(;it!=a.end();it++)
1149     {
1150       if(cm!=&((*it)->getCellModel()))
1151         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1SGTUMeshes impossible !");
1152       nbOfCells+=(*it)->getNumberOfCells();
1153     }
1154   std::vector<const MEDCouplingPointSet *> aps(a.size());
1155   std::copy(a.begin(),a.end(),aps.begin());
1156   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> pts=MergeNodesArray(aps);
1157   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh("merge",*cm));
1158   ret->setCoords(pts);
1159   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=DataArrayInt::New();
1160   c->alloc(nbOfCells*nbNodesPerCell,1);
1161   int *cPtr=c->getPointer();
1162   int offset=0;
1163   for(it=a.begin();it!=a.end();it++)
1164     {
1165       int curConnLgth=(*it)->getNodalConnectivityLength();
1166       const int *curC=(*it)->_conn->begin();
1167       cPtr=std::transform(curC,curC+curConnLgth,cPtr,std::bind2nd(std::plus<int>(),offset));
1168       offset+=(*it)->getNumberOfNodes();
1169     }
1170   //
1171   ret->setNodalConnectivity(c);
1172   return ret.retn();
1173 }
1174
1175 MEDCouplingPointSet *MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords(const int *begin, const int *end) const
1176 {
1177   int ncell=getNumberOfCells();
1178   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh(getName(),*_cm));
1179   ret->setCoords(_coords);
1180   std::size_t nbOfElemsRet=std::distance(begin,end);
1181   const int *inConn=_conn->getConstPointer();
1182   int sz=getNumberOfNodesPerCell();
1183   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> connRet=DataArrayInt::New(); connRet->alloc((int)nbOfElemsRet*sz,1);
1184   int *connPtr=connRet->getPointer();
1185   for(const int *work=begin;work!=end;work++,connPtr+=sz)
1186     {
1187       if(*work>=0 && *work<ncell)
1188         std::copy(inConn+(work[0])*sz,inConn+(work[0]+1)*sz,connPtr);
1189       else
1190         {
1191           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords : On pos #" << std::distance(begin,work) << " input cell id =" << *work << " should be in [0," << ncell << ") !";
1192           throw INTERP_KERNEL::Exception(oss.str().c_str());
1193         }
1194     }
1195   ret->_conn=connRet;
1196   ret->copyTinyInfoFrom(this);
1197   return ret.retn();
1198 }
1199
1200 MEDCouplingPointSet *MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords2(int start, int end, int step) const
1201 {
1202   int ncell=getNumberOfCells();
1203   int nbOfElemsRet=DataArray::GetNumberOfItemGivenBESRelative(start,end,step,"MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords2 : ");
1204   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh(getName(),*_cm));
1205   ret->setCoords(_coords);
1206   const int *inConn=_conn->getConstPointer();
1207   int sz=getNumberOfNodesPerCell();
1208   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> connRet=DataArrayInt::New(); connRet->alloc((int)nbOfElemsRet*sz,1);
1209   int *connPtr=connRet->getPointer();
1210   int curId=start;
1211   for(int i=0;i<nbOfElemsRet;i++,connPtr+=sz,curId+=step)
1212     {
1213       if(curId>=0 && curId<ncell)
1214         std::copy(inConn+curId*sz,inConn+(curId+1)*sz,connPtr);
1215       else
1216         {
1217           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords2 : On pos #" << i << " input cell id =" << curId  << " should be in [0," << ncell << ") !";
1218           throw INTERP_KERNEL::Exception(oss.str().c_str());
1219         }
1220     }
1221   ret->_conn=connRet;
1222   ret->copyTinyInfoFrom(this);
1223   return ret.retn();
1224 }
1225
1226 void MEDCoupling1SGTUMesh::computeNodeIdsAlg(std::vector<bool>& nodeIdsInUse) const
1227 {
1228   int sz((int)nodeIdsInUse.size());
1229   int nbCells(getNumberOfCells());
1230   int nbOfNodesPerCell(getNumberOfNodesPerCell());
1231   const int *w(_conn->begin());
1232   for(int i=0;i<nbCells;i++)
1233     for(int j=0;j<nbOfNodesPerCell;j++,w++)
1234       {
1235         if(*w>=0 && *w<sz)
1236           nodeIdsInUse[*w]=true;
1237         else
1238           {
1239             std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeNodeIdsAlg : At cell #" << i << " presence of node id #" << *w << " should be in [0," << sz << ") !";
1240             throw INTERP_KERNEL::Exception(oss.str().c_str());
1241           }
1242       }
1243 }
1244
1245 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::buildSetInstanceFromThis(int spaceDim) const
1246 {
1247   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh(getName(),*_cm));
1248   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp1;
1249   const DataArrayInt *nodalConn(_conn);
1250   if(!nodalConn)
1251     {
1252       tmp1=DataArrayInt::New(); tmp1->alloc(0,1);
1253     }
1254   else
1255     tmp1=_conn;
1256   ret->_conn=tmp1;
1257   if(!_coords)
1258     {
1259       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords=DataArrayDouble::New(); coords->alloc(0,spaceDim);
1260       ret->setCoords(coords);
1261     }
1262   else
1263     ret->setCoords(_coords);
1264   return ret.retn();
1265 }
1266
1267 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePol0()
1268 {
1269   int nbOfCells=getNumberOfCells();
1270   if(getCellModelEnum()!=INTERP_KERNEL::NORM_QUAD4)
1271     return DataArrayInt::Range(0,nbOfCells,1);
1272   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(2*3*nbOfCells,1);
1273   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(2*nbOfCells,1);
1274   const int *c(_conn->begin());
1275   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1276   for(int i=0;i<nbOfCells;i++,c+=4,newConnPtr+=6,retPtr+=2)
1277     {
1278       newConnPtr[0]=c[0]; newConnPtr[1]=c[1]; newConnPtr[2]=c[2];
1279       newConnPtr[3]=c[0]; newConnPtr[4]=c[2]; newConnPtr[5]=c[3];
1280       retPtr[0]=i; retPtr[1]=i;
1281     }
1282   _conn=newConn;
1283   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TRI3);
1284   updateTime();
1285   return ret.retn();
1286 }
1287
1288 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePol1()
1289 {
1290   int nbOfCells=getNumberOfCells();
1291   if(getCellModelEnum()!=INTERP_KERNEL::NORM_QUAD4)
1292     return DataArrayInt::Range(0,nbOfCells,1);
1293   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(2*3*nbOfCells,1);
1294   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(2*nbOfCells,1);
1295   const int *c(_conn->begin());
1296   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1297   for(int i=0;i<nbOfCells;i++,c+=4,newConnPtr+=6,retPtr+=2)
1298     {
1299       newConnPtr[0]=c[0]; newConnPtr[1]=c[1]; newConnPtr[2]=c[3];
1300       newConnPtr[3]=c[1]; newConnPtr[4]=c[2]; newConnPtr[5]=c[3];
1301       retPtr[0]=i; retPtr[1]=i;
1302     }
1303   _conn=newConn;
1304   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TRI3);
1305   updateTime();
1306   return ret.retn();
1307 }
1308
1309 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePlanarFace5()
1310 {
1311   int nbOfCells=getNumberOfCells();
1312   if(getCellModelEnum()!=INTERP_KERNEL::NORM_HEXA8)
1313     return DataArrayInt::Range(0,nbOfCells,1);
1314   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(5*4*nbOfCells,1);
1315   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(5*nbOfCells,1);
1316   const int *c(_conn->begin());
1317   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1318   for(int i=0;i<nbOfCells;i++,c+=8,newConnPtr+=20,retPtr+=5)
1319     {
1320       for(int j=0;j<20;j++)
1321         newConnPtr[j]=c[INTERP_KERNEL::SPLIT_NODES_5_WO[j]];
1322       retPtr[0]=i; retPtr[1]=i; retPtr[2]=i; retPtr[3]=i; retPtr[4]=i;
1323     }
1324   _conn=newConn;
1325   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TETRA4);
1326   updateTime();
1327   return ret.retn();
1328 }
1329
1330 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePlanarFace6()
1331 {
1332   int nbOfCells=getNumberOfCells();
1333   if(getCellModelEnum()!=INTERP_KERNEL::NORM_HEXA8)
1334     return DataArrayInt::Range(0,nbOfCells,1);
1335   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(6*4*nbOfCells,1);
1336   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(6*nbOfCells,1);
1337   const int *c(_conn->begin());
1338   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1339   for(int i=0;i<nbOfCells;i++,c+=8,newConnPtr+=24,retPtr+=6)
1340     {
1341       for(int j=0;j<24;j++)
1342         newConnPtr[j]=c[INTERP_KERNEL::SPLIT_NODES_6_WO[j]];
1343       retPtr[0]=i; retPtr[1]=i; retPtr[2]=i; retPtr[3]=i; retPtr[4]=i; retPtr[5]=i;
1344     }
1345   _conn=newConn;
1346   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TETRA4);
1347   updateTime();
1348   return ret.retn();
1349 }
1350
1351 void MEDCoupling1SGTUMesh::reprQuickOverview(std::ostream& stream) const
1352 {
1353   stream << "MEDCoupling1SGTUMesh C++ instance at " << this << ". Type=" << _cm->getRepr() << ". Name : \"" << getName() << "\".";
1354   stream << " Mesh dimension : " << getMeshDimension() << ".";
1355   if(!_coords)
1356     { stream << " No coordinates set !"; return ; }
1357   if(!_coords->isAllocated())
1358     { stream << " Coordinates set but not allocated !"; return ; }
1359   stream << " Space dimension : " << _coords->getNumberOfComponents() << "." << std::endl;
1360   stream << "Number of nodes : " << _coords->getNumberOfTuples() << ".";
1361   if(!(const DataArrayInt *)_conn)
1362     { stream << std::endl << "Nodal connectivity NOT set !"; return ; }
1363   if(_conn->isAllocated())
1364     {
1365       if(_conn->getNumberOfComponents()==1)
1366         stream << std::endl << "Number of cells : " << getNumberOfCells() << ".";
1367     }
1368 }
1369
1370 void MEDCoupling1SGTUMesh::checkFullyDefined() const
1371 {
1372   if(!((const DataArrayInt *)_conn) || !((const DataArrayDouble *)_coords))
1373     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFullyDefined : part of this is not fully defined.");
1374 }
1375
1376 /*!
1377  * First step of unserialization process.
1378  */
1379 bool MEDCoupling1SGTUMesh::isEmptyMesh(const std::vector<int>& tinyInfo) const
1380 {
1381   throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::isEmptyMesh : not implemented yet !");
1382 }
1383
1384 void MEDCoupling1SGTUMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
1385 {
1386   int it,order;
1387   double time=getTime(it,order);
1388   tinyInfo.clear(); tinyInfoD.clear(); littleStrings.clear();
1389   //
1390   littleStrings.push_back(getName());
1391   littleStrings.push_back(getDescription());
1392   littleStrings.push_back(getTimeUnit());
1393   //
1394   std::vector<std::string> littleStrings2,littleStrings3;
1395   if((const DataArrayDouble *)_coords)
1396     _coords->getTinySerializationStrInformation(littleStrings2);
1397   if((const DataArrayInt *)_conn)
1398     _conn->getTinySerializationStrInformation(littleStrings3);
1399   int sz0((int)littleStrings2.size()),sz1((int)littleStrings3.size());
1400   littleStrings.insert(littleStrings.end(),littleStrings2.begin(),littleStrings2.end());
1401   littleStrings.insert(littleStrings.end(),littleStrings3.begin(),littleStrings3.end());
1402   //
1403   tinyInfo.push_back(getCellModelEnum());
1404   tinyInfo.push_back(it);
1405   tinyInfo.push_back(order);
1406   std::vector<int> tinyInfo2,tinyInfo3;
1407   if((const DataArrayDouble *)_coords)
1408     _coords->getTinySerializationIntInformation(tinyInfo2);
1409   if((const DataArrayInt *)_conn)
1410     _conn->getTinySerializationIntInformation(tinyInfo3);
1411   int sz2((int)tinyInfo2.size()),sz3((int)tinyInfo3.size());
1412   tinyInfo.push_back(sz0); tinyInfo.push_back(sz1); tinyInfo.push_back(sz2); tinyInfo.push_back(sz3);
1413   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
1414   tinyInfo.insert(tinyInfo.end(),tinyInfo3.begin(),tinyInfo3.end());
1415   //
1416   tinyInfoD.push_back(time);
1417 }
1418
1419 void MEDCoupling1SGTUMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
1420 {
1421   std::vector<int> tinyInfo2(tinyInfo.begin()+7,tinyInfo.begin()+7+tinyInfo[5]);
1422   std::vector<int> tinyInfo1(tinyInfo.begin()+7+tinyInfo[5],tinyInfo.begin()+7+tinyInfo[5]+tinyInfo[6]);
1423   a1->resizeForUnserialization(tinyInfo1);
1424   a2->resizeForUnserialization(tinyInfo2);
1425 }
1426
1427 void MEDCoupling1SGTUMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
1428 {
1429   int sz(0);
1430   if((const DataArrayInt *)_conn)
1431     if(_conn->isAllocated())
1432       sz=_conn->getNbOfElems();
1433   a1=DataArrayInt::New();
1434   a1->alloc(sz,1);
1435   if(sz!=0 && (const DataArrayInt *)_conn)
1436     std::copy(_conn->begin(),_conn->end(),a1->getPointer());
1437   sz=0;
1438   if((const DataArrayDouble *)_coords)
1439     if(_coords->isAllocated())
1440       sz=_coords->getNbOfElems();
1441   a2=DataArrayDouble::New();
1442   a2->alloc(sz,1);
1443   if(sz!=0 && (const DataArrayDouble *)_coords)
1444     std::copy(_coords->begin(),_coords->end(),a2->getPointer());
1445 }
1446
1447 void MEDCoupling1SGTUMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
1448                                            const std::vector<std::string>& littleStrings)
1449 {
1450   INTERP_KERNEL::NormalizedCellType gt((INTERP_KERNEL::NormalizedCellType)tinyInfo[0]);
1451   _cm=&INTERP_KERNEL::CellModel::GetCellModel(gt);
1452   setName(littleStrings[0]);
1453   setDescription(littleStrings[1]);
1454   setTimeUnit(littleStrings[2]);
1455   setTime(tinyInfoD[0],tinyInfo[1],tinyInfo[2]);
1456   int sz0(tinyInfo[3]),sz1(tinyInfo[4]),sz2(tinyInfo[5]),sz3(tinyInfo[6]);
1457   //
1458   _coords=DataArrayDouble::New();
1459   std::vector<int> tinyInfo2(tinyInfo.begin()+7,tinyInfo.begin()+7+sz2);
1460   _coords->resizeForUnserialization(tinyInfo2);
1461   std::copy(a2->begin(),a2->end(),_coords->getPointer());
1462   _conn=DataArrayInt::New();
1463   std::vector<int> tinyInfo3(tinyInfo.begin()+7+sz2,tinyInfo.begin()+7+sz2+sz3);
1464   _conn->resizeForUnserialization(tinyInfo3);
1465   std::copy(a1->begin(),a1->end(),_conn->getPointer());
1466   std::vector<std::string> littleStrings2(littleStrings.begin()+3,littleStrings.begin()+3+sz0);
1467   _coords->finishUnserialization(tinyInfo2,littleStrings2);
1468   std::vector<std::string> littleStrings3(littleStrings.begin()+3+sz0,littleStrings.begin()+3+sz0+sz1);
1469   _conn->finishUnserialization(tinyInfo3,littleStrings3);
1470 }
1471
1472 /*!
1473  * Checks if \a this and \a other meshes are geometrically equivalent with high
1474  * probability, else an exception is thrown. The meshes are considered equivalent if
1475  * (1) meshes contain the same number of nodes and the same number of elements of the
1476  * same types (2) three cells of the two meshes (first, last and middle) are based
1477  * on coincident nodes (with a specified precision).
1478  *  \param [in] other - the mesh to compare with.
1479  *  \param [in] prec - the precision used to compare nodes of the two meshes.
1480  *  \throw If the two meshes do not match.
1481  */
1482 void MEDCoupling1SGTUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double prec) const
1483 {
1484   MEDCouplingPointSet::checkFastEquivalWith(other,prec);
1485   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
1486   if(!otherC)
1487     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : Two meshes are not unstructured with single static geometric type !");
1488   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
1489   if(c1==c2)
1490     return;
1491   if(!c1 || !c2)
1492     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : presence of nodal connectivity only in one of the 2 meshes !");
1493   if((c1->isAllocated() && !c2->isAllocated()) || (!c1->isAllocated() && c2->isAllocated()))
1494     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : in nodal connectivity, only one is allocated !");
1495   if(c1->getNumberOfComponents()!=1 || c1->getNumberOfComponents()!=1)
1496     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : in nodal connectivity, must have 1 and only 1 component !");
1497   if(c1->getHashCode()!=c2->getHashCode())
1498     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : nodal connectivity differs");
1499 }
1500
1501 MEDCouplingPointSet *MEDCoupling1SGTUMesh::mergeMyselfWithOnSameCoords(const MEDCouplingPointSet *other) const
1502 {
1503   if(!other)
1504     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::mergeMyselfWithOnSameCoords : input other is null !");
1505   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
1506   if(!otherC)
1507     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::mergeMyselfWithOnSameCoords : the input other mesh is not of type single statuc geo type unstructured !");
1508   std::vector<const MEDCoupling1SGTUMesh *> ms(2);
1509   ms[0]=this;
1510   ms[1]=otherC;
1511   return Merge1SGTUMeshesOnSameCoords(ms);
1512 }
1513
1514 void MEDCoupling1SGTUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const
1515 {
1516   checkFullyDefined();
1517   int nbOfNodes=getNumberOfNodes();
1518   int *revNodalIndxPtr=(int *)malloc((nbOfNodes+1)*sizeof(int));
1519   revNodalIndx->useArray(revNodalIndxPtr,true,C_DEALLOC,nbOfNodes+1,1);
1520   std::fill(revNodalIndxPtr,revNodalIndxPtr+nbOfNodes+1,0);
1521   const int *conn=_conn->begin();
1522   int nbOfCells=getNumberOfCells();
1523   int nbOfEltsInRevNodal=0;
1524   int nbOfNodesPerCell=getNumberOfNodesPerCell();
1525   for(int eltId=0;eltId<nbOfCells;eltId++)
1526     {
1527       for(int j=0;j<nbOfNodesPerCell;j++,conn++)
1528         {
1529           if(conn[0]>=0 && conn[0]<nbOfNodes)
1530             {
1531               nbOfEltsInRevNodal++;
1532               revNodalIndxPtr[conn[0]+1]++;
1533             }
1534           else
1535             {
1536               std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getReverseNodalConnectivity : At cell #" << eltId << " presence of nodeId #" << conn[0] << " should be in [0," << nbOfNodes << ") !";
1537               throw INTERP_KERNEL::Exception(oss.str().c_str());
1538             }
1539         }
1540     }
1541   std::transform(revNodalIndxPtr+1,revNodalIndxPtr+nbOfNodes+1,revNodalIndxPtr,revNodalIndxPtr+1,std::plus<int>());
1542   conn=_conn->begin();
1543   int *revNodalPtr=(int *)malloc((nbOfEltsInRevNodal)*sizeof(int));
1544   revNodal->useArray(revNodalPtr,true,C_DEALLOC,nbOfEltsInRevNodal,1);
1545   std::fill(revNodalPtr,revNodalPtr+nbOfEltsInRevNodal,-1);
1546   for(int eltId=0;eltId<nbOfCells;eltId++)
1547     {
1548       for(int j=0;j<nbOfNodesPerCell;j++,conn++)
1549         {
1550           *std::find_if(revNodalPtr+revNodalIndxPtr[*conn],revNodalPtr+revNodalIndxPtr[*conn+1],std::bind2nd(std::equal_to<int>(),-1))=eltId;
1551         }
1552     }
1553 }
1554
1555 /*!
1556  * Use \a nodalConn array as nodal connectivity of \a this. The input \a nodalConn pointer can be null.
1557  */
1558 void MEDCoupling1SGTUMesh::setNodalConnectivity(DataArrayInt *nodalConn)
1559 {
1560   if(nodalConn)
1561     nodalConn->incrRef();
1562   _conn=nodalConn;
1563   declareAsNew();
1564 }
1565
1566 /*!
1567  * \return DataArrayInt * - the internal reference to the nodal connectivity. The caller is not reponsible to deallocate it.
1568  */
1569 DataArrayInt *MEDCoupling1SGTUMesh::getNodalConnectivity() const
1570 {
1571   const DataArrayInt *ret(_conn);
1572   return const_cast<DataArrayInt *>(ret);
1573 }
1574
1575 /*!
1576  * Allocates memory to store an estimation of the given number of cells. Closer is the estimation to the number of cells effectively inserted,
1577  * less will be the needs to realloc. If the number of cells to be inserted is not known simply put 0 to this parameter.
1578  * If a nodal connectivity previouly existed before the call of this method, it will be reset.
1579  *
1580  *  \param [in] nbOfCells - estimation of the number of cell \a this mesh will contain.
1581  */
1582 void MEDCoupling1SGTUMesh::allocateCells(int nbOfCells)
1583 {
1584   if(nbOfCells<0)
1585     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::allocateCells : the input number of cells should be >= 0 !");
1586   _conn=DataArrayInt::New();
1587   _conn->reserve(getNumberOfNodesPerCell()*nbOfCells);
1588   declareAsNew();
1589 }
1590
1591 /*!
1592  * Appends at the end of \a this a cell having nodal connectivity array defined in [ \a nodalConnOfCellBg, \a nodalConnOfCellEnd ).
1593  *
1594  * \param [in] nodalConnOfCellBg - the begin (included) of nodal connectivity of the cell to add.
1595  * \param [in] nodalConnOfCellEnd - the end (excluded) of nodal connectivity of the cell to add.
1596  * \throw If the length of the input nodal connectivity array of the cell to add is not equal to number of nodes per cell relative to the unique geometric type
1597  *        attached to \a this.
1598  * \thow If the nodal connectivity array in \a this is null (call MEDCoupling1SGTUMesh::allocateCells before).
1599  */
1600 void MEDCoupling1SGTUMesh::insertNextCell(const int *nodalConnOfCellBg, const int *nodalConnOfCellEnd)
1601 {
1602   int sz=(int)std::distance(nodalConnOfCellBg,nodalConnOfCellEnd);
1603   int ref=getNumberOfNodesPerCell();
1604   if(sz==ref)
1605     {
1606       DataArrayInt *c(_conn);
1607       if(c)
1608         c->pushBackValsSilent(nodalConnOfCellBg,nodalConnOfCellEnd);
1609       else
1610         throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::insertNextCell : nodal connectivity array is null ! Call MEDCoupling1SGTUMesh::allocateCells before !");
1611     }
1612   else
1613     {
1614       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::insertNextCell : input nodal size (" << sz << ") does not match number of nodes per cell of this (";
1615       oss << ref << ") !";
1616       throw INTERP_KERNEL::Exception(oss.str().c_str());
1617     }
1618 }
1619
1620 /*!
1621  * This method builds the dual mesh of \a this and returns it.
1622  * 
1623  * \return MEDCoupling1SGTUMesh * - newly object created to be managed by the caller.
1624  * \throw If \a this is not a mesh containing only simplex cells.
1625  * \throw If \a this is not correctly allocated (coordinates and connectivities have to be correctly set !).
1626  * \throw If at least one node in \a this is orphan (without any simplex cell lying on it !)
1627  */
1628 MEDCoupling1GTUMesh *MEDCoupling1SGTUMesh::computeDualMesh() const
1629 {
1630   const INTERP_KERNEL::CellModel& cm(getCellModel());
1631   if(!cm.isSimplex())
1632     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh : this mesh is not a simplex mesh ! Please invoke simplexize of tetrahedrize on this before calling this method !");
1633   switch(getMeshDimension())
1634   {
1635     case 3:
1636       return computeDualMesh3D();
1637     case 2:
1638       return computeDualMesh2D();
1639     default:
1640       throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh : meshdimension must be in [2,3] !");
1641   }
1642 }
1643
1644 /*!
1645  * This method explode each NORM_HEXA8 cells in \a this into 6 NORM_QUAD4 cells and put the result into the MEDCoupling1SGTUMesh returned instance.
1646  * 
1647  * \return MEDCoupling1SGTUMesh * - a newly allocated instances (to be managed by the caller) storing the result of the explosion.
1648  * \throw If \a this is not a mesh containing only NORM_HEXA8 cells.
1649  * \throw If \a this is not properly allocated.
1650  */
1651 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::explodeEachHexa8To6Quad4() const
1652 {
1653   const INTERP_KERNEL::CellModel& cm(getCellModel());
1654   if(cm.getEnum()!=INTERP_KERNEL::NORM_HEXA8)
1655     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::explodeEachHexa8To6Quad4 : this method can be applied only on HEXA8 mesh !");
1656   int nbHexa8(getNumberOfCells());
1657   const int *inConnPtr(getNodalConnectivity()->begin());
1658   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(MEDCoupling1SGTUMesh::New(getName(),INTERP_KERNEL::NORM_QUAD4));
1659   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(DataArrayInt::New()); c->alloc(nbHexa8*6*4,1);
1660   int *cPtr(c->getPointer());
1661   for(int i=0;i<nbHexa8;i++,inConnPtr+=8)
1662     {
1663       for(int j=0;j<6;j++,cPtr+=4)
1664         cm.fillSonCellNodalConnectivity(j,inConnPtr,cPtr);
1665     }
1666   ret->setCoords(getCoords());
1667   ret->setNodalConnectivity(c);
1668   return ret.retn();
1669 }
1670
1671 /*!
1672  * This method starts from an unstructured mesh that hides in reality a cartesian mesh.
1673  * If it is not the case, an exception will be thrown.
1674  * This method returns three objects : The cartesian mesh geometrically equivalent to \a this (within a precision of \a eps) and a permutation of cells
1675  * and a permutation of nodes.
1676  *
1677  * \param [out] cellPerm the permutation array of size \c this->getNumberOfCells()
1678  * \param [out] nodePerm the permutation array of size \c this->getNumberOfNodes()
1679  * \return MEDCouplingCMesh * - a newly allocated mesh that is the result of the structurization of \a this.
1680  */
1681 MEDCouplingCMesh *MEDCoupling1SGTUMesh::structurizeMe(DataArrayInt *& cellPerm, DataArrayInt *& nodePerm, double eps) const
1682 {
1683   checkCoherency();
1684   int spaceDim(getSpaceDimension()),meshDim(getMeshDimension()),nbNodes(getNumberOfNodes());
1685   if(MEDCouplingStructuredMesh::GetGeoTypeGivenMeshDimension(meshDim)!=getCellModelEnum())
1686     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::structurizeMe : the unique geo type in this is not compatible with the geometric type regarding mesh dimension !");
1687   MEDCouplingAutoRefCountObjectPtr<MEDCouplingCMesh> cm(MEDCouplingCMesh::New());
1688   for(int i=0;i<spaceDim;i++)
1689     {
1690       std::vector<int> tmp(1,i);
1691       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> elt(static_cast<DataArrayDouble*>(getCoords()->keepSelectedComponents(tmp)));
1692       elt=elt->getDifferentValues(eps);
1693       elt->sort(true);
1694       cm->setCoordsAt(i,elt);
1695     }
1696   if(nbNodes!=cm->getNumberOfNodes())
1697     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::structurizeMe : considering the number of nodes after split per components in space this can't be a cartesian mesh ! Maybe your epsilon parameter is invalid ?");
1698   try
1699   { cm->copyTinyInfoFrom(this); }
1700   catch(INTERP_KERNEL::Exception&) { }
1701   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> um(cm->buildUnstructured()),self(buildUnstructured());
1702   self->checkGeoEquivalWith(um,12,eps,cellPerm,nodePerm);
1703   return cm.retn();
1704 }
1705
1706 /// @cond INTERNAL
1707
1708 bool UpdateHexa8Cell(int validAxis, int neighId, const int *validConnQuad4NeighSide, int *allFacesNodalConn, int *myNeighbours)
1709 {
1710   static const int TAB[48]={
1711     0,1,2,3,4,5,6,7,//0
1712     4,7,6,5,0,3,2,1,//1
1713     0,3,7,4,1,2,6,5,//2
1714     4,0,3,7,5,1,2,6,//3
1715     5,1,0,4,6,2,3,7,//4
1716     3,7,4,0,2,6,5,1 //5
1717   };
1718   static const int TAB2[6]={0,0,3,3,3,3};
1719   if(myNeighbours[validAxis]==neighId && allFacesNodalConn[4*validAxis+0]==validConnQuad4NeighSide[TAB2[validAxis]])
1720     return true;
1721   int oldAxis((int)std::distance(myNeighbours,std::find(myNeighbours,myNeighbours+6,neighId)));
1722   std::size_t pos(std::distance(MEDCoupling1SGTUMesh::HEXA8_FACE_PAIRS,std::find(MEDCoupling1SGTUMesh::HEXA8_FACE_PAIRS,MEDCoupling1SGTUMesh::HEXA8_FACE_PAIRS+6,oldAxis)));
1723   std::size_t pos0(pos/2),pos1(pos%2);
1724   int oldAxisOpp(MEDCoupling1SGTUMesh::HEXA8_FACE_PAIRS[2*pos0+(pos1+1)%2]);
1725   int oldConn[8],myConn2[8]={-1,-1,-1,-1,-1,-1,-1,-1},myConn[8],edgeConn[2],allFacesTmp[24],neighTmp[6];
1726   oldConn[0]=allFacesNodalConn[0]; oldConn[1]=allFacesNodalConn[1]; oldConn[2]=allFacesNodalConn[2]; oldConn[3]=allFacesNodalConn[3];
1727   oldConn[4]=allFacesNodalConn[4]; oldConn[5]=allFacesNodalConn[7]; oldConn[6]=allFacesNodalConn[6]; oldConn[7]=allFacesNodalConn[5];
1728   const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_HEXA8));
1729   for(int i=0;i<4;i++)
1730     myConn2[i]=validConnQuad4NeighSide[(4-i+TAB2[validAxis])%4];
1731   for(int i=0;i<4;i++)
1732     {
1733       int nodeId(myConn2[i]);//the node id for which the opposite one will be found
1734       bool found(false);
1735       INTERP_KERNEL::NormalizedCellType typeOfSon;
1736       for(int j=0;j<12 && !found;j++)
1737         {
1738           cm.fillSonEdgesNodalConnectivity3D(j,oldConn,-1,edgeConn,typeOfSon);
1739           if(edgeConn[0]==nodeId || edgeConn[1]==nodeId)
1740             {
1741               if(std::find(allFacesNodalConn+4*oldAxisOpp,allFacesNodalConn+4*oldAxisOpp+4,edgeConn[0]==nodeId?edgeConn[1]:edgeConn[0])!=allFacesNodalConn+4*oldAxisOpp+4)
1742                 {
1743                   myConn2[i+4]=edgeConn[0]==nodeId?edgeConn[1]:edgeConn[0];
1744                   found=true;
1745                 }
1746             }
1747         }
1748       if(!found)
1749         throw INTERP_KERNEL::Exception("UpdateHexa8Cell : Internal Error !");
1750     }
1751   const int *myTab(TAB+8*validAxis);
1752   for(int i=0;i<8;i++)
1753     myConn[i]=myConn2[myTab[i]];
1754   for(int i=0;i<6;i++)
1755     {
1756       cm.fillSonCellNodalConnectivity(i,myConn,allFacesTmp+4*i);
1757       std::set<int> s(allFacesTmp+4*i,allFacesTmp+4*i+4);
1758       bool found(false);
1759       for(int j=0;j<6 && !found;j++)
1760         {
1761           std::set<int> s1(allFacesNodalConn+4*j,allFacesNodalConn+4*j+4);
1762           if(s==s1)
1763             {
1764               neighTmp[i]=myNeighbours[j];
1765               found=true;
1766             }
1767         }
1768       if(!found)
1769         throw INTERP_KERNEL::Exception("UpdateHexa8Cell : Internal Error #2 !");
1770     }
1771   std::copy(allFacesTmp,allFacesTmp+24,allFacesNodalConn);
1772   std::copy(neighTmp,neighTmp+6,myNeighbours);
1773   return false;
1774 }
1775
1776 /// @endcond
1777
1778 /*!
1779  * This method expects the \a this contains NORM_HEXA8 cells only. This method will sort each cells in \a this so that their numbering was
1780  * homogeneous. If it succeeds the result of MEDCouplingUMesh::tetrahedrize will return a conform mesh.
1781  *
1782  * \return DataArrayInt * - a newly allocated array (to be managed by the caller) containing renumbered cell ids.
1783  *
1784  * \throw If \a this is not a mesh containing only NORM_HEXA8 cells.
1785  * \throw If \a this is not properly allocated.
1786  * \sa MEDCouplingUMesh::tetrahedrize, MEDCouplingUMesh::simplexize.
1787  */
1788 DataArrayInt *MEDCoupling1SGTUMesh::sortHexa8EachOther()
1789 {
1790   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> quads(explodeEachHexa8To6Quad4());//checks that only hexa8
1791   int nbHexa8(getNumberOfCells()),*cQuads(quads->getNodalConnectivity()->getPointer());
1792   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> neighOfQuads(DataArrayInt::New()); neighOfQuads->alloc(nbHexa8*6,1); neighOfQuads->fillWithValue(-1);
1793   int *ptNeigh(neighOfQuads->getPointer());
1794   {//neighOfQuads tells for each face of each Quad8 which cell (if!=-1) is connected to this face.
1795     MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> quadsTmp(quads->buildUnstructured());
1796     MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ccSafe,cciSafe;
1797     DataArrayInt *cc(0),*cci(0);
1798     quadsTmp->findCommonCells(3,0,cc,cci);
1799     ccSafe=cc; cciSafe=cci;
1800     const int *ccPtr(ccSafe->begin()),nbOfPair(cci->getNumberOfTuples()-1);
1801     for(int i=0;i<nbOfPair;i++)
1802       { ptNeigh[ccPtr[2*i+0]]=ccPtr[2*i+1]/6; ptNeigh[ccPtr[2*i+1]]=ccPtr[2*i+0]/6; }
1803   }
1804   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(0,1);
1805   std::vector<bool> fetched(nbHexa8,false);
1806   std::vector<bool>::iterator it(std::find(fetched.begin(),fetched.end(),false));
1807   while(it!=fetched.end())//it will turns as time as number of connected zones
1808     {
1809       int cellId((int)std::distance(fetched.begin(),it));//it is the seed of the connected zone.
1810       std::set<int> s; s.insert(cellId);//s contains already organized.
1811       while(!s.empty())
1812         {
1813           std::set<int> sNext;
1814           for(std::set<int>::const_iterator it0=s.begin();it0!=s.end();it0++)
1815             {
1816               fetched[*it0]=true;
1817               int *myNeighb(ptNeigh+6*(*it0));
1818               for(int i=0;i<6;i++)
1819                 {
1820                   if(myNeighb[i]!=-1 && !fetched[myNeighb[i]])
1821                     {
1822                       std::size_t pos(std::distance(HEXA8_FACE_PAIRS,std::find(HEXA8_FACE_PAIRS,HEXA8_FACE_PAIRS+6,i)));
1823                       std::size_t pos0(pos/2),pos1(pos%2);
1824                       if(!UpdateHexa8Cell(HEXA8_FACE_PAIRS[2*pos0+(pos1+1)%2],*it0,cQuads+6*4*(*it0)+4*i,cQuads+6*4*myNeighb[i],ptNeigh+6*myNeighb[i]))
1825                         ret->pushBackSilent(myNeighb[i]);
1826                       fetched[myNeighb[i]]=true;
1827                       sNext.insert(myNeighb[i]);
1828                     }
1829                 }
1830             }
1831           s=sNext;
1832         }
1833       it=std::find(fetched.begin(),fetched.end(),false);
1834     }
1835   if(!ret->empty())
1836     {
1837       int *conn(getNodalConnectivity()->getPointer());
1838       for(const int *pt=ret->begin();pt!=ret->end();pt++)
1839         {
1840           int cellId(*pt);
1841           conn[8*cellId+0]=cQuads[24*cellId+0]; conn[8*cellId+1]=cQuads[24*cellId+1]; conn[8*cellId+2]=cQuads[24*cellId+2]; conn[8*cellId+3]=cQuads[24*cellId+3];
1842           conn[8*cellId+4]=cQuads[24*cellId+4]; conn[8*cellId+5]=cQuads[24*cellId+7]; conn[8*cellId+6]=cQuads[24*cellId+6]; conn[8*cellId+7]=cQuads[24*cellId+5];
1843         }
1844       declareAsNew();
1845     }
1846   return ret.retn();
1847 }
1848
1849 MEDCoupling1DGTUMesh *MEDCoupling1SGTUMesh::computeDualMesh3D() const
1850 {
1851   static const int DUAL_TETRA_0[36]={
1852     4,1,0, 6,0,3, 7,3,1,
1853     4,0,1, 5,2,0, 8,1,2,
1854     6,3,0, 5,0,2, 9,2,3,
1855     7,1,3, 9,3,2, 8,2,1
1856   };
1857   static const int DUAL_TETRA_1[36]={
1858     8,4,10, 11,5,8, 10,7,11,
1859     9,4,8, 8,5,12, 12,6,9,
1860     10,4,9, 9,6,13, 13,7,10,
1861     12,5,11, 13,6,12, 11,7,13
1862   };
1863   static const int FACEID_NOT_SH_NODE[4]={2,3,1,0};
1864   if(getCellModelEnum()!=INTERP_KERNEL::NORM_TETRA4)
1865     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh3D : only TETRA4 supported !");
1866   checkFullyDefined();
1867   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> thisu(buildUnstructured());
1868   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> revNodArr(DataArrayInt::New()),revNodIArr(DataArrayInt::New());
1869   thisu->getReverseNodalConnectivity(revNodArr,revNodIArr);
1870   const int *revNod(revNodArr->begin()),*revNodI(revNodIArr->begin()),*nodal(_conn->begin());
1871   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> d1Arr(DataArrayInt::New()),di1Arr(DataArrayInt::New()),rd1Arr(DataArrayInt::New()),rdi1Arr(DataArrayInt::New());
1872   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> edges(thisu->explode3DMeshTo1D(d1Arr,di1Arr,rd1Arr,rdi1Arr));
1873   const int *d1(d1Arr->begin());
1874   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> d2Arr(DataArrayInt::New()),di2Arr(DataArrayInt::New()),rd2Arr(DataArrayInt::New()),rdi2Arr(DataArrayInt::New());
1875   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> faces(thisu->buildDescendingConnectivity(d2Arr,di2Arr,rd2Arr,rdi2Arr));  thisu=0;
1876   const int *d2(d2Arr->begin()),*rdi2(rdi2Arr->begin());
1877   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> edgesBaryArr(edges->getBarycenterAndOwner()),facesBaryArr(faces->getBarycenterAndOwner()),baryArr(getBarycenterAndOwner());
1878   const int nbOfNodes(getNumberOfNodes()),offset0(nbOfNodes+faces->getNumberOfCells()),offset1(offset0+edges->getNumberOfCells());
1879   edges=0; faces=0;
1880   std::vector<const DataArrayDouble *> v(4); v[0]=getCoords(); v[1]=facesBaryArr; v[2]=edgesBaryArr; v[3]=baryArr;
1881   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> zeArr(DataArrayDouble::Aggregate(v)); baryArr=0; edgesBaryArr=0; facesBaryArr=0;
1882   std::string name("DualOf_"); name+=getName();
1883   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(MEDCoupling1DGTUMesh::New(name,INTERP_KERNEL::NORM_POLYHED)); ret->setCoords(zeArr);
1884   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cArr(DataArrayInt::New()),ciArr(DataArrayInt::New()); ciArr->alloc(nbOfNodes+1,1); ciArr->setIJ(0,0,0); cArr->alloc(0,1);
1885   for(int i=0;i<nbOfNodes;i++,revNodI++)
1886     {
1887       int nbOfCellsSharingNode(revNodI[1]-revNodI[0]);
1888       if(nbOfCellsSharingNode==0)
1889         {
1890           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeDualMesh3D : Node #" << i << " is orphan !"; 
1891           throw INTERP_KERNEL::Exception(oss.str().c_str());
1892         }
1893       for(int j=0;j<nbOfCellsSharingNode;j++)
1894         {
1895           int curCellId(revNod[revNodI[0]+j]);
1896           const int *connOfCurCell(nodal+4*curCellId);
1897           std::size_t nodePosInCurCell(std::distance(connOfCurCell,std::find(connOfCurCell,connOfCurCell+4,i)));
1898           if(j!=0) cArr->pushBackSilent(-1);
1899           int tmp[14];
1900           //
1901           tmp[0]=d1[6*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+0]-4]+offset0; tmp[1]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+1]]+nbOfNodes;
1902           tmp[2]=curCellId+offset1; tmp[3]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+2]]+nbOfNodes;
1903           tmp[4]=-1;
1904           tmp[5]=d1[6*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+3]-4]+offset0; tmp[6]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+4]]+nbOfNodes;
1905           tmp[7]=curCellId+offset1; tmp[8]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+5]]+nbOfNodes;
1906           tmp[9]=-1;
1907           tmp[10]=d1[6*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+6]-4]+offset0; tmp[11]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+7]]+nbOfNodes;
1908           tmp[12]=curCellId+offset1; tmp[13]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+8]]+nbOfNodes;
1909           cArr->insertAtTheEnd(tmp,tmp+14);
1910           int kk(0);
1911           for(int k=0;k<4;k++)
1912             {
1913               if(FACEID_NOT_SH_NODE[nodePosInCurCell]!=k)
1914                 {
1915                   const int *faceId(d2+4*curCellId+k);
1916                   if(rdi2[*faceId+1]-rdi2[*faceId]==1)
1917                     {
1918                       int tmp2[5]; tmp2[0]=-1; tmp2[1]=i;
1919                       tmp2[2]=d1[6*curCellId+DUAL_TETRA_1[9*nodePosInCurCell+3*kk+0]-8]+offset0;
1920                       tmp2[3]=d2[4*curCellId+DUAL_TETRA_1[9*nodePosInCurCell+3*kk+1]-4]+nbOfNodes;
1921                       tmp2[4]=d1[6*curCellId+DUAL_TETRA_1[9*nodePosInCurCell+3*kk+2]-8]+offset0;
1922                       cArr->insertAtTheEnd(tmp2,tmp2+5);
1923                     }
1924                   kk++;
1925                 }
1926             }
1927         }
1928       ciArr->setIJ(i+1,0,cArr->getNumberOfTuples());
1929     }
1930   ret->setNodalConnectivity(cArr,ciArr);
1931   return ret.retn();
1932 }
1933
1934 MEDCoupling1DGTUMesh *MEDCoupling1SGTUMesh::computeDualMesh2D() const
1935 {
1936   static const int DUAL_TRI_0[6]={0,2, 1,0, 2,1};
1937   static const int DUAL_TRI_1[6]={-3,+5, +3,-4, +4,-5};
1938   static const int FACEID_NOT_SH_NODE[3]={1,2,0};
1939   if(getCellModelEnum()!=INTERP_KERNEL::NORM_TRI3)
1940     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh2D : only TRI3 supported !");
1941   checkFullyDefined();
1942   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> thisu(buildUnstructured());
1943   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> revNodArr(DataArrayInt::New()),revNodIArr(DataArrayInt::New());
1944   thisu->getReverseNodalConnectivity(revNodArr,revNodIArr);
1945   const int *revNod(revNodArr->begin()),*revNodI(revNodIArr->begin()),*nodal(_conn->begin());
1946   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> d2Arr(DataArrayInt::New()),di2Arr(DataArrayInt::New()),rd2Arr(DataArrayInt::New()),rdi2Arr(DataArrayInt::New());
1947   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> edges(thisu->buildDescendingConnectivity(d2Arr,di2Arr,rd2Arr,rdi2Arr));  thisu=0;
1948   const int *d2(d2Arr->begin()),*rdi2(rdi2Arr->begin());
1949   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> edgesBaryArr(edges->getBarycenterAndOwner()),baryArr(getBarycenterAndOwner());
1950   const int nbOfNodes(getNumberOfNodes()),offset0(nbOfNodes+edges->getNumberOfCells());
1951   edges=0;
1952   std::vector<const DataArrayDouble *> v(3); v[0]=getCoords(); v[1]=edgesBaryArr; v[2]=baryArr;
1953   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> zeArr(DataArrayDouble::Aggregate(v)); baryArr=0; edgesBaryArr=0;
1954   std::string name("DualOf_"); name+=getName();
1955   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(MEDCoupling1DGTUMesh::New(name,INTERP_KERNEL::NORM_POLYGON)); ret->setCoords(zeArr);
1956   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cArr(DataArrayInt::New()),ciArr(DataArrayInt::New()); ciArr->alloc(nbOfNodes+1,1); ciArr->setIJ(0,0,0); cArr->alloc(0,1);
1957   for(int i=0;i<nbOfNodes;i++,revNodI++)
1958     {
1959       int nbOfCellsSharingNode(revNodI[1]-revNodI[0]);
1960       if(nbOfCellsSharingNode==0)
1961         {
1962           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeDualMesh2D : Node #" << i << " is orphan !"; 
1963           throw INTERP_KERNEL::Exception(oss.str().c_str());
1964         }
1965       std::vector< std::vector<int> > polyg;
1966       for(int j=0;j<nbOfCellsSharingNode;j++)
1967         {
1968           int curCellId(revNod[revNodI[0]+j]);
1969           const int *connOfCurCell(nodal+3*curCellId);
1970           std::size_t nodePosInCurCell(std::distance(connOfCurCell,std::find(connOfCurCell,connOfCurCell+4,i)));
1971           std::vector<int> locV(3);
1972           locV[0]=d2[3*curCellId+DUAL_TRI_0[2*nodePosInCurCell+0]]+nbOfNodes; locV[1]=curCellId+offset0; locV[2]=d2[3*curCellId+DUAL_TRI_0[2*nodePosInCurCell+1]]+nbOfNodes;
1973           polyg.push_back(locV);
1974           int kk(0);
1975           for(int k=0;k<3;k++)
1976             {
1977               if(FACEID_NOT_SH_NODE[nodePosInCurCell]!=k)
1978                 {
1979                   const int *edgeId(d2+3*curCellId+k);
1980                   if(rdi2[*edgeId+1]-rdi2[*edgeId]==1)
1981                     {
1982                       std::vector<int> locV2(2);
1983                       int zeLocEdgeIdRel(DUAL_TRI_1[2*nodePosInCurCell+kk]);
1984                       if(zeLocEdgeIdRel>0)
1985                         {  locV2[0]=d2[3*curCellId+zeLocEdgeIdRel-3]+nbOfNodes;  locV2[1]=i; }
1986                       else
1987                         {  locV2[0]=i; locV2[1]=d2[3*curCellId-zeLocEdgeIdRel-3]+nbOfNodes; }
1988                       polyg.push_back(locV2);
1989                     }
1990                   kk++;
1991                 }
1992             }
1993         }
1994       std::vector<int> zePolyg(MEDCoupling1DGTUMesh::BuildAPolygonFromParts(polyg));
1995       cArr->insertAtTheEnd(zePolyg.begin(),zePolyg.end());
1996       ciArr->setIJ(i+1,0,cArr->getNumberOfTuples());
1997     }
1998   ret->setNodalConnectivity(cArr,ciArr);
1999   return ret.retn();
2000 }
2001
2002 /*!
2003  * This method aggregate the bbox of each cell and put it into bbox 
2004  *
2005  * \param [in] arcDetEps - a parameter specifying in case of 2D quadratic polygon cell the detection limit between linear and arc circle. (By default 1e-12)
2006  *                         For all other cases this input parameter is ignored.
2007  * \return DataArrayDouble * - newly created object (to be managed by the caller) \a this number of cells tuples and 2*spacedim components.
2008  * 
2009  * \throw If \a this is not fully set (coordinates and connectivity).
2010  * \throw If a cell in \a this has no valid nodeId.
2011  */
2012 DataArrayDouble *MEDCoupling1SGTUMesh::getBoundingBoxForBBTree(double arcDetEps) const
2013 {
2014   int spaceDim(getSpaceDimension()),nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes()),nbOfNodesPerCell(getNumberOfNodesPerCell());
2015   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
2016   double *bbox(ret->getPointer());
2017   for(int i=0;i<nbOfCells*spaceDim;i++)
2018     {
2019       bbox[2*i]=std::numeric_limits<double>::max();
2020       bbox[2*i+1]=-std::numeric_limits<double>::max();
2021     }
2022   const double *coordsPtr(_coords->getConstPointer());
2023   const int *conn(_conn->getConstPointer());
2024   for(int i=0;i<nbOfCells;i++)
2025     {
2026       int kk(0);
2027       for(int j=0;j<nbOfNodesPerCell;j++,conn++)
2028         {
2029           int nodeId(*conn);
2030           if(nodeId>=0 && nodeId<nbOfNodes)
2031             {
2032               for(int k=0;k<spaceDim;k++)
2033                 {
2034                   bbox[2*spaceDim*i+2*k]=std::min(bbox[2*spaceDim*i+2*k],coordsPtr[spaceDim*nodeId+k]);
2035                   bbox[2*spaceDim*i+2*k+1]=std::max(bbox[2*spaceDim*i+2*k+1],coordsPtr[spaceDim*nodeId+k]);
2036                 }
2037               kk++;
2038             }
2039         }
2040       if(kk==0)
2041         {
2042           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getBoundingBoxForBBTree : cell #" << i << " contains no valid nodeId !";
2043           throw INTERP_KERNEL::Exception(oss.str().c_str());
2044         }
2045     }
2046   return ret.retn();
2047 }
2048
2049 //== 
2050
2051 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::New()
2052 {
2053   return new MEDCoupling1DGTUMesh;
2054 }
2055
2056 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::New(const std::string& name, INTERP_KERNEL::NormalizedCellType type)
2057 {
2058   if(type==INTERP_KERNEL::NORM_ERROR)
2059     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::New : NORM_ERROR is not a valid type to be used as base geometric type for a mesh !");
2060   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
2061   if(!cm.isDynamic())
2062     {
2063       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::New : the input geometric type " << cm.getRepr() << " is static ! Only dynamic types are allowed here !";
2064       throw INTERP_KERNEL::Exception(oss.str().c_str());
2065     }
2066   return new MEDCoupling1DGTUMesh(name,cm);
2067 }
2068
2069 MEDCoupling1DGTUMesh::MEDCoupling1DGTUMesh()
2070 {
2071 }
2072
2073 MEDCoupling1DGTUMesh::MEDCoupling1DGTUMesh(const std::string& name, const INTERP_KERNEL::CellModel& cm):MEDCoupling1GTUMesh(name,cm)
2074 {
2075 }
2076
2077 MEDCoupling1DGTUMesh::MEDCoupling1DGTUMesh(const MEDCoupling1DGTUMesh& other, bool recDeepCpy):MEDCoupling1GTUMesh(other,recDeepCpy),_conn(other._conn)
2078 {
2079   if(recDeepCpy)
2080     {
2081       const DataArrayInt *c(other._conn);
2082       if(c)
2083         _conn=c->deepCpy();
2084       c=other._conn_indx;
2085       if(c)
2086         _conn_indx=c->deepCpy();
2087     }
2088 }
2089
2090 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::clone(bool recDeepCpy) const
2091 {
2092   return new MEDCoupling1DGTUMesh(*this,recDeepCpy);
2093 }
2094
2095 /*!
2096  * This method behaves mostly like MEDCoupling1DGTUMesh::deepCpy method, except that only nodal connectivity arrays are deeply copied.
2097  * The coordinates are shared between \a this and the returned instance.
2098  * 
2099  * \return MEDCouplingUMesh * - A new object instance holding the copy of \a this (deep for connectivity, shallow for coordiantes)
2100  * \sa MEDCoupling1DGTUMesh::deepCpy
2101  */
2102 MEDCouplingPointSet *MEDCoupling1DGTUMesh::deepCpyConnectivityOnly() const
2103 {
2104   checkCoherency();
2105   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(clone(false));
2106   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(_conn->deepCpy()),ci(_conn_indx->deepCpy());
2107   ret->setNodalConnectivity(c,ci);
2108   return ret.retn();
2109 }
2110
2111 void MEDCoupling1DGTUMesh::updateTime() const
2112 {
2113   MEDCoupling1GTUMesh::updateTime();
2114   const DataArrayInt *c(_conn);
2115   if(c)
2116     updateTimeWith(*c);
2117   c=_conn_indx;
2118   if(c)
2119     updateTimeWith(*c);
2120 }
2121
2122 std::size_t MEDCoupling1DGTUMesh::getHeapMemorySizeWithoutChildren() const
2123 {
2124   return MEDCoupling1GTUMesh::getHeapMemorySizeWithoutChildren();
2125 }
2126
2127 std::vector<const BigMemoryObject *> MEDCoupling1DGTUMesh::getDirectChildren() const
2128 {
2129   std::vector<const BigMemoryObject *> ret(MEDCoupling1GTUMesh::getDirectChildren());
2130   const DataArrayInt *c(_conn);
2131   if(c)
2132     ret.push_back(c);
2133   c=_conn_indx;
2134   if(c)
2135     ret.push_back(c);
2136   return ret;
2137 }
2138
2139 MEDCouplingMesh *MEDCoupling1DGTUMesh::deepCpy() const
2140 {
2141   return clone(true);
2142 }
2143
2144 bool MEDCoupling1DGTUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
2145 {
2146   if(!other)
2147     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::isEqualIfNotWhy : input other pointer is null !");
2148   std::ostringstream oss; oss.precision(15);
2149   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2150   if(!otherC)
2151     {
2152       reason="mesh given in input is not castable in MEDCoupling1DGTUMesh !";
2153       return false;
2154     }
2155   if(!MEDCoupling1GTUMesh::isEqualIfNotWhy(other,prec,reason))
2156     return false;
2157   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
2158   if(c1==c2)
2159     return true;
2160   if(!c1 || !c2)
2161     {
2162       reason="in connectivity of single dynamic geometric type exactly one among this and other is null !";
2163       return false;
2164     }
2165   if(!c1->isEqualIfNotWhy(*c2,reason))
2166     {
2167       reason.insert(0,"Nodal connectivity DataArrayInt differs : ");
2168       return false;
2169     }
2170   c1=_conn_indx; c2=otherC->_conn_indx;
2171   if(c1==c2)
2172     return true;
2173   if(!c1 || !c2)
2174     {
2175       reason="in connectivity index of single dynamic geometric type exactly one among this and other is null !";
2176       return false;
2177     }
2178   if(!c1->isEqualIfNotWhy(*c2,reason))
2179     {
2180       reason.insert(0,"Nodal connectivity index DataArrayInt differs : ");
2181       return false;
2182     }
2183   return true;
2184 }
2185
2186 bool MEDCoupling1DGTUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
2187 {
2188   if(!other)
2189     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::isEqualWithoutConsideringStr : input other pointer is null !");
2190   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2191   if(!otherC)
2192     return false;
2193   if(!MEDCoupling1GTUMesh::isEqualWithoutConsideringStr(other,prec))
2194     return false;
2195   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
2196   if(c1==c2)
2197     return true;
2198   if(!c1 || !c2)
2199     return false;
2200   if(!c1->isEqualWithoutConsideringStr(*c2))
2201     return false;
2202   return true;
2203   c1=_conn_indx; c2=otherC->_conn_indx;
2204   if(c1==c2)
2205     return true;
2206   if(!c1 || !c2)
2207     return false;
2208   if(!c1->isEqualWithoutConsideringStr(*c2))
2209     return false;
2210   return true;
2211 }
2212
2213 /*!
2214  * Checks if \a this and \a other meshes are geometrically equivalent with high
2215  * probability, else an exception is thrown. The meshes are considered equivalent if
2216  * (1) meshes contain the same number of nodes and the same number of elements of the
2217  * same types (2) three cells of the two meshes (first, last and middle) are based
2218  * on coincident nodes (with a specified precision).
2219  *  \param [in] other - the mesh to compare with.
2220  *  \param [in] prec - the precision used to compare nodes of the two meshes.
2221  *  \throw If the two meshes do not match.
2222  */
2223 void MEDCoupling1DGTUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double prec) const
2224 {
2225   MEDCouplingPointSet::checkFastEquivalWith(other,prec);
2226   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2227   if(!otherC)
2228     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : Two meshes are not unstructured with single dynamic geometric type !");
2229   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
2230   if(c1!=c2)
2231     {
2232       if(!c1 || !c2)
2233         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : presence of nodal connectivity only in one of the 2 meshes !");
2234       if((c1->isAllocated() && !c2->isAllocated()) || (!c1->isAllocated() && c2->isAllocated()))
2235         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity, only one is allocated !");
2236       if(c1->getNumberOfComponents()!=1 || c1->getNumberOfComponents()!=1)
2237         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity, must have 1 and only 1 component !");
2238       if(c1->getHashCode()!=c2->getHashCode())
2239         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : nodal connectivity differs");
2240     }
2241   c1=_conn_indx; c2=otherC->_conn_indx;
2242   if(c1!=c2)
2243     {
2244       if(!c1 || !c2)
2245         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : presence of nodal connectivity index only in one of the 2 meshes !");
2246       if((c1->isAllocated() && !c2->isAllocated()) || (!c1->isAllocated() && c2->isAllocated()))
2247         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity index, only one is allocated !");
2248       if(c1->getNumberOfComponents()!=1 || c1->getNumberOfComponents()!=1)
2249         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity index, must have 1 and only 1 component !");
2250       if(c1->getHashCode()!=c2->getHashCode())
2251         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : nodal connectivity index differs");
2252     }
2253 }
2254
2255 void MEDCoupling1DGTUMesh::checkCoherencyOfConnectivity() const
2256 {
2257   const DataArrayInt *c1(_conn);
2258   if(c1)
2259     {
2260       if(c1->getNumberOfComponents()!=1)
2261         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to be with number of components set to one !");
2262       if(c1->getInfoOnComponent(0)!="")
2263         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to have no info on its single component !");
2264       c1->checkAllocated();
2265     }
2266   else
2267     throw INTERP_KERNEL::Exception("Nodal connectivity array not defined !");
2268   //
2269   int sz2=_conn->getNumberOfTuples();
2270   c1=_conn_indx;
2271   if(c1)
2272     {
2273       if(c1->getNumberOfComponents()!=1)
2274         throw INTERP_KERNEL::Exception("Nodal connectivity index array is expected to be with number of components set to one !");
2275       c1->checkAllocated();
2276       if(c1->getNumberOfTuples()<1)
2277         throw INTERP_KERNEL::Exception("Nodal connectivity index array is expected to have a a size of 1 at least !");
2278       if(c1->getInfoOnComponent(0)!="")
2279         throw INTERP_KERNEL::Exception("Nodal connectivity index array is expected to have no info on its single component !");
2280       int f=c1->front(),ll=c1->back();
2281       if(f<0 || f>=sz2)
2282         {
2283           std::ostringstream oss; oss << "Nodal connectivity index array first value (" << f << ") is expected to be exactly in [0," << sz2 << ") !";
2284           throw INTERP_KERNEL::Exception(oss.str().c_str());
2285         }
2286       if(ll<0 || ll>sz2)
2287         {
2288           std::ostringstream oss; oss << "Nodal connectivity index array last value (" << ll << ") is expected to be exactly in [0," << sz2 << "] !";
2289           throw INTERP_KERNEL::Exception(oss.str().c_str());
2290         }
2291       if(f>ll)
2292         {
2293           std::ostringstream oss; oss << "Nodal connectivity index array looks very bad (not increasing monotonic) because front (" << f << ") is greater that back (" << ll << ") !";
2294           throw INTERP_KERNEL::Exception(oss.str().c_str());
2295         }
2296     }
2297   else
2298     throw INTERP_KERNEL::Exception("Nodal connectivity index array not defined !");
2299   int szOfC1Exp=_conn_indx->back();
2300   if(sz2<szOfC1Exp)
2301     {
2302       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::checkCoherencyOfConnectivity : The expected length of nodal connectivity array regarding index is " << szOfC1Exp << " but the actual size of it is " << c1->getNumberOfTuples() << " !";
2303       throw INTERP_KERNEL::Exception(oss.str().c_str());
2304     }
2305 }
2306
2307 /*!
2308  * If \a this pass this method, you are sure that connectivity arrays are not null, with exactly one component, no name, no component name, allocated.
2309  * In addition you are sure that the length of nodal connectivity index array is bigger than or equal to one.
2310  * In addition you are also sure that length of nodal connectivity is coherent with the content of the last value in the index array.
2311  */
2312 void MEDCoupling1DGTUMesh::checkCoherency() const
2313 {
2314   MEDCouplingPointSet::checkCoherency();
2315   checkCoherencyOfConnectivity();
2316 }
2317
2318 void MEDCoupling1DGTUMesh::checkCoherency1(double eps) const
2319 {
2320   checkCoherency();
2321   const DataArrayInt *c1(_conn),*c2(_conn_indx);
2322   if(!c2->isMonotonic(true))
2323     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkCoherency1 : the nodal connectivity index is expected to be increasing monotinic !");
2324   //
2325   int nbOfTuples=c1->getNumberOfTuples();
2326   int nbOfNodes=getNumberOfNodes();
2327   const int *w(c1->begin());
2328   for(int i=0;i<nbOfTuples;i++,w++)
2329     {
2330       if(*w==-1) continue;
2331       if(*w<0 || *w>=nbOfNodes)
2332         {
2333           std::ostringstream oss; oss << "At pos #" << i << " of nodal connectivity array references to node id #" << *w << " must be in [0," << nbOfNodes << ") !";
2334           throw INTERP_KERNEL::Exception(oss.str().c_str());
2335         }
2336     }
2337 }
2338
2339 void MEDCoupling1DGTUMesh::checkCoherency2(double eps) const
2340 {
2341   checkCoherency1(eps);
2342 }
2343
2344 int MEDCoupling1DGTUMesh::getNumberOfCells() const
2345 {
2346   checkCoherencyOfConnectivity();//do not remove
2347   return _conn_indx->getNumberOfTuples()-1;
2348 }
2349
2350 /*!
2351  * This method returns a newly allocated array containing this->getNumberOfCells() tuples and 1 component.
2352  * For each cell in \b this the number of nodes constituting cell is computed.
2353  * For each polyhedron cell, the sum of the number of nodes of each face constituting polyhedron cell is returned.
2354  * So for pohyhedrons some nodes can be counted several times in the returned result.
2355  * 
2356  * \return a newly allocated array
2357  */
2358 DataArrayInt *MEDCoupling1DGTUMesh::computeNbOfNodesPerCell() const
2359 {
2360   checkCoherency();
2361   _conn_indx->checkMonotonic(true);
2362   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED)
2363     return _conn_indx->deltaShiftIndex();
2364   // for polyhedrons
2365   int nbOfCells=_conn_indx->getNumberOfTuples()-1;
2366   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2367   ret->alloc(nbOfCells,1);
2368   int *retPtr=ret->getPointer();
2369   const int *ci=_conn_indx->begin(),*c=_conn->begin();
2370   for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2371     *retPtr=ci[1]-ci[0]-std::count(c+ci[0],c+ci[1],-1);
2372   return ret.retn();
2373 }
2374
2375 /*!
2376  * This method returns a newly allocated array containing this->getNumberOfCells() tuples and 1 component.
2377  * For each cell in \b this the number of faces constituting (entity of dimension this->getMeshDimension()-1) cell is computed.
2378  * 
2379  * \return a newly allocated array
2380  */
2381 DataArrayInt *MEDCoupling1DGTUMesh::computeNbOfFacesPerCell() const
2382 {
2383   checkCoherency();
2384   _conn_indx->checkMonotonic(true);
2385   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED && getCellModelEnum()!=INTERP_KERNEL::NORM_QPOLYG)
2386     return _conn_indx->deltaShiftIndex();
2387   if(getCellModelEnum()==INTERP_KERNEL::NORM_QPOLYG)
2388     {
2389       MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=_conn_indx->deltaShiftIndex();
2390       ret->applyDivideBy(2);
2391       return ret.retn();
2392     }
2393   // for polyhedrons
2394   int nbOfCells=_conn_indx->getNumberOfTuples()-1;
2395   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2396   ret->alloc(nbOfCells,1);
2397   int *retPtr=ret->getPointer();
2398   const int *ci=_conn_indx->begin(),*c=_conn->begin();
2399   for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2400     *retPtr=std::count(c+ci[0],c+ci[1],-1)+1;
2401   return ret.retn();
2402 }
2403
2404 /*!
2405  * This method computes effective number of nodes per cell. That is to say nodes appearing several times in nodal connectivity of a cell,
2406  * will be counted only once here whereas it will be counted several times in MEDCoupling1DGTUMesh::computeNbOfNodesPerCell method.
2407  *
2408  * \return DataArrayInt * - new object to be deallocated by the caller.
2409  * \sa MEDCoupling1DGTUMesh::computeNbOfNodesPerCell
2410  */
2411 DataArrayInt *MEDCoupling1DGTUMesh::computeEffectiveNbOfNodesPerCell() const
2412 {
2413   checkCoherency();
2414   _conn_indx->checkMonotonic(true);
2415   int nbOfCells(_conn_indx->getNumberOfTuples()-1);
2416   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2417   ret->alloc(nbOfCells,1);
2418   int *retPtr(ret->getPointer());
2419   const int *ci(_conn_indx->begin()),*c(_conn->begin());
2420   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED)
2421     {
2422       for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2423         {
2424           std::set<int> s(c+ci[0],c+ci[1]);
2425           *retPtr=(int)s.size();
2426         }
2427     }
2428   else
2429     {
2430       for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2431         {
2432           std::set<int> s(c+ci[0],c+ci[1]); s.erase(-1);
2433           *retPtr=(int)s.size();
2434         }
2435     }
2436   return ret.retn();
2437 }
2438
2439 void MEDCoupling1DGTUMesh::getNodeIdsOfCell(int cellId, std::vector<int>& conn) const
2440 {
2441   int nbOfCells(getNumberOfCells());//performs checks
2442   if(cellId>=0 && cellId<nbOfCells)
2443     {
2444       int strt=_conn_indx->getIJ(cellId,0),stp=_conn_indx->getIJ(cellId+1,0);
2445       int nbOfNodes=stp-strt;
2446       if(nbOfNodes<0)
2447         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::getNodeIdsOfCell : the index array is invalid ! Should be increasing monotonic !");
2448       conn.resize(nbOfNodes);
2449       std::copy(_conn->begin()+strt,_conn->begin()+stp,conn.begin());
2450     }
2451   else
2452     {
2453       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getNodeIdsOfCell : request for cellId #" << cellId << " must be in [0," << nbOfCells << ") !";
2454       throw INTERP_KERNEL::Exception(oss.str().c_str());
2455     }
2456 }
2457
2458 int MEDCoupling1DGTUMesh::getNumberOfNodesInCell(int cellId) const
2459 {
2460   int nbOfCells(getNumberOfCells());//performs checks
2461   if(cellId>=0 && cellId<nbOfCells)
2462     {
2463       const int *conn(_conn->begin());
2464       int strt=_conn_indx->getIJ(cellId,0),stp=_conn_indx->getIJ(cellId+1,0);
2465       return stp-strt-std::count(conn+strt,conn+stp,-1);
2466     }
2467   else
2468     {
2469       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getNumberOfNodesInCell : request for cellId #" << cellId << " must be in [0," << nbOfCells << ") !";
2470       throw INTERP_KERNEL::Exception(oss.str().c_str());
2471     }
2472 }
2473
2474 std::string MEDCoupling1DGTUMesh::simpleRepr() const
2475 {
2476   static const char msg0[]="No coordinates specified !";
2477   std::ostringstream ret;
2478   ret << "Single dynamic geometic type (" << _cm->getRepr() << ") unstructured mesh with name : \"" << getName() << "\"\n";
2479   ret << "Description of mesh : \"" << getDescription() << "\"\n";
2480   int tmpp1,tmpp2;
2481   double tt=getTime(tmpp1,tmpp2);
2482   ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
2483   ret << "Iteration : " << tmpp1  << " Order : " << tmpp2 << "\n";
2484   ret << "Mesh dimension : " << getMeshDimension() << "\nSpace dimension : ";
2485   if(_coords!=0)
2486     {
2487       const int spaceDim=getSpaceDimension();
2488       ret << spaceDim << "\nInfo attached on space dimension : ";
2489       for(int i=0;i<spaceDim;i++)
2490         ret << "\"" << _coords->getInfoOnComponent(i) << "\" ";
2491       ret << "\n";
2492     }
2493   else
2494     ret << msg0 << "\n";
2495   ret << "Number of nodes : ";
2496   if(_coords!=0)
2497     ret << getNumberOfNodes() << "\n";
2498   else
2499     ret << msg0 << "\n";
2500   ret << "Number of cells : ";
2501   bool isOK=true;
2502   try { checkCoherency(); } catch(INTERP_KERNEL::Exception& /* e */)
2503   {
2504       ret << "Nodal connectivity arrays are not set or badly set !\n";
2505       isOK=false;
2506   }
2507   if(isOK)
2508     ret << getNumberOfCells() << "\n";
2509   ret << "Cell type : " << _cm->getRepr() << "\n";
2510   return ret.str();
2511 }
2512
2513 std::string MEDCoupling1DGTUMesh::advancedRepr() const
2514 {
2515   std::ostringstream ret;
2516   ret << simpleRepr();
2517   ret << "\nCoordinates array : \n___________________\n\n";
2518   if(_coords)
2519     _coords->reprWithoutNameStream(ret);
2520   else
2521     ret << "No array set !\n";
2522   ret << "\n\nNodal Connectivity : \n____________________\n\n";
2523   //
2524   bool isOK=true;
2525   try { checkCoherency1(); } catch(INTERP_KERNEL::Exception& /* e */)
2526   {
2527       ret << "Nodal connectivity arrays are not set or badly set !\n";
2528       isOK=false;
2529   }
2530   if(!isOK)
2531     return ret.str();
2532   int nbOfCells=getNumberOfCells();
2533   const int *ci=_conn_indx->begin(),*c=_conn->begin();
2534   for(int i=0;i<nbOfCells;i++,ci++)
2535     {
2536       ret << "Cell #" << i << " : ";
2537       std::copy(c+ci[0],c+ci[1],std::ostream_iterator<int>(ret," "));
2538       ret << "\n";
2539     }
2540   return ret.str();
2541 }
2542
2543 DataArrayDouble *MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell() const
2544 {
2545   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
2546   int spaceDim=getSpaceDimension();
2547   int nbOfCells=getNumberOfCells();//checkCoherency()
2548   int nbOfNodes=getNumberOfNodes();
2549   ret->alloc(nbOfCells,spaceDim);
2550   double *ptToFill=ret->getPointer();
2551   const double *coor=_coords->begin();
2552   const int *nodal=_conn->begin(),*nodali=_conn_indx->begin();
2553   nodal+=nodali[0];
2554   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED)
2555     {
2556       for(int i=0;i<nbOfCells;i++,ptToFill+=spaceDim,nodali++)
2557         {
2558           std::fill(ptToFill,ptToFill+spaceDim,0.);
2559           if(nodali[0]<nodali[1])// >= to avoid division by 0.
2560             {
2561               for(int j=nodali[0];j<nodali[1];j++,nodal++)
2562                 {
2563                   if(*nodal>=0 && *nodal<nbOfNodes)
2564                     std::transform(coor+spaceDim*nodal[0],coor+spaceDim*(nodal[0]+1),ptToFill,ptToFill,std::plus<double>());
2565                   else
2566                     {
2567                       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell : on cell #" << i << " presence of nodeId #" << *nodal << " should be in [0," <<   nbOfNodes << ") !";
2568                       throw INTERP_KERNEL::Exception(oss.str().c_str());
2569                     }
2570                   std::transform(ptToFill,ptToFill+spaceDim,ptToFill,std::bind2nd(std::multiplies<double>(),1./(nodali[1]-nodali[0])));
2571                 }
2572             }
2573           else
2574             {
2575               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell : at cell #" << i << " the nodal index array is invalid !";
2576               throw INTERP_KERNEL::Exception(oss.str().c_str());
2577             }
2578         }
2579     }
2580   else
2581     {
2582       for(int i=0;i<nbOfCells;i++,ptToFill+=spaceDim,nodali++)
2583         {
2584           std::fill(ptToFill,ptToFill+spaceDim,0.);
2585           if(nodali[0]<nodali[1])// >= to avoid division by 0.
2586             {
2587               int nbOfNod=0;
2588               for(int j=nodali[0];j<nodali[1];j++,nodal++)
2589                 {
2590                   if(*nodal==-1) continue;
2591                   if(*nodal>=0 && *nodal<nbOfNodes)
2592                     {
2593                       std::transform(coor+spaceDim*nodal[0],coor+spaceDim*(nodal[0]+1),ptToFill,ptToFill,std::plus<double>());
2594                       nbOfNod++;
2595                     }
2596                   else
2597                     {
2598                       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell (polyhedron) : on cell #" << i << " presence of nodeId #" << *nodal << " should be in [0," <<   nbOfNodes << ") !";
2599                       throw INTERP_KERNEL::Exception(oss.str().c_str());
2600                     }
2601                 }
2602               if(nbOfNod!=0)
2603                 std::transform(ptToFill,ptToFill+spaceDim,ptToFill,std::bind2nd(std::multiplies<double>(),1./nbOfNod));
2604               else
2605                 {
2606                   std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell (polyhedron) : no nodes in cell #" << i << " !";
2607                   throw INTERP_KERNEL::Exception(oss.str().c_str());
2608                 }
2609             }
2610           else
2611             {
2612               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell (polyhedron)  : at cell #" << i << " the nodal index array is invalid !";
2613               throw INTERP_KERNEL::Exception(oss.str().c_str());
2614             }
2615         }
2616     }
2617   return ret.retn();
2618 }
2619
2620 void MEDCoupling1DGTUMesh::renumberCells(const int *old2NewBg, bool check)
2621 {
2622   int nbCells=getNumberOfCells();
2623   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> o2n=DataArrayInt::New();
2624   o2n->useArray(old2NewBg,false,C_DEALLOC,nbCells,1);
2625   if(check)
2626     o2n=o2n->checkAndPreparePermutation();
2627   //
2628   const int *o2nPtr=o2n->getPointer();
2629   const int *conn=_conn->begin(),*conni=_conn_indx->begin();
2630   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New();
2631   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConnI=DataArrayInt::New();
2632   newConn->alloc(_conn->getNumberOfTuples(),1); newConnI->alloc(nbCells,1);
2633   newConn->copyStringInfoFrom(*_conn); newConnI->copyStringInfoFrom(*_conn_indx);
2634   //
2635   int *newC=newConn->getPointer(),*newCI=newConnI->getPointer();
2636   for(int i=0;i<nbCells;i++)
2637     {
2638       int newPos=o2nPtr[i];
2639       int sz=conni[i+1]-conni[i];
2640       if(sz>=0)
2641         newCI[newPos]=sz;
2642       else
2643         {
2644           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::renumberCells : the index nodal array is invalid for cell #" << i << " !";
2645           throw INTERP_KERNEL::Exception(oss.str().c_str());
2646         }
2647     }
2648   newConnI->computeOffsets2(); newCI=newConnI->getPointer();
2649   //
2650   for(int i=0;i<nbCells;i++,conni++)
2651     {
2652       int newp=o2nPtr[i];
2653       std::copy(conn+conni[0],conn+conni[1],newC+newCI[newp]);
2654     }
2655   _conn=newConn;
2656   _conn_indx=newConnI;
2657 }
2658
2659 MEDCouplingMesh *MEDCoupling1DGTUMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
2660 {
2661   if(other->getType()!=SINGLE_DYNAMIC_GEO_TYPE_UNSTRUCTURED)
2662     throw INTERP_KERNEL::Exception("Merge of umesh only available with umesh single dynamic geo type each other !");
2663   const MEDCoupling1DGTUMesh *otherC=static_cast<const MEDCoupling1DGTUMesh *>(other);
2664   return Merge1DGTUMeshes(this,otherC);
2665 }
2666
2667 MEDCouplingUMesh *MEDCoupling1DGTUMesh::buildUnstructured() const
2668 {
2669   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName(),getMeshDimension());
2670   ret->setCoords(getCoords());
2671   const int *nodalConn=_conn->begin(),*nodalConnI=_conn_indx->begin();
2672   int nbCells=getNumberOfCells();//checkCoherency
2673   int geoType=(int)getCellModelEnum();
2674   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=DataArrayInt::New(); c->alloc(nbCells+_conn->getNumberOfTuples(),1);
2675   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cI=DataArrayInt::New(); cI->alloc(nbCells+1);
2676   int *cPtr=c->getPointer(),*ciPtr=cI->getPointer();
2677   ciPtr[0]=0;
2678   for(int i=0;i<nbCells;i++,ciPtr++)
2679     {
2680       int sz=nodalConnI[i+1]-nodalConnI[i];
2681       if(sz>=0)
2682         {
2683           *cPtr++=geoType;
2684           cPtr=std::copy(nodalConn+nodalConnI[i],nodalConn+nodalConnI[i+1],cPtr);
2685           ciPtr[1]=ciPtr[0]+sz+1;
2686         }
2687       else
2688         {
2689           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::buildUnstructured : Invalid for nodal index for cell #" << i << " !";
2690           throw INTERP_KERNEL::Exception(oss.str().c_str());
2691         }
2692     }
2693   ret->setConnectivity(c,cI,true);
2694   try
2695   { ret->copyTinyInfoFrom(this); }
2696   catch(INTERP_KERNEL::Exception&) { }
2697   return ret.retn();
2698 }
2699
2700 /*!
2701  * Do nothing for the moment, because there is no policy that allows to split polygons, polyhedrons ... into simplexes
2702  */
2703 DataArrayInt *MEDCoupling1DGTUMesh::simplexize(int policy)
2704 {
2705   int nbOfCells=getNumberOfCells();
2706   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2707   ret->alloc(nbOfCells,1);
2708   ret->iota(0);
2709   return ret.retn();
2710 }
2711
2712 void MEDCoupling1DGTUMesh::reprQuickOverview(std::ostream& stream) const
2713 {
2714   stream << "MEDCoupling1DGTUMesh C++ instance at " << this << ". Type=" << _cm->getRepr() << ". Name : \"" << getName() << "\".";
2715   stream << " Mesh dimension : " << getMeshDimension() << ".";
2716   if(!_coords)
2717     { stream << " No coordinates set !"; return ; }
2718   if(!_coords->isAllocated())
2719     { stream << " Coordinates set but not allocated !"; return ; }
2720   stream << " Space dimension : " << _coords->getNumberOfComponents() << "." << std::endl;
2721   stream << "Number of nodes : " << _coords->getNumberOfTuples() << ".";
2722   bool isOK=true;
2723   try { checkCoherency(); } catch(INTERP_KERNEL::Exception&  /* e */)
2724   {
2725       stream << std::endl << "Nodal connectivity NOT set properly !\n";
2726       isOK=false;
2727   }
2728   if(isOK)
2729     stream << std::endl << "Number of cells : " << getNumberOfCells() << ".";
2730 }
2731
2732 void MEDCoupling1DGTUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other)
2733 {
2734   if(!other)
2735     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::shallowCopyConnectivityFrom : input pointer is null !");
2736   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2737   if(!otherC)
2738     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::shallowCopyConnectivityFrom : input pointer is not an MEDCoupling1DGTUMesh instance !");
2739   setNodalConnectivity(otherC->getNodalConnectivity(),otherC->getNodalConnectivityIndex());
2740 }
2741
2742 MEDCouplingPointSet *MEDCoupling1DGTUMesh::mergeMyselfWithOnSameCoords(const MEDCouplingPointSet *other) const
2743 {
2744   if(!other)
2745     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::mergeMyselfWithOnSameCoords : input other is null !");
2746   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2747   if(!otherC)
2748     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::mergeMyselfWithOnSameCoords : the input other mesh is not of type single statuc geo type unstructured !");
2749   std::vector<const MEDCoupling1DGTUMesh *> ms(2);
2750   ms[0]=this;
2751   ms[1]=otherC;
2752   return Merge1DGTUMeshesOnSameCoords(ms);
2753 }
2754
2755 MEDCouplingPointSet *MEDCoupling1DGTUMesh::buildPartOfMySelfKeepCoords(const int *begin, const int *end) const
2756 {
2757   checkCoherency();
2758   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName(),*_cm));
2759   ret->setCoords(_coords);
2760   DataArrayInt *c=0,*ci=0;
2761   MEDCouplingUMesh::ExtractFromIndexedArrays(begin,end,_conn,_conn_indx,c,ci);
2762   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cSafe(c),ciSafe(ci);
2763   ret->setNodalConnectivity(c,ci);
2764   return ret.retn();
2765 }
2766
2767 MEDCouplingPointSet *MEDCoupling1DGTUMesh::buildPartOfMySelfKeepCoords2(int start, int end, int step) const
2768 {
2769   checkCoherency();
2770   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName(),*_cm));
2771   ret->setCoords(_coords);
2772   DataArrayInt *c=0,*ci=0;
2773   MEDCouplingUMesh::ExtractFromIndexedArrays2(start,end,step,_conn,_conn_indx,c,ci);
2774   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cSafe(c),ciSafe(ci);
2775   ret->setNodalConnectivity(c,ci);
2776   return ret.retn();
2777 }
2778
2779 void MEDCoupling1DGTUMesh::computeNodeIdsAlg(std::vector<bool>& nodeIdsInUse) const
2780 {
2781   int sz((int)nodeIdsInUse.size());
2782   int nbCells(getNumberOfCells());
2783   const int *w(_conn->begin()),*wi(_conn_indx->begin());
2784   for(int i=0;i<nbCells;i++,wi++)
2785     for(const int *pt=w+wi[0];pt!=w+wi[1];pt++)
2786       if(*pt!=-1)
2787         {
2788           if(*pt>=0 && *pt<sz)
2789             nodeIdsInUse[*pt]=true;
2790           else
2791             {
2792               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeNodeIdsAlg : At cell #" << i << " presence of node id #" << *pt << " should be in [0," << sz << ") !";
2793               throw INTERP_KERNEL::Exception(oss.str().c_str());
2794             }
2795         }
2796 }
2797
2798 void MEDCoupling1DGTUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const
2799 {
2800   checkFullyDefined();
2801   int nbOfNodes=getNumberOfNodes();
2802   int *revNodalIndxPtr=(int *)malloc((nbOfNodes+1)*sizeof(int));
2803   revNodalIndx->useArray(revNodalIndxPtr,true,C_DEALLOC,nbOfNodes+1,1);
2804   std::fill(revNodalIndxPtr,revNodalIndxPtr+nbOfNodes+1,0);
2805   const int *conn=_conn->begin(),*conni=_conn_indx->begin();
2806   int nbOfCells=getNumberOfCells();
2807   int nbOfEltsInRevNodal=0;
2808   for(int eltId=0;eltId<nbOfCells;eltId++)
2809     {
2810       int nbOfNodesPerCell=conni[eltId+1]-conni[eltId];
2811       if(nbOfNodesPerCell>=0)
2812         {
2813           for(int j=0;j<nbOfNodesPerCell;j++)
2814             {
2815               int nodeId=conn[conni[eltId]+j];
2816               if(nodeId==-1) continue;            
2817               if(nodeId>=0 && nodeId<nbOfNodes)
2818                 {
2819                   nbOfEltsInRevNodal++;
2820                   revNodalIndxPtr[nodeId+1]++;
2821                 }
2822               else
2823                 {
2824                   std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getReverseNodalConnectivity : At cell #" << eltId << " presence of nodeId #" << conn[0] << " should be in [0," << nbOfNodes << ") !";
2825                   throw INTERP_KERNEL::Exception(oss.str().c_str());
2826                 }
2827             }
2828         }
2829       else
2830         {
2831           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getReverseNodalConnectivity : At cell #" << eltId << "nodal connectivity is invalid !";
2832           throw INTERP_KERNEL::Exception(oss.str().c_str());
2833         }
2834     }
2835   std::transform(revNodalIndxPtr+1,revNodalIndxPtr+nbOfNodes+1,revNodalIndxPtr,revNodalIndxPtr+1,std::plus<int>());
2836   conn=_conn->begin();
2837   int *revNodalPtr=(int *)malloc((nbOfEltsInRevNodal)*sizeof(int));
2838   revNodal->useArray(revNodalPtr,true,C_DEALLOC,nbOfEltsInRevNodal,1);
2839   std::fill(revNodalPtr,revNodalPtr+nbOfEltsInRevNodal,-1);
2840   for(int eltId=0;eltId<nbOfCells;eltId++)
2841     {
2842       int nbOfNodesPerCell=conni[eltId+1]-conni[eltId];
2843       for(int j=0;j<nbOfNodesPerCell;j++)
2844         {
2845           int nodeId=conn[conni[eltId]+j];
2846           if(nodeId!=-1)
2847             *std::find_if(revNodalPtr+revNodalIndxPtr[nodeId],revNodalPtr+revNodalIndxPtr[nodeId+1],std::bind2nd(std::equal_to<int>(),-1))=eltId;
2848         }
2849     }
2850 }
2851
2852 void MEDCoupling1DGTUMesh::checkFullyDefined() const
2853 {
2854   if(!((const DataArrayInt *)_conn) || !((const DataArrayInt *)_conn_indx) || !((const DataArrayDouble *)_coords))
2855     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFullyDefined : part of this is not fully defined.");
2856 }
2857
2858 bool MEDCoupling1DGTUMesh::isEmptyMesh(const std::vector<int>& tinyInfo) const
2859 {
2860   throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::isEmptyMesh : not implemented yet !");
2861 }
2862
2863 void MEDCoupling1DGTUMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
2864 {
2865   int it,order;
2866   double time=getTime(it,order);
2867   tinyInfo.clear(); tinyInfoD.clear(); littleStrings.clear();
2868   //
2869   littleStrings.push_back(getName());
2870   littleStrings.push_back(getDescription());
2871   littleStrings.push_back(getTimeUnit());
2872   //
2873   std::vector<std::string> littleStrings2,littleStrings3,littleStrings4;
2874   if((const DataArrayDouble *)_coords)
2875     _coords->getTinySerializationStrInformation(littleStrings2);
2876   if((const DataArrayInt *)_conn)
2877     _conn->getTinySerializationStrInformation(littleStrings3);
2878   if((const DataArrayInt *)_conn_indx)
2879     _conn_indx->getTinySerializationStrInformation(littleStrings4);
2880   int sz0((int)littleStrings2.size()),sz1((int)littleStrings3.size()),sz2((int)littleStrings4.size());
2881   littleStrings.insert(littleStrings.end(),littleStrings2.begin(),littleStrings2.end());
2882   littleStrings.insert(littleStrings.end(),littleStrings3.begin(),littleStrings3.end());
2883   littleStrings.insert(littleStrings.end(),littleStrings4.begin(),littleStrings4.end());
2884   //
2885   tinyInfo.push_back(getCellModelEnum());
2886   tinyInfo.push_back(it);
2887   tinyInfo.push_back(order);
2888   std::vector<int> tinyInfo2,tinyInfo3,tinyInfo4;
2889   if((const DataArrayDouble *)_coords)
2890     _coords->getTinySerializationIntInformation(tinyInfo2);
2891   if((const DataArrayInt *)_conn)
2892     _conn->getTinySerializationIntInformation(tinyInfo3);
2893   if((const DataArrayInt *)_conn_indx)
2894     _conn_indx->getTinySerializationIntInformation(tinyInfo4);
2895   int sz3((int)tinyInfo2.size()),sz4((int)tinyInfo3.size()),sz5((int)tinyInfo4.size());
2896   tinyInfo.push_back(sz0); tinyInfo.push_back(sz1); tinyInfo.push_back(sz2); tinyInfo.push_back(sz3); tinyInfo.push_back(sz4);  tinyInfo.push_back(sz5);
2897   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
2898   tinyInfo.insert(tinyInfo.end(),tinyInfo3.begin(),tinyInfo3.end());
2899   tinyInfo.insert(tinyInfo.end(),tinyInfo4.begin(),tinyInfo4.end());
2900   //
2901   tinyInfoD.push_back(time);
2902 }
2903
2904 void MEDCoupling1DGTUMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
2905 {
2906   std::vector<int> tinyInfo2(tinyInfo.begin()+9,tinyInfo.begin()+9+tinyInfo[6]);
2907   std::vector<int> tinyInfo1(tinyInfo.begin()+9+tinyInfo[6],tinyInfo.begin()+9+tinyInfo[6]+tinyInfo[7]);
2908   std::vector<int> tinyInfo12(tinyInfo.begin()+9+tinyInfo[6]+tinyInfo[7],tinyInfo.begin()+9+tinyInfo[6]+tinyInfo[7]+tinyInfo[8]);
2909   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> p1(DataArrayInt::New()); p1->resizeForUnserialization(tinyInfo1);
2910   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> p2(DataArrayInt::New()); p2->resizeForUnserialization(tinyInfo12);
2911   std::vector<const DataArrayInt *> v(2); v[0]=p1; v[1]=p2;
2912   p2=DataArrayInt::Aggregate(v);
2913   a2->resizeForUnserialization(tinyInfo2);
2914   a1->alloc(p2->getNbOfElems(),1);
2915 }
2916
2917 void MEDCoupling1DGTUMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
2918 {
2919   int sz(0);
2920   if((const DataArrayInt *)_conn)
2921     if(_conn->isAllocated())
2922       sz=_conn->getNbOfElems();
2923   if((const DataArrayInt *)_conn_indx)
2924     if(_conn_indx->isAllocated())
2925       sz+=_conn_indx->getNbOfElems();
2926   a1=DataArrayInt::New();
2927   a1->alloc(sz,1);
2928   int *work(a1->getPointer());
2929   if(sz!=0 && (const DataArrayInt *)_conn)
2930     work=std::copy(_conn->begin(),_conn->end(),a1->getPointer());
2931   if(sz!=0 && (const DataArrayInt *)_conn_indx)
2932     std::copy(_conn_indx->begin(),_conn_indx->end(),work);
2933   sz=0;
2934   if((const DataArrayDouble *)_coords)
2935     if(_coords->isAllocated())
2936       sz=_coords->getNbOfElems();
2937   a2=DataArrayDouble::New();
2938   a2->alloc(sz,1);
2939   if(sz!=0 && (const DataArrayDouble *)_coords)
2940     std::copy(_coords->begin(),_coords->end(),a2->getPointer());
2941 }
2942
2943 void MEDCoupling1DGTUMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
2944                                            const std::vector<std::string>& littleStrings)
2945 {
2946   INTERP_KERNEL::NormalizedCellType gt((INTERP_KERNEL::NormalizedCellType)tinyInfo[0]);
2947   _cm=&INTERP_KERNEL::CellModel::GetCellModel(gt);
2948   setName(littleStrings[0]);
2949   setDescription(littleStrings[1]);
2950   setTimeUnit(littleStrings[2]);
2951   setTime(tinyInfoD[0],tinyInfo[1],tinyInfo[2]);
2952   int sz0(tinyInfo[3]),sz1(tinyInfo[4]),sz2(tinyInfo[5]),sz3(tinyInfo[6]),sz4(tinyInfo[7]),sz5(tinyInfo[8]);
2953   //
2954   _coords=DataArrayDouble::New();
2955   std::vector<int> tinyInfo2(tinyInfo.begin()+9,tinyInfo.begin()+9+sz3);
2956   _coords->resizeForUnserialization(tinyInfo2);
2957   std::copy(a2->begin(),a2->end(),_coords->getPointer());
2958   _conn=DataArrayInt::New();
2959   std::vector<int> tinyInfo3(tinyInfo.begin()+9+sz3,tinyInfo.begin()+9+sz3+sz4);
2960   _conn->resizeForUnserialization(tinyInfo3);
2961   std::copy(a1->begin(),a1->begin()+_conn->getNbOfElems(),_conn->getPointer());
2962   _conn_indx=DataArrayInt::New();
2963   std::vector<int> tinyInfo4(tinyInfo.begin()+9+sz3+sz4,tinyInfo.begin()+9+sz3+sz4+sz5);
2964   _conn_indx->resizeForUnserialization(tinyInfo4);
2965   std::copy(a1->begin()+_conn->getNbOfElems(),a1->end(),_conn_indx->getPointer());
2966   std::vector<std::string> littleStrings2(littleStrings.begin()+3,littleStrings.begin()+3+sz0);
2967   _coords->finishUnserialization(tinyInfo2,littleStrings2);
2968   std::vector<std::string> littleStrings3(littleStrings.begin()+3+sz0,littleStrings.begin()+3+sz0+sz1);
2969   _conn->finishUnserialization(tinyInfo3,littleStrings3);
2970   std::vector<std::string> littleStrings4(littleStrings.begin()+3+sz0+sz1,littleStrings.begin()+3+sz0+sz1+sz2);
2971   _conn_indx->finishUnserialization(tinyInfo4,littleStrings4);
2972 }
2973
2974 /*!
2975  * Finds nodes not used in any cell and returns an array giving a new id to every node
2976  * by excluding the unused nodes, for which the array holds -1. The result array is
2977  * a mapping in "Old to New" mode. 
2978  *  \param [out] nbrOfNodesInUse - number of node ids present in the nodal connectivity.
2979  *  \return DataArrayInt * - a new instance of DataArrayInt. Its length is \a
2980  *          this->getNumberOfNodes(). It holds for each node of \a this mesh either -1
2981  *          if the node is unused or a new id else. The caller is to delete this
2982  *          array using decrRef() as it is no more needed.  
2983  *  \throw If the coordinates array is not set.
2984  *  \throw If the nodal connectivity of cells is not defined.
2985  *  \throw If the nodal connectivity includes an invalid id.
2986  */
2987 DataArrayInt *MEDCoupling1DGTUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const
2988 {
2989   nbrOfNodesInUse=-1;
2990   int nbOfNodes=getNumberOfNodes();
2991   int nbOfCells=getNumberOfCells();//checkCoherency
2992   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2993   ret->alloc(nbOfNodes,1);
2994   int *traducer=ret->getPointer();
2995   std::fill(traducer,traducer+nbOfNodes,-1);
2996   const int *conn=_conn->begin(),*conni(_conn_indx->begin());
2997   for(int i=0;i<nbOfCells;i++,conni++)
2998     {
2999       int nbNodesPerCell=conni[1]-conni[0];
3000       for(int j=0;j<nbNodesPerCell;j++)
3001         {
3002           int nodeId=conn[conni[0]+j];
3003           if(nodeId==-1) continue;
3004           if(nodeId>=0 && nodeId<nbOfNodes)
3005             traducer[nodeId]=1;
3006           else
3007             {
3008               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getNodeIdsInUse : In cell #" << i  << " presence of node id " <<  nodeId << " not in [0," << nbOfNodes << ") !";
3009               throw INTERP_KERNEL::Exception(oss.str().c_str());
3010             }
3011         }
3012     }
3013   nbrOfNodesInUse=(int)std::count(traducer,traducer+nbOfNodes,1);
3014   std::transform(traducer,traducer+nbOfNodes,traducer,MEDCouplingAccVisit());
3015   return ret.retn();
3016 }
3017
3018 /*!
3019  * Changes ids of nodes within the nodal connectivity arrays according to a permutation
3020  * array in "Old to New" mode. The node coordinates array is \b not changed by this method.
3021  * This method is a generalization of shiftNodeNumbersInConn().
3022  *  \warning This method performs no check of validity of new ids. **Use it with care !**
3023  *  \param [in] newNodeNumbersO2N - a permutation array, of length \a
3024  *         this->getNumberOfNodes(), in "Old to New" mode. 
3025  *         See \ref MEDCouplingArrayRenumbering for more info on renumbering modes.
3026  *  \throw If the nodal connectivity of cells is not defined.
3027  */
3028 void MEDCoupling1DGTUMesh::renumberNodesInConn(const int *newNodeNumbersO2N)
3029 {
3030   getNumberOfCells();//only to check that all is well defined.
3031   //
3032   int nbElemsIn=getNumberOfNodes();
3033   int nbOfTuples=_conn->getNumberOfTuples();
3034   int *pt=_conn->getPointer();
3035   for(int i=0;i<nbOfTuples;i++,pt++)
3036     {
3037       if(*pt==-1) continue;
3038       if(*pt>=0 && *pt<nbElemsIn)
3039         *pt=newNodeNumbersO2N[*pt];
3040       else
3041         {
3042           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::renumberNodesInConn : error on tuple #" << i << " value is " << *pt << " and indirectionnal array as a size equal to " << nbElemsIn;
3043           throw INTERP_KERNEL::Exception(oss.str().c_str());
3044         }
3045     }
3046   _conn->declareAsNew();
3047   //
3048   updateTime();
3049 }
3050
3051 /*!
3052  * Keeps from \a this only cells which constituing point id are in the ids specified by [\a begin,\a end).
3053  * The resulting cell ids are stored at the end of the 'cellIdsKept' parameter.
3054  * Parameter \a fullyIn specifies if a cell that has part of its nodes in ids array is kept or not.
3055  * If \a fullyIn is true only cells whose ids are \b fully contained in [\a begin,\a end) tab will be kept.
3056  *
3057  * \param [in] begin input start of array of node ids.
3058  * \param [in] end input end of array of node ids.
3059  * \param [in] fullyIn input that specifies if all node ids must be in [\a begin,\a end) array to consider cell to be in.
3060  * \param [in,out] cellIdsKeptArr array where all candidate cell ids are put at the end.
3061  */
3062 void MEDCoupling1DGTUMesh::fillCellIdsToKeepFromNodeIds(const int *begin, const int *end, bool fullyIn, DataArrayInt *&cellIdsKeptArr) const
3063 {
3064   int nbOfCells=getNumberOfCells();
3065   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellIdsKept=DataArrayInt::New(); cellIdsKept->alloc(0,1);
3066   int tmp=-1;
3067   int sz=_conn->getMaxValue(tmp); sz=std::max(sz,0)+1;
3068   std::vector<bool> fastFinder(sz,false);
3069   for(const int *work=begin;work!=end;work++)
3070     if(*work>=0 && *work<sz)
3071       fastFinder[*work]=true;
3072   const int *conn=_conn->begin(),*conni=_conn_indx->begin();
3073   for(int i=0;i<nbOfCells;i++,conni++)
3074     {
3075       int ref=0,nbOfHit=0;
3076       int nbNodesPerCell=conni[1]-conni[0];
3077       if(nbNodesPerCell>=0)
3078         {
3079           for(int j=0;j<nbNodesPerCell;j++)
3080             {
3081               int nodeId=conn[conni[0]+j];
3082               if(nodeId>=0)
3083                 {
3084                   ref++;
3085                   if(fastFinder[nodeId])
3086                     nbOfHit++;
3087                 }
3088             }
3089         }
3090       else
3091         {
3092           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::fillCellIdsToKeepFromNodeIds : invalid index array for cell #" << i << " !";
3093           throw INTERP_KERNEL::Exception(oss.str().c_str());
3094         }
3095       if((ref==nbOfHit && fullyIn) || (nbOfHit!=0 && !fullyIn))
3096         cellIdsKept->pushBackSilent(i);
3097     }
3098   cellIdsKeptArr=cellIdsKept.retn();
3099 }
3100
3101 void MEDCoupling1DGTUMesh::allocateCells(int nbOfCells)
3102 {
3103   if(nbOfCells<0)
3104     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::allocateCells : the input number of cells should be >= 0 !");
3105   _conn=DataArrayInt::New();
3106   _conn->reserve(nbOfCells*3);
3107   _conn_indx=DataArrayInt::New();
3108   _conn_indx->reserve(nbOfCells+1); _conn_indx->pushBackSilent(0);
3109   declareAsNew();
3110 }
3111
3112 /*!
3113  * Appends at the end of \a this a cell having nodal connectivity array defined in [ \a nodalConnOfCellBg, \a nodalConnOfCellEnd ).
3114  *
3115  * \param [in] nodalConnOfCellBg - the begin (included) of nodal connectivity of the cell to add.
3116  * \param [in] nodalConnOfCellEnd - the end (excluded) of nodal connectivity of the cell to add.
3117  * \throw If the length of the input nodal connectivity array of the cell to add is not equal to number of nodes per cell relative to the unique geometric type
3118  *        attached to \a this.
3119  * \thow If the nodal connectivity array in \a this is null (call MEDCoupling1SGTUMesh::allocateCells before).
3120  */
3121 void MEDCoupling1DGTUMesh::insertNextCell(const int *nodalConnOfCellBg, const int *nodalConnOfCellEnd)
3122 {
3123   int sz=(int)std::distance(nodalConnOfCellBg,nodalConnOfCellEnd);
3124   DataArrayInt *c(_conn),*c2(_conn_indx);
3125   if(c && c2)
3126     {
3127       int pos=c2->back();
3128       if(pos==c->getNumberOfTuples())
3129         {
3130           c->pushBackValsSilent(nodalConnOfCellBg,nodalConnOfCellEnd);
3131           c2->pushBackSilent(pos+sz);
3132         }
3133       else
3134         {
3135           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::insertNextCell : The nodal index array (end=" << pos << ") mismatches with nodal array (length=" << c->getNumberOfTuples() << ") !";
3136           throw INTERP_KERNEL::Exception(oss.str().c_str());
3137         }
3138     }
3139   else
3140     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::insertNextCell : nodal connectivity array is null ! Call MEDCoupling1DGTUMesh::allocateCells before !");
3141 }
3142
3143 void MEDCoupling1DGTUMesh::setNodalConnectivity(DataArrayInt *nodalConn, DataArrayInt *nodalConnIndex)
3144 {
3145   if(nodalConn)
3146     nodalConn->incrRef();
3147   _conn=nodalConn;
3148   if(nodalConnIndex)
3149     nodalConnIndex->incrRef();
3150   _conn_indx=nodalConnIndex;
3151   declareAsNew();
3152 }
3153
3154 /*!
3155  * \return DataArrayInt * - the internal reference to the nodal connectivity. The caller is not reponsible to deallocate it.
3156  */
3157 DataArrayInt *MEDCoupling1DGTUMesh::getNodalConnectivity() const
3158 {
3159   const DataArrayInt *ret(_conn);
3160   return const_cast<DataArrayInt *>(ret);
3161 }
3162
3163 /*!
3164  * \return DataArrayInt * - the internal reference to the nodal connectivity index. The caller is not reponsible to deallocate it.
3165  */
3166 DataArrayInt *MEDCoupling1DGTUMesh::getNodalConnectivityIndex() const
3167 {
3168   const DataArrayInt *ret(_conn_indx);
3169   return const_cast<DataArrayInt *>(ret);
3170 }
3171
3172 /*!
3173  * See the definition of the nodal connectivity pack \ref MEDCoupling1DGTUMesh::isPacked "here".
3174  * This method tries to build a new instance geometrically equivalent to \a this, by limiting at most the number of new object (nodal connectivity).
3175  * Geometrically the returned mesh is equal to \a this. So if \a this is already packed, the return value is a shallow copy of \a this.
3176  *
3177  * Whatever the status of pack of \a this, the coordinates array of the returned newly created instance is the same than those in \a this.
3178  * 
3179  * \param [out] isShallowCpyOfNodalConnn - tells if the returned instance share the same pair of nodal connectivity arrays (true) or if nodal
3180  *              connectivity arrays are different (false)
3181  * \return a new object to be managed by the caller.
3182  * 
3183  * \sa MEDCoupling1DGTUMesh::retrievePackedNodalConnectivity, MEDCoupling1DGTUMesh::isPacked
3184  */
3185 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::copyWithNodalConnectivityPacked(bool& isShallowCpyOfNodalConnn) const
3186 {
3187   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName(),*_cm));
3188   DataArrayInt *nc=0,*nci=0;
3189   isShallowCpyOfNodalConnn=retrievePackedNodalConnectivity(nc,nci);
3190   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ncs(nc),ncis(nci);
3191   ret->_conn=ncs; ret->_conn_indx=ncis;
3192   ret->setCoords(getCoords());
3193   return ret.retn();
3194 }
3195
3196 /*!
3197  * This method allows to compute, if needed, the packed nodal connectivity pair.
3198  * Indeed, it is possible to store in \a this a nodal connectivity array bigger than ranges convered by nodal connectivity index array.
3199  * It is typically the case when nodalConnIndx starts with an id greater than 0, and finishes with id less than number of tuples in \c this->_conn.
3200  * 
3201  * If \a this looks packed (the front of nodal connectivity index equal to 0 and back of connectivity index equal to number of tuple of nodal connectivity array)
3202  * true will be returned and respectively \a this->_conn and \a this->_conn_indx (with ref counter incremented). This is the classical case.
3203  *
3204  * If nodal connectivity index points to a subpart of nodal connectivity index the packed pair of arrays will be computed (new objects) and returned and false
3205  * will be returned.
3206  * 
3207  * This method return 3 elements.
3208  * \param [out] nodalConn - a pointer that can be equal to \a this->_conn if true is returned (general case). Whatever the value of return parameter
3209  *                          this pointer can be seen as a new object, that is to managed by the caller.
3210  * \param [out] nodalConnIndx - a pointer that can be equal to \a this->_conn_indx if true is returned (general case). Whatever the value of return parameter
3211  *                              this pointer can be seen as a new object, that is to managed by the caller.
3212  * \return bool - an indication of the content of the 2 output parameters. If true, \a this looks packed (general case), if true, \a this is not packed then
3213  * output parameters are newly created objects.
3214  *
3215  * \throw if \a this does not pass MEDCoupling1DGTUMesh::checkCoherency test
3216  */
3217 bool MEDCoupling1DGTUMesh::retrievePackedNodalConnectivity(DataArrayInt *&nodalConn, DataArrayInt *&nodalConnIndx) const
3218 {
3219   if(isPacked())//performs the checkCoherency
3220     {
3221       const DataArrayInt *c0(_conn),*c1(_conn_indx);
3222       nodalConn=const_cast<DataArrayInt *>(c0); nodalConnIndx=const_cast<DataArrayInt *>(c1);
3223       nodalConn->incrRef(); nodalConnIndx->incrRef();
3224       return true;
3225     }
3226   int bg=_conn_indx->front(),end=_conn_indx->back();
3227   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nc(_conn->selectByTupleId2(bg,end,1));
3228   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nci(_conn_indx->deepCpy());
3229   nci->applyLin(1,-bg);
3230   nodalConn=nc.retn(); nodalConnIndx=nci.retn();
3231   return false;
3232 }
3233
3234 /*
3235  * If \a this looks packed (the front of nodal connectivity index equal to 0 and back of connectivity index equal to number of tuple of nodal connectivity array)
3236  * true will be returned and respectively \a this->_conn and \a this->_conn_indx (with ref counter incremented). This is the classical case.
3237  * If nodal connectivity index points to a subpart of nodal connectivity index false will be returned.
3238  * \return bool - true if \a this looks packed, false is not.
3239  *
3240  * \throw if \a this does not pass MEDCoupling1DGTUMesh::checkCoherency test
3241  */
3242 bool MEDCoupling1DGTUMesh::isPacked() const
3243 {
3244   checkCoherency();
3245   return _conn_indx->front()==0 && _conn_indx->back()==_conn->getNumberOfTuples();
3246 }
3247
3248 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshes(const MEDCoupling1DGTUMesh *mesh1, const MEDCoupling1DGTUMesh *mesh2)
3249 {
3250   std::vector<const MEDCoupling1DGTUMesh *> tmp(2);
3251   tmp[0]=const_cast<MEDCoupling1DGTUMesh *>(mesh1); tmp[1]=const_cast<MEDCoupling1DGTUMesh *>(mesh2);
3252   return Merge1DGTUMeshes(tmp);
3253 }
3254
3255 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshes(std::vector<const MEDCoupling1DGTUMesh *>& a)
3256 {
3257   std::size_t sz=a.size();
3258   if(sz==0)
3259     return Merge1DGTUMeshesLL(a);
3260   for(std::size_t ii=0;ii<sz;ii++)
3261     if(!a[ii])
3262       {
3263         std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::Merge1DGTUMeshes : item #" << ii << " in input array of size "<< sz << " is empty !";
3264         throw INTERP_KERNEL::Exception(oss.str().c_str());
3265       }
3266   const INTERP_KERNEL::CellModel *cm=&(a[0]->getCellModel());
3267   for(std::size_t ii=0;ii<sz;ii++)
3268     if(&(a[ii]->getCellModel())!=cm)
3269       throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshes : all items must have the same geo type !");
3270   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> > bb(sz);
3271   std::vector< const MEDCoupling1DGTUMesh * > aa(sz);
3272   int spaceDim=-3;
3273   for(std::size_t i=0;i<sz && spaceDim==-3;i++)
3274     {
3275       const MEDCoupling1DGTUMesh *cur=a[i];
3276       const DataArrayDouble *coo=cur->getCoords();
3277       if(coo)
3278         spaceDim=coo->getNumberOfComponents();
3279     }
3280   if(spaceDim==-3)
3281     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshes : no spaceDim specified ! unable to perform merge !");
3282   for(std::size_t i=0;i<sz;i++)
3283     {
3284       bb[i]=a[i]->buildSetInstanceFromThis(spaceDim);
3285       aa[i]=bb[i];
3286     }
3287   return Merge1DGTUMeshesLL(aa);
3288 }
3289
3290 /*!
3291  * \throw If presence of a null instance in the input vector \a a.
3292  * \throw If a is empty
3293  */
3294 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords(std::vector<const MEDCoupling1DGTUMesh *>& a)
3295 {
3296   if(a.empty())
3297     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : input array must be NON EMPTY !");
3298   std::vector<const MEDCoupling1DGTUMesh *>::const_iterator it=a.begin();
3299   if(!(*it))
3300     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : null instance in the first element of input vector !");
3301   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> > objs(a.size());
3302   std::vector<const DataArrayInt *> ncs(a.size()),ncis(a.size());
3303   (*it)->getNumberOfCells();//to check that all is OK
3304   const DataArrayDouble *coords=(*it)->getCoords();
3305   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
3306   bool tmp;
3307   objs[0]=(*it)->copyWithNodalConnectivityPacked(tmp);
3308   ncs[0]=objs[0]->getNodalConnectivity(); ncis[0]=objs[0]->getNodalConnectivityIndex();
3309   it++;
3310   for(int i=1;it!=a.end();i++,it++)
3311     {
3312       if(!(*it))
3313         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : presence of null instance !");
3314       if(cm!=&((*it)->getCellModel()))
3315         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1DGTUMeshes impossible !");
3316       (*it)->getNumberOfCells();//to check that all is OK
3317       objs[i]=(*it)->copyWithNodalConnectivityPacked(tmp);
3318       ncs[i]=objs[i]->getNodalConnectivity(); ncis[i]=objs[i]->getNodalConnectivityIndex();
3319       if(coords!=(*it)->getCoords())
3320         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : not lying on same coords !");
3321     }
3322   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh("merge",*cm));
3323   ret->setCoords(coords);
3324   ret->_conn=DataArrayInt::Aggregate(ncs);
3325   ret->_conn_indx=DataArrayInt::AggregateIndexes(ncis);
3326   return ret.retn();
3327 }
3328
3329 /*!
3330  * Assume that all instances in \a a are non null. If null it leads to a crash. That's why this method is assigned to be low level (LL)
3331  */
3332 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshesLL(std::vector<const MEDCoupling1DGTUMesh *>& a)
3333 {
3334   if(a.empty())
3335     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshes : input array must be NON EMPTY !");
3336   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> > objs(a.size());
3337   std::vector<const DataArrayInt *> ncs(a.size()),ncis(a.size());
3338   std::vector<const MEDCoupling1DGTUMesh *>::const_iterator it=a.begin();
3339   std::vector<int> nbNodesPerElt(a.size());
3340   int nbOfCells=(*it)->getNumberOfCells();
3341   bool tmp;
3342   objs[0]=(*it)->copyWithNodalConnectivityPacked(tmp);
3343   ncs[0]=objs[0]->getNodalConnectivity(); ncis[0]=objs[0]->getNodalConnectivityIndex();
3344   nbNodesPerElt[0]=0;
3345   int prevNbOfNodes=(*it)->getNumberOfNodes();
3346   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
3347   it++;
3348   for(int i=1;it!=a.end();i++,it++)
3349     {
3350       if(cm!=&((*it)->getCellModel()))
3351         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1DGTUMeshes impossible !");
3352       objs[i]=(*it)->copyWithNodalConnectivityPacked(tmp);
3353       ncs[i]=objs[i]->getNodalConnectivity(); ncis[i]=objs[i]->getNodalConnectivityIndex();
3354       nbOfCells+=(*it)->getNumberOfCells();
3355       nbNodesPerElt[i]=nbNodesPerElt[i-1]+prevNbOfNodes;
3356       prevNbOfNodes=(*it)->getNumberOfNodes();
3357     }
3358   std::vector<const MEDCouplingPointSet *> aps(a.size());
3359   std::copy(a.begin(),a.end(),aps.begin());
3360   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> pts=MergeNodesArray(aps);
3361   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh("merge",*cm));
3362   ret->setCoords(pts);
3363   ret->_conn=AggregateNodalConnAndShiftNodeIds(ncs,nbNodesPerElt);
3364   ret->_conn_indx=DataArrayInt::AggregateIndexes(ncis);
3365   return ret.retn();
3366 }
3367
3368 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::buildSetInstanceFromThis(int spaceDim) const
3369 {
3370   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName(),*_cm));
3371   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp1,tmp2;
3372   const DataArrayInt *nodalConn(_conn),*nodalConnI(_conn_indx);
3373   if(!nodalConn)
3374     {
3375       tmp1=DataArrayInt::New(); tmp1->alloc(0,1);
3376     }
3377   else
3378     tmp1=_conn;
3379   ret->_conn=tmp1;
3380   //
3381   if(!nodalConnI)
3382     {
3383       tmp2=DataArrayInt::New(); tmp2->alloc(1,1); tmp2->setIJ(0,0,0);
3384     }
3385   else
3386     tmp2=_conn_indx;
3387   ret->_conn_indx=tmp2;
3388   //
3389   if(!_coords)
3390     {
3391       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords=DataArrayDouble::New(); coords->alloc(0,spaceDim);
3392       ret->setCoords(coords);
3393     }
3394   else
3395     ret->setCoords(_coords);
3396   return ret.retn();
3397 }
3398
3399 /*!
3400  * This method aggregate the bbox of each cell and put it into bbox parameter.
3401  * 
3402  * \param [in] arcDetEps - a parameter specifying in case of 2D quadratic polygon cell the detection limit between linear and arc circle. (By default 1e-12)
3403  *                         For all other cases this input parameter is ignored.
3404  * \return DataArrayDouble * - newly created object (to be managed by the caller) \a this number of cells tuples and 2*spacedim components.
3405  * 
3406  * \throw If \a this is not fully set (coordinates and connectivity).
3407  * \throw If a cell in \a this has no valid nodeId.
3408  */
3409 DataArrayDouble *MEDCoupling1DGTUMesh::getBoundingBoxForBBTree(double arcDetEps) const
3410 {
3411   checkFullyDefined();
3412   int spaceDim(getSpaceDimension()),nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes());
3413   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
3414   double *bbox(ret->getPointer());
3415   for(int i=0;i<nbOfCells*spaceDim;i++)
3416     {
3417       bbox[2*i]=std::numeric_limits<double>::max();
3418       bbox[2*i+1]=-std::numeric_limits<double>::max();
3419     }
3420   const double *coordsPtr(_coords->getConstPointer());
3421   const int *conn(_conn->getConstPointer()),*connI(_conn_indx->getConstPointer());
3422   for(int i=0;i<nbOfCells;i++)
3423     {
3424       int offset=connI[i];
3425       int nbOfNodesForCell(connI[i+1]-offset),kk(0);
3426       for(int j=0;j<nbOfNodesForCell;j++)
3427         {
3428           int nodeId=conn[offset+j];
3429           if(nodeId>=0 && nodeId<nbOfNodes)
3430             {
3431               for(int k=0;k<spaceDim;k++)
3432                 {
3433                   bbox[2*spaceDim*i+2*k]=std::min(bbox[2*spaceDim*i+2*k],coordsPtr[spaceDim*nodeId+k]);
3434                   bbox[2*spaceDim*i+2*k+1]=std::max(bbox[2*spaceDim*i+2*k+1],coordsPtr[spaceDim*nodeId+k]);
3435                 }
3436               kk++;
3437             }
3438         }
3439       if(kk==0)
3440         {
3441           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getBoundingBoxForBBTree : cell #" << i << " contains no valid nodeId !";
3442           throw INTERP_KERNEL::Exception(oss.str().c_str());
3443         }
3444     }
3445   return ret.retn();
3446 }
3447
3448 std::vector<int> MEDCoupling1DGTUMesh::BuildAPolygonFromParts(const std::vector< std::vector<int> >& parts)
3449 {
3450   std::vector<int> ret;
3451   if(parts.empty())
3452     return ret;
3453   ret.insert(ret.end(),parts[0].begin(),parts[0].end());
3454   int ref(ret.back());
3455   std::size_t sz(parts.size()),nbh(1);
3456   std::vector<bool> b(sz,true); b[0]=false;
3457   while(nbh<sz)
3458     {
3459       std::size_t i(0);
3460       for(;i<sz;i++) if(b[i] && parts[i].front()==ref) { ret.insert(ret.end(),parts[i].begin()+1,parts[i].end()); nbh++; break; }
3461       if(i<sz)
3462         ref=ret.back();
3463       else
3464         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::BuildAPolygonFromParts : the input vector is not a part of a single polygon !");
3465     }
3466   if(ret.back()==ret.front())
3467     ret.pop_back();
3468   return ret;
3469 }
3470
3471 /*!
3472  * This method performs an aggregation of \a nodalConns (as DataArrayInt::Aggregate does) but in addition of that a shift is applied on the 
3473  * values contained in \a nodalConns using corresponding offset specified in input \a offsetInNodeIdsPerElt.
3474  * But it also manage the values -1, that have a semantic in MEDCoupling1DGTUMesh class (separator for polyhedron).
3475  *
3476  * \param [in] nodalConns - a list of nodal connectivity arrays same size than \a offsetInNodeIdsPerElt.
3477  * \param [in] offsetInNodeIdsPerElt - a list of offsets to apply.
3478  * \return DataArrayInt * - A new object (to be managed by the caller) that is the result of the aggregation.
3479  * \throw If \a nodalConns or \a offsetInNodeIdsPerElt are empty.
3480  * \throw If \a nodalConns and \a offsetInNodeIdsPerElt have not the same size.
3481  * \throw If presence of null pointer in \a nodalConns.
3482  * \throw If presence of not allocated or array with not exactly one component in \a nodalConns.
3483  */
3484 DataArrayInt *MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds(const std::vector<const DataArrayInt *>& nodalConns, const std::vector<int>& offsetInNodeIdsPerElt)
3485 {
3486   std::size_t sz1(nodalConns.size()),sz2(offsetInNodeIdsPerElt.size());
3487   if(sz1!=sz2)
3488     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : input vectors do not have the same size !");
3489   if(sz1==0)
3490     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : empty vectors in input !");
3491   int nbOfTuples=0;
3492   for(std::vector<const DataArrayInt *>::const_iterator it=nodalConns.begin();it!=nodalConns.end();it++)
3493     {
3494       if(!(*it))
3495         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : presence of null pointer in input vector !");
3496       if(!(*it)->isAllocated())
3497         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : presence of non allocated array in input vector !");
3498       if((*it)->getNumberOfComponents()!=1)
3499         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : presence of array with not exactly one component !");
3500       nbOfTuples+=(*it)->getNumberOfTuples();
3501     }
3502   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(nbOfTuples,1);
3503   int *pt=ret->getPointer();
3504   int i=0;
3505   for(std::vector<const DataArrayInt *>::const_iterator it=nodalConns.begin();it!=nodalConns.end();it++,i++)
3506     {
3507       int curNbt=(*it)->getNumberOfTuples();
3508       const int *inPt=(*it)->begin();
3509       int offset=offsetInNodeIdsPerElt[i];
3510       for(int j=0;j<curNbt;j++,pt++)
3511         {
3512           if(inPt[j]!=-1)
3513             *pt=inPt[j]+offset;
3514           else
3515             *pt=-1;
3516         }
3517     }
3518   return ret.retn();
3519 }
3520
3521 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::New(const MEDCouplingUMesh *m)
3522 {
3523   if(!m)
3524     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::New : input mesh is null !");
3525   std::set<INTERP_KERNEL::NormalizedCellType> gts(m->getAllGeoTypes());
3526   if(gts.size()!=1)
3527     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::New : input mesh must have exactly one geometric type !");
3528   int geoType((int)*gts.begin());
3529   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(MEDCoupling1DGTUMesh::New(m->getName(),*gts.begin()));
3530   ret->setCoords(m->getCoords()); ret->setDescription(m->getDescription());
3531   int nbCells(m->getNumberOfCells());
3532   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()),connI(DataArrayInt::New());
3533   conn->alloc(m->getMeshLength()-nbCells,1); connI->alloc(nbCells+1,1);
3534   int *c(conn->getPointer()),*ci(connI->getPointer()); *ci=0;
3535   const int *cin(m->getNodalConnectivity()->begin()),*ciin(m->getNodalConnectivityIndex()->begin());
3536   for(int i=0;i<nbCells;i++,ciin++,ci++)
3537     {
3538       if(cin[ciin[0]]==geoType)
3539         {
3540           if(ciin[1]-ciin[0]>=1)
3541             {
3542               c=std::copy(cin+ciin[0]+1,cin+ciin[1],c);
3543               ci[1]=ci[0]+ciin[1]-ciin[0]-1;
3544             }
3545           else
3546             {
3547               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The size of cell is not >=0 !";
3548               throw INTERP_KERNEL::Exception(oss.str().c_str());
3549             }
3550         }
3551       else
3552         {
3553           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The geometric type is not those expected !";
3554           throw INTERP_KERNEL::Exception(oss.str().c_str());
3555         }
3556     }
3557   ret->setNodalConnectivity(conn,connI);
3558   return ret.retn();
3559 }