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