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