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