]> SALOME platform Git repositories - modules/med.git/blob - medtool/src/MEDCoupling/MEDCouplingIMesh.cxx
Salome HOME
Move medtool folder to MED base repository
[modules/med.git] / medtool / src / MEDCoupling / MEDCouplingIMesh.cxx
1 // Copyright (C) 2007-2015  CEA/DEN, EDF R&D
2 //
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
7 //
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
11 // Lesser General Public License for more details.
12 //
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
16 //
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
18 //
19 // Author : Anthony Geay (CEA/DEN)
20
21 #include "MEDCouplingIMesh.hxx"
22 #include "MEDCouplingCMesh.hxx"
23 #include "MEDCouplingMemArray.hxx"
24 #include "MEDCouplingFieldDouble.hxx"
25
26 #include <functional>
27 #include <algorithm>
28 #include <sstream>
29 #include <numeric>
30
31 using namespace ParaMEDMEM;
32
33 MEDCouplingIMesh::MEDCouplingIMesh():_space_dim(-1)
34 {
35   _origin[0]=0.; _origin[1]=0.; _origin[2]=0.;
36   _dxyz[0]=0.; _dxyz[1]=0.; _dxyz[2]=0.;
37   _structure[0]=0; _structure[1]=0; _structure[2]=0;
38 }
39
40 MEDCouplingIMesh::MEDCouplingIMesh(const MEDCouplingIMesh& other, bool deepCopy):MEDCouplingStructuredMesh(other,deepCopy),_space_dim(other._space_dim),_axis_unit(other._axis_unit)
41 {
42   _origin[0]=other._origin[0]; _origin[1]=other._origin[1]; _origin[2]=other._origin[2];
43   _dxyz[0]=other._dxyz[0]; _dxyz[1]=other._dxyz[1]; _dxyz[2]=other._dxyz[2];
44   _structure[0]=other._structure[0]; _structure[1]=other._structure[1]; _structure[2]=other._structure[2];
45 }
46
47 MEDCouplingIMesh::~MEDCouplingIMesh()
48 {
49 }
50
51 MEDCouplingIMesh *MEDCouplingIMesh::New()
52 {
53   return new MEDCouplingIMesh;
54 }
55
56 MEDCouplingIMesh *MEDCouplingIMesh::New(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
57                                         const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
58 {
59   MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(new MEDCouplingIMesh);
60   ret->setName(meshName);
61   ret->setSpaceDimension(spaceDim);
62   ret->setNodeStruct(nodeStrctStart,nodeStrctStop);
63   ret->setOrigin(originStart,originStop);
64   ret->setDXYZ(dxyzStart,dxyzStop);
65   return ret.retn();
66 }
67
68 MEDCouplingMesh *MEDCouplingIMesh::deepCpy() const
69 {
70   return clone(true);
71 }
72
73 MEDCouplingIMesh *MEDCouplingIMesh::clone(bool recDeepCpy) const
74 {
75   return new MEDCouplingIMesh(*this,recDeepCpy);
76 }
77
78 /*!
79  * This method creates a copy of \a this enlarged by \a ghostLev cells on each axis.
80  * If \a ghostLev equal to 0 this method behaves as MEDCouplingIMesh::clone.
81  *
82  * \param [in] ghostLev - the ghost level expected
83  * \return MEDCouplingIMesh * - a newly alloacted object to be managed by the caller.
84  * \throw if \a ghostLev < 0.
85  */
86 MEDCouplingIMesh *MEDCouplingIMesh::buildWithGhost(int ghostLev) const
87 {
88   if(ghostLev<0)
89     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::buildWithGhost : the ghostLev must be >= 0 !");
90   checkCoherency();
91   int spaceDim(getSpaceDimension());
92   double origin[3],dxyz[3];
93   int structure[3];
94   for(int i=0;i<spaceDim;i++)
95     {
96       origin[i]=_origin[i]-ghostLev*_dxyz[i];
97       dxyz[i]=_dxyz[i];
98       structure[i]=_structure[i]+2*ghostLev;
99     }
100   MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(MEDCouplingIMesh::New(getName(),spaceDim,structure,structure+spaceDim,origin,origin+spaceDim,dxyz,dxyz+spaceDim));
101   ret->copyTinyInfoFrom(this);
102   return ret.retn();
103 }
104
105 void MEDCouplingIMesh::setNodeStruct(const int *nodeStrctStart, const int *nodeStrctStop)
106 {
107   checkSpaceDimension();
108   int sz((int)std::distance(nodeStrctStart,nodeStrctStop));
109   if(sz!=_space_dim)
110     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setNodeStruct : input vector of node structure has not the right size ! Or change space dimension before calling it !");
111   std::copy(nodeStrctStart,nodeStrctStop,_structure);
112   declareAsNew();
113 }
114
115 std::vector<int> MEDCouplingIMesh::getNodeStruct() const
116 {
117   checkSpaceDimension();
118   return std::vector<int>(_structure,_structure+_space_dim);
119 }
120
121 void MEDCouplingIMesh::setOrigin(const double *originStart, const double *originStop)
122 {
123   checkSpaceDimension();
124   int sz((int)std::distance(originStart,originStop));
125   if(sz!=_space_dim)
126     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setOrigin : input vector of origin vector has not the right size ! Or change space dimension before calling it !");
127   std::copy(originStart,originStop,_origin);
128   declareAsNew();
129 }
130
131 std::vector<double> MEDCouplingIMesh::getOrigin() const
132 {
133   checkSpaceDimension();
134   return std::vector<double>(_origin,_origin+_space_dim);
135 }
136
137 void MEDCouplingIMesh::setDXYZ(const double *dxyzStart, const double *dxyzStop)
138 {
139   checkSpaceDimension();
140   int sz((int)std::distance(dxyzStart,dxyzStop));
141   if(sz!=_space_dim)
142     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setDXYZ : input vector of dxyz vector has not the right size ! Or change space dimension before calling it !");
143   std::copy(dxyzStart,dxyzStop,_dxyz);
144   declareAsNew();
145 }
146
147 std::vector<double> MEDCouplingIMesh::getDXYZ() const
148 {
149   checkSpaceDimension();
150   return std::vector<double>(_dxyz,_dxyz+_space_dim);
151 }
152
153 void MEDCouplingIMesh::setAxisUnit(const std::string& unitName)
154 {
155   _axis_unit=unitName;
156   declareAsNew();
157 }
158
159 std::string MEDCouplingIMesh::getAxisUnit() const
160 {
161   return _axis_unit;
162 }
163
164 /*!
165  * This method returns the measure of any cell in \a this.
166  * This specific method of image grid mesh utilizes the fact that any cell in \a this have the same measure.
167  * The value returned by this method is those used to feed the returned field in the MEDCouplingIMesh::getMeasureField.
168  *
169  * \sa getMeasureField
170  */
171 double MEDCouplingIMesh::getMeasureOfAnyCell() const
172 {
173   checkCoherency();
174   int dim(getSpaceDimension());
175   double ret(1.);
176   for(int i=0;i<dim;i++)
177     ret*=fabs(_dxyz[i]);
178   return ret;
179 }
180
181 /*!
182  * This method is allows to convert \a this into MEDCouplingCMesh instance.
183  * This method is the middle level between MEDCouplingIMesh and the most general MEDCouplingUMesh.
184  * This method is useful for MED writers that do not have still the image grid support.
185  *
186  * \sa MEDCouplingMesh::buildUnstructured
187  */
188 MEDCouplingCMesh *MEDCouplingIMesh::convertToCartesian() const
189 {
190   checkCoherency();
191   MEDCouplingAutoRefCountObjectPtr<MEDCouplingCMesh> ret(MEDCouplingCMesh::New());
192   try
193   { ret->copyTinyInfoFrom(this); }
194   catch(INTERP_KERNEL::Exception& ) { }
195   int spaceDim(getSpaceDimension());
196   std::vector<std::string> infos(buildInfoOnComponents());
197   for(int i=0;i<spaceDim;i++)
198     {
199       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr(DataArrayDouble::New()); arr->alloc(_structure[i],1); arr->setInfoOnComponent(0,infos[i]);
200       arr->iota(); arr->applyLin(_dxyz[i],_origin[i]);
201       ret->setCoordsAt(i,arr);
202     }
203   return ret.retn();
204 }
205
206 /*!
207  * This method refines \a this uniformaly along all of its dimensions. In case of success the space covered by \a this will remain
208  * the same before the invocation except that the number of cells will be multiplied by \a factor ^ this->getMeshDimension().
209  * The origin of \a this will be not touched only spacing and node structure will be changed.
210  * This method can be useful for AMR users.
211  */
212 void MEDCouplingIMesh::refineWithFactor(const std::vector<int>& factors)
213 {
214   if((int)factors.size()!=_space_dim)
215     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::refineWithFactor : refinement factors must have size equal to spaceDim !");
216   checkCoherency();
217   std::vector<int> structure(_structure,_structure+3);
218   std::vector<double> dxyz(_dxyz,_dxyz+3);
219   for(int i=0;i<_space_dim;i++)
220     {
221       if(factors[i]<=0)
222         {
223           std::ostringstream oss; oss << "MEDCouplingIMesh::refineWithFactor : factor for axis #" << i << " (" << factors[i] << ")is invalid ! Must be > 0 !";
224           throw INTERP_KERNEL::Exception(oss.str().c_str());
225         }
226       int factAbs(std::abs(factors[i]));
227       double fact2(1./(double)factors[i]);
228       structure[i]=(_structure[i]-1)*factAbs+1;
229       dxyz[i]=fact2*_dxyz[i];
230     }
231   std::copy(structure.begin(),structure.end(),_structure);
232   std::copy(dxyz.begin(),dxyz.end(),_dxyz);
233   declareAsNew();
234 }
235
236 /*!
237  * This method returns a newly created mesh containing a single cell in it. This returned cell covers exactly the space covered by \a this.
238  *
239  * \return MEDCouplingIMesh * - A newly created object (to be managed by the caller with decrRef) containing simply one cell.
240  *
241  * \throw if \a this does not pass the \c checkCoherency test.
242  */
243 MEDCouplingIMesh *MEDCouplingIMesh::asSingleCell() const
244 {
245   checkCoherency();
246   int spaceDim(getSpaceDimension()),nodeSt[3];
247   double dxyz[3];
248   for(int i=0;i<spaceDim;i++)
249     {
250       if(_structure[i]>=2)
251         {
252           nodeSt[i]=2;
253           dxyz[i]=(_structure[i]-1)*_dxyz[i];
254         }
255       else
256         {
257           nodeSt[i]=_structure[i];
258           dxyz[i]=_dxyz[i];
259         }
260     }
261   MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(MEDCouplingIMesh::New(getName(),getSpaceDimension(),nodeSt,nodeSt+spaceDim,_origin,_origin+spaceDim,dxyz,dxyz+spaceDim));
262   ret->copyTinyInfoFrom(this);
263   return ret.retn();
264 }
265
266 /*!
267  * This static method is useful to condense field on cells of a MEDCouplingIMesh instance coming from a refinement ( MEDCouplingIMesh::refineWithFactor for example)
268  * to a coarse MEDCouplingIMesh instance. So this method can be seen as a specialization in P0P0 conservative interpolation non overlaping from fine image mesh
269  * to a coarse image mesh. Only tuples ( deduced from \a fineLocInCoarse ) of \a coarseDA will be modified. Other tuples of \a coarseDA will be let unchanged.
270  *
271  * \param [in] coarseSt The cell structure of coarse mesh.
272  * \param [in] fineDA The DataArray containing the cell field on uniformly refined mesh
273  * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
274  * \param [in] facts The refinement coefficient per axis.
275  * \param [in,out] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
276  *
277  * \sa CondenseFineToCoarseGhost,SpreadCoarseToFine
278  */
279 void MEDCouplingIMesh::CondenseFineToCoarse(const std::vector<int>& coarseSt, const DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, DataArrayDouble *coarseDA)
280 {
281   if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
282     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : All input vectors (dimension) must have the same size !");
283   if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
284     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the parameters 1 or 3 are NULL or not allocated !");
285   int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
286   int nbCompo(fineDA->getNumberOfComponents());
287   if(coarseDA->getNumberOfComponents()!=nbCompo)
288     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the number of components of fine DA and coarse one mismatches !");
289   if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
290     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
291   if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
292     {
293       std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
294       throw INTERP_KERNEL::Exception(oss.str().c_str());
295     }
296   int nbTuplesFine(fineDA->getNumberOfTuples());
297   if(nbOfTuplesInFineExp==0)
298     {
299       if(nbTuplesFine==0)
300         return ;
301       else
302         throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : Nothing to condense considering the range specified ! But DataArray is not empty !");
303     }
304   if(nbTuplesFine%nbOfTuplesInFineExp!=0)
305     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : Invalid nb of tuples in fine DataArray regarding its structure !");
306   int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
307   if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
308     {
309       std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Invalid number of tuples ("  << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
310       throw INTERP_KERNEL::Exception(oss.str().c_str());
311     }
312   // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
313   double *outPtr(coarseDA->getPointer());
314   const double *inPtr(fineDA->begin());
315   //
316   std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
317   switch(meshDim)
318   {
319     case 1:
320       {
321         int offset(fineLocInCoarse[0].first),fact0(facts[0]);
322         for(int i=0;i<dims[0];i++)
323           {
324             double *loc(outPtr+(offset+i)*nbCompo);
325             for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
326               {
327                 if(ifact!=0)
328                   std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
329                 else
330                   std::copy(inPtr,inPtr+nbCompo,loc);
331               }
332           }
333         break;
334       }
335     case 2:
336       {
337         int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact1(facts[1]),fact0(facts[0]);
338         for(int j=0;j<dims[1];j++)
339           {
340             for(int jfact=0;jfact<fact1;jfact++)
341               {
342                 for(int i=0;i<dims[0];i++)
343                   {
344                     double *loc(outPtr+(kk+i)*nbCompo);
345                     for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
346                       {
347                         if(jfact!=0 || ifact!=0)
348                           std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
349                         else
350                           std::copy(inPtr,inPtr+nbCompo,loc);
351                       }
352                   }
353               }
354             kk+=coarseSt[0];
355           }
356         break;
357       }
358     case 3:
359       {
360         int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first+coarseSt[0]*coarseSt[1]*fineLocInCoarse[2].first),fact2(facts[2]),fact1(facts[1]),fact0(facts[0]);
361         for(int k=0;k<dims[2];k++)
362           {
363             for(int kfact=0;kfact<fact2;kfact++)
364               {
365                 for(int j=0;j<dims[1];j++)
366                   {
367                     for(int jfact=0;jfact<fact1;jfact++)
368                       {
369                         for(int i=0;i<dims[0];i++)
370                           {
371                             double *loc(outPtr+(kk+i+j*coarseSt[0])*nbCompo);
372                             for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
373                               {
374                                 if(kfact!=0 || jfact!=0 || ifact!=0)
375                                   std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
376                                 else
377                                   std::copy(inPtr,inPtr+nbCompo,loc);
378                               }
379                           }
380                       }
381                   }
382               }
383             kk+=coarseSt[0]*coarseSt[1];
384           }
385         break;
386       }
387     default:
388       throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : only dimensions 1, 2 and 3 supported !");
389   }
390 }
391
392 /*!
393  * This static method is useful to condense field on cells of a MEDCouplingIMesh instance coming from a refinement ( MEDCouplingIMesh::refineWithFactor for example)
394  * to a coarse MEDCouplingIMesh instance. So this method can be seen as a specialization in P0P0 conservative interpolation non overlaping from fine image mesh
395  * to a coarse image mesh. Only tuples ( deduced from \a fineLocInCoarse ) of \a coarseDA will be modified. Other tuples of \a coarseDA will be let unchanged.
396  *
397  * \param [in] coarseSt The cell structure of coarse mesh.
398  * \param [in] fineDA The DataArray containing the cell field on uniformly refined mesh
399  * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
400  * \param [in] facts The refinement coefficient per axis.
401  * \param [in,out] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
402  * \param [in] ghostSize - The size of the ghost zone. The ghost zone is expected to be the same for all axis and both for coarse and fine meshes.
403  *
404  * \sa CondenseFineToCoarse,SpreadCoarseToFineGhost
405  */
406 void MEDCouplingIMesh::CondenseFineToCoarseGhost(const std::vector<int>& coarseSt, const DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, DataArrayDouble *coarseDA, int ghostSize)
407 {
408   if(ghostSize<0)
409     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : ghost level has to be >= 0 !");
410   if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
411     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : All input vectors (dimension) must have the same size !");
412   if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
413     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the parameters 1 or 3 are NULL or not allocated !");
414   std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
415   int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
416   int nbCompo(fineDA->getNumberOfComponents());
417   if(coarseDA->getNumberOfComponents()!=nbCompo)
418     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the number of components of fine DA and coarse one mismatches !");
419   if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
420     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
421   if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
422     {
423       std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarseGhost : Expecting " << nbOfTuplesInCoarseExp << " tuples in coarse DataArray having " << coarseDA->getNumberOfTuples() << " !";
424       throw INTERP_KERNEL::Exception(oss.str().c_str());
425     }
426   //
427   std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
428   std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
429   std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
430   int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
431   if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
432     {
433       std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarseGhost : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
434       throw INTERP_KERNEL::Exception(oss.str().c_str());
435     }
436   //
437   double *outPtr(coarseDA->getPointer());
438   const double *inPtr(fineDA->begin());
439   //
440   std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
441   switch(meshDim)
442   {
443     case 1:
444       {
445         int offset(fineLocInCoarse[0].first+ghostSize),fact0(facts[0]);
446         inPtr+=ghostSize*nbCompo;
447         for(int i=0;i<dims[0];i++)
448           {
449             double *loc(outPtr+(offset+i)*nbCompo);
450             for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
451               {
452                 if(ifact!=0)
453                   std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
454                 else
455                   std::copy(inPtr,inPtr+nbCompo,loc);
456               }
457           }
458         break;
459       }
460     case 2:
461       {
462         int nxwg(coarseSt[0]+2*ghostSize);
463         int kk(fineLocInCoarse[0].first+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize)),fact1(facts[1]),fact0(facts[0]);
464         inPtr+=(dims[0]*fact0+2*ghostSize)*ghostSize*nbCompo;
465         for(int j=0;j<dims[1];j++)
466           {
467              for(int jfact=0;jfact<fact1;jfact++)
468               {
469                 inPtr+=ghostSize*nbCompo;
470                 for(int i=0;i<dims[0];i++)
471                   {
472                     double *loc(outPtr+(kk+i)*nbCompo);
473                     for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
474                       {
475                         if(jfact!=0 || ifact!=0)
476                           std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
477                         else
478                           std::copy(inPtr,inPtr+nbCompo,loc);
479                       }
480                   }
481                 inPtr+=ghostSize*nbCompo;
482               }
483             kk+=nxwg;
484           }
485         break;
486       }
487     case 3:
488       {
489         int nxwg(coarseSt[0]+2*ghostSize),nxywg((coarseSt[0]+2*ghostSize)*(coarseSt[1]+2*ghostSize));
490         int kk(fineLocInCoarse[0].first+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize)+nxywg*(fineLocInCoarse[2].first+ghostSize)),fact2(facts[2]),fact1(facts[1]),fact0(facts[0]);
491         inPtr+=(dims[0]*fact0+2*ghostSize)*(dims[1]*fact1+2*ghostSize)*ghostSize*nbCompo;
492         for(int k=0;k<dims[2];k++)
493           {
494             for(int kfact=0;kfact<fact2;kfact++)
495               {
496                 inPtr+=ghostSize*(dims[0]*fact0+2*ghostSize)*nbCompo;
497                 for(int j=0;j<dims[1];j++)
498                   {
499                     int kky(j*nxwg);
500                     for(int jfact=0;jfact<fact1;jfact++)
501                       {
502                         inPtr+=ghostSize*nbCompo;
503                         for(int i=0;i<dims[0];i++)
504                           {
505                             double *loc(outPtr+(kky+kk+i)*nbCompo);
506                             for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
507                               {
508                                 if(kfact!=0 || jfact!=0 || ifact!=0)
509                                   std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
510                                 else
511                                   std::copy(inPtr,inPtr+nbCompo,loc);
512                               }
513                           }
514                         inPtr+=ghostSize*nbCompo;
515                       }
516                   }
517                 inPtr+=ghostSize*(dims[0]*fact0+2*ghostSize)*nbCompo;
518               }
519             kk+=nxywg;
520           }
521         break;
522       }
523     default:
524       throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : only dimensions 1, 2, 3 supported !");
525   }
526 }
527
528 /*!
529  * This method spreads the values of coarse data \a coarseDA into \a fineDA.
530  *
531  * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
532  * \param [in] coarseSt The cell structure of coarse mesh.
533  * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
534  * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
535  * \param [in] facts The refinement coefficient per axis.
536  * \sa SpreadCoarseToFineGhost, CondenseFineToCoarse
537  */
538 void MEDCouplingIMesh::SpreadCoarseToFine(const DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts)
539 {
540   if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
541       throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : All input vectors (dimension) must have the same size !");
542   if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
543     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the parameters 1 or 3 are NULL or not allocated !");
544   int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
545   int nbCompo(fineDA->getNumberOfComponents());
546   if(coarseDA->getNumberOfComponents()!=nbCompo)
547     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the number of components of fine DA and coarse one mismatches !");
548   if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
549     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
550   if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
551     {
552       std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
553       throw INTERP_KERNEL::Exception(oss.str().c_str());
554     }
555   int nbTuplesFine(fineDA->getNumberOfTuples());
556   if(nbTuplesFine%nbOfTuplesInFineExp!=0)
557     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : Invalid nb of tuples in fine DataArray regarding its structure !");
558   int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
559   if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
560     {
561       std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Invalid number of tuples ("  << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
562       throw INTERP_KERNEL::Exception(oss.str().c_str());
563     }
564   // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
565   double *outPtr(fineDA->getPointer());
566   const double *inPtr(coarseDA->begin());
567   //
568   std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
569   switch(meshDim)
570   {
571     case 1:
572       {
573         int offset(fineLocInCoarse[0].first),fact0(facts[0]);
574         for(int i=0;i<dims[0];i++)
575           {
576             const double *loc(inPtr+(offset+i)*nbCompo);
577             for(int ifact=0;ifact<fact0;ifact++)
578               outPtr=std::copy(loc,loc+nbCompo,outPtr);
579           }
580         break;
581       }
582     case 2:
583       {
584         int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact0(facts[0]),fact1(facts[1]);
585         for(int j=0;j<dims[1];j++)
586           {
587             for(int jfact=0;jfact<fact1;jfact++)
588               {
589                 for(int i=0;i<dims[0];i++)
590                   {
591                     const double *loc(inPtr+(kk+i)*nbCompo);
592                     for(int ifact=0;ifact<fact0;ifact++)
593                       outPtr=std::copy(loc,loc+nbCompo,outPtr);
594                   }
595               }
596             kk+=coarseSt[0];
597           }
598         break;
599       }
600     case 3:
601       {
602         int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first+coarseSt[0]*coarseSt[1]*fineLocInCoarse[2].first),fact0(facts[0]),fact1(facts[2]),fact2(facts[2]);
603         for(int k=0;k<dims[2];k++)
604           {
605             for(int kfact=0;kfact<fact2;kfact++)
606               {
607                 for(int j=0;j<dims[1];j++)
608                   {
609                     for(int jfact=0;jfact<fact1;jfact++)
610                       {
611                         for(int i=0;i<dims[0];i++)
612                           {
613                             const double *loc(inPtr+(kk+i+j*coarseSt[0])*nbCompo);
614                             for(int ifact=0;ifact<fact0;ifact++)
615                               outPtr=std::copy(loc,loc+nbCompo,outPtr);
616                           }
617                       }
618                   }
619               }
620             kk+=coarseSt[0]*coarseSt[1];
621           }
622         break;
623       }
624     default:
625       throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : only dimensions 1, 2 and 3 supported !");
626   }
627 }
628
629 /*!
630  * This method spreads the values of coarse data \a coarseDA into \a fineDA.
631  *
632  * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
633  * \param [in] coarseSt The cell structure of coarse mesh.
634  * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
635  * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
636  * \param [in] facts The refinement coefficient per axis.
637  * \param [in] ghostSize - The size of the ghost zone. The ghost zone is expected to be the same for all axis and both for coarse and fine meshes.
638  *
639  * \sa CondenseFineToCoarse, SpreadCoarseToFineGhostZone
640  */
641 void MEDCouplingIMesh::SpreadCoarseToFineGhost(const DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, int ghostSize)
642 {
643   if(ghostSize<0)
644     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : ghost level has to be >= 0 !");
645   if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
646     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : All input vectors (dimension) must have the same size !");
647   if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
648     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the parameters 1 or 3 are NULL or not allocated !");
649   std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
650   int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
651   int nbCompo(fineDA->getNumberOfComponents());
652   if(coarseDA->getNumberOfComponents()!=nbCompo)
653     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the number of components of fine DA and coarse one mismatches !");
654   if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
655     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
656   if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
657     {
658       std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhost : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
659       throw INTERP_KERNEL::Exception(oss.str().c_str());
660     }
661   //
662   std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
663   std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
664   std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
665   int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
666   if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
667     {
668       std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhost : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
669       throw INTERP_KERNEL::Exception(oss.str().c_str());
670     }
671   //
672   double *outPtr(fineDA->getPointer());
673   const double *inPtr(coarseDA->begin());
674   //
675   switch(meshDim)
676   {
677     case 1:
678       {
679         std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
680         int offset(fineLocInCoarse[0].first+ghostSize-1),fact0(facts[0]);//offset is always >=0 thanks to the fact that ghostSize>=1 !
681         for(int i=0;i<ghostSize;i++)
682           outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
683         offset=fineLocInCoarse[0].first+ghostSize;
684         for(int i=0;i<dims[0];i++)
685           {
686             const double *loc(inPtr+(offset+i)*nbCompo);
687             for(int ifact=0;ifact<fact0;ifact++)
688               outPtr=std::copy(loc,loc+nbCompo,outPtr);
689           }
690         offset=fineLocInCoarse[0].second+ghostSize;
691         for(int i=0;i<ghostSize;i++)
692           outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
693         break;
694       }
695     case 2:
696       {
697         SpreadCoarseToFineGhost2D(inPtr,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
698         break;
699       }
700     case 3:
701       {
702         std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
703         int fact0(facts[0]),fact1(facts[1]),fact2(facts[2]);
704         int nxyWgCoarse((coarseSt[0]+2*ghostSize)*(coarseSt[1]+2*ghostSize)),nxyWgFine((dims[0]*fact0+2*ghostSize)*(dims[1]*fact1+2*ghostSize));
705         int offset((fineLocInCoarse[2].first+ghostSize-1)*nxyWgCoarse);//offset is always >=0 thanks to the fact that ghostSize>=1 !
706         for(int i=0;i<ghostSize;i++,outPtr+=nxyWgFine*nbCompo)
707           SpreadCoarseToFineGhost2D(inPtr+offset*nbCompo,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
708         offset+=nxyWgCoarse;
709         for(int i=0;i<dims[2];i++,offset+=nxyWgCoarse)
710           for(int j=0;j<fact2;j++,outPtr+=nxyWgFine*nbCompo)
711             SpreadCoarseToFineGhost2D(inPtr+offset*nbCompo,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
712         for(int i=0;i<ghostSize;i++,outPtr+=nxyWgFine*nbCompo)
713           SpreadCoarseToFineGhost2D(inPtr+offset*nbCompo,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
714         break;
715       }
716     default:
717       throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : only dimensions 1, 2, 3 supported !");
718   }
719 }
720
721 /*!
722  * This method spreads the values of coarse data \a coarseDA into \a fineDA \b ONLY \b in \b the \b ghost \b zone (contrary to SpreadCoarseToFineGhost that spread the values everywhere).
723  *
724  * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
725  * \param [in] coarseSt The cell structure of coarse mesh.
726  * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
727  * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
728  * \param [in] facts The refinement coefficient per axis.
729  * \param [in] ghostSize - The size of the ghost zone. The ghost zone is expected to be the same for all axis and both for coarse and fine meshes.
730  *
731  * \sa SpreadCoarseToFineGhost
732  */
733 void MEDCouplingIMesh::SpreadCoarseToFineGhostZone(const DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, int ghostSize)
734 {
735   if(ghostSize<0)
736     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : ghost level has to be >= 0 !");
737   if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
738     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : All input vectors (dimension) must have the same size !");
739   if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
740     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the parameters 1 or 3 are NULL or not allocated !");
741   std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
742   int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
743   int nbCompo(fineDA->getNumberOfComponents());
744   if(coarseDA->getNumberOfComponents()!=nbCompo)
745     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the number of components of fine DA and coarse one mismatches !");
746   if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
747     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
748   if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
749     {
750       std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhostZone : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
751       throw INTERP_KERNEL::Exception(oss.str().c_str());
752     }
753   //
754   std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
755   std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
756   std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
757   int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
758   if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
759     {
760       std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhostZone : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
761       throw INTERP_KERNEL::Exception(oss.str().c_str());
762     }
763   //
764   double *outPtr(fineDA->getPointer());
765   const double *inPtr(coarseDA->begin());
766   //
767   std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
768   switch(meshDim)
769   {
770     case 1:
771       {
772         int offset(fineLocInCoarse[0].first+ghostSize-1),fact0(facts[0]);//offset is always >=0 thanks to the fact that ghostSize>=1 !
773         for(int i=0;i<ghostSize;i++)
774           outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
775         outPtr+=nbCompo*fact0*dims[0];
776         offset=fineLocInCoarse[0].second+ghostSize;
777         for(int i=0;i<ghostSize;i++)
778           outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
779         break;
780       }
781     case 2:
782       {
783         SpreadCoarseToFineGhostZone2D(inPtr,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
784         break;
785       }
786     case 3:
787       {
788         int fact0(facts[0]),fact1(facts[1]),fact2(facts[2]);
789         int nxyWgCoarse((coarseSt[0]+2*ghostSize)*(coarseSt[1]+2*ghostSize)),nxyWgFine((dims[0]*fact0+2*ghostSize)*(dims[1]*fact1+2*ghostSize));
790         int offset((fineLocInCoarse[2].first+ghostSize-1)*nxyWgCoarse);//offset is always >=0 thanks to the fact that ghostSize>=1 !
791         for(int i=0;i<ghostSize;i++,outPtr+=nxyWgFine*nbCompo)
792           SpreadCoarseToFineGhost2D(inPtr+offset*nbCompo,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
793         offset+=nxyWgCoarse;
794         for(int i=0;i<dims[2];i++,offset+=nxyWgCoarse)
795           for(int j=0;j<fact2;j++,outPtr+=nxyWgFine*nbCompo)
796             SpreadCoarseToFineGhostZone2D(inPtr+offset*nbCompo,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
797         for(int i=0;i<ghostSize;i++,outPtr+=nxyWgFine*nbCompo)
798           SpreadCoarseToFineGhost2D(inPtr+offset*nbCompo,outPtr,nbCompo,coarseSt,fineLocInCoarse,facts,ghostSize);
799         break;
800       }
801     default:
802       throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : only dimensions 1, 2, 3 supported !");
803   }
804 }
805
806 void MEDCouplingIMesh::setSpaceDimension(int spaceDim)
807 {
808   if(spaceDim==_space_dim)
809     return ;
810   CheckSpaceDimension(spaceDim);
811   _space_dim=spaceDim;
812   declareAsNew();
813 }
814
815 void MEDCouplingIMesh::updateTime() const
816 {
817 }
818
819 std::size_t MEDCouplingIMesh::getHeapMemorySizeWithoutChildren() const
820 {
821   return MEDCouplingStructuredMesh::getHeapMemorySizeWithoutChildren();
822 }
823
824 std::vector<const BigMemoryObject *> MEDCouplingIMesh::getDirectChildrenWithNull() const
825 {
826   return std::vector<const BigMemoryObject *>();
827 }
828
829 /*!
830  * This method copyies all tiny strings from other (name and components name).
831  * @throw if other and this have not same mesh type.
832  */
833 void MEDCouplingIMesh::copyTinyStringsFrom(const MEDCouplingMesh *other)
834
835   const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
836   if(!otherC)
837     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::copyTinyStringsFrom : meshes have not same type !");
838   MEDCouplingStructuredMesh::copyTinyStringsFrom(other);
839   declareAsNew();
840 }
841
842 bool MEDCouplingIMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
843 {
844   if(!other)
845     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::isEqualIfNotWhy : input other pointer is null !");
846   const MEDCouplingIMesh *otherC(dynamic_cast<const MEDCouplingIMesh *>(other));
847   if(!otherC)
848     {
849       reason="mesh given in input is not castable in MEDCouplingIMesh !";
850       return false;
851     }
852   if(!MEDCouplingStructuredMesh::isEqualIfNotWhy(other,prec,reason))
853     return false;
854   if(!isEqualWithoutConsideringStrInternal(otherC,prec,reason))
855     return false;
856   if(_axis_unit!=otherC->_axis_unit)
857     {
858       reason="The units of axis are not the same !";
859       return false;
860     }
861   return true;
862 }
863
864 bool MEDCouplingIMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
865 {
866   const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
867   if(!otherC)
868     return false;
869   std::string tmp;
870   return isEqualWithoutConsideringStrInternal(other,prec,tmp);
871 }
872
873 bool MEDCouplingIMesh::isEqualWithoutConsideringStrInternal(const MEDCouplingMesh *other, double prec, std::string& reason) const
874 {
875   const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
876   if(!otherC)
877     return false;
878   if(_space_dim!=otherC->_space_dim)
879     {
880       std::ostringstream oss;
881       oss << "The spaceDimension of this (" << _space_dim << ") is not equal to those of other (" << otherC->_space_dim << ") !";
882       return false;
883     }
884   checkSpaceDimension();
885   for(int i=0;i<_space_dim;i++)
886     {
887       if(fabs(_origin[i]-otherC->_origin[i])>prec)
888         {
889           std::ostringstream oss;
890           oss << "The origin of this and other differs at " << i << " !";
891           reason=oss.str();
892           return false;
893         }
894     }
895   for(int i=0;i<_space_dim;i++)
896     {
897       if(fabs(_dxyz[i]-otherC->_dxyz[i])>prec)
898         {
899           std::ostringstream oss;
900           oss << "The delta of this and other differs at " << i << " !";
901           reason=oss.str();
902           return false;
903         }
904     }
905   for(int i=0;i<_space_dim;i++)
906     {
907       if(_structure[i]!=otherC->_structure[i])
908         {
909           std::ostringstream oss;
910           oss << "The structure of this and other differs at " << i << " !";
911           reason=oss.str();
912           return false;
913         }
914     }
915   return true;
916 }
917
918 void MEDCouplingIMesh::checkDeepEquivalWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
919                                             DataArrayInt *&cellCor, DataArrayInt *&nodeCor) const
920 {
921   if(!isEqualWithoutConsideringStr(other,prec))
922     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalWith : Meshes are not the same !");
923 }
924
925 /*!
926  * Nothing is done here (except to check that the other is a ParaMEDMEM::MEDCouplingIMesh instance too).
927  * The user intend that the nodes are the same, so by construction of ParaMEDMEM::MEDCouplingIMesh, \a this and \a other are the same !
928  */
929 void MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
930                                                        DataArrayInt *&cellCor) const
931 {
932   if(!isEqualWithoutConsideringStr(other,prec))
933     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith : Meshes are not the same !");
934 }
935
936 void MEDCouplingIMesh::checkCoherency() const
937 {
938   checkSpaceDimension();
939   for(int i=0;i<_space_dim;i++)
940     if(_structure[i]<1)
941       {
942         std::ostringstream oss; oss << "MEDCouplingIMesh::checkCoherency : On axis " << i << "/" << _space_dim << ", number of nodes is equal to " << _structure[i] << " ! must be >=1 !";
943         throw INTERP_KERNEL::Exception(oss.str().c_str());
944       }
945 }
946
947 void MEDCouplingIMesh::checkCoherency1(double eps) const
948 {
949   checkCoherency();
950 }
951
952 void MEDCouplingIMesh::checkCoherency2(double eps) const
953 {
954   checkCoherency1(eps);
955 }
956
957 void MEDCouplingIMesh::getNodeGridStructure(int *res) const
958 {
959   checkSpaceDimension();
960   std::copy(_structure,_structure+_space_dim,res);
961 }
962
963 std::vector<int> MEDCouplingIMesh::getNodeGridStructure() const
964 {
965   checkSpaceDimension();
966   std::vector<int> ret(_structure,_structure+_space_dim);
967   return ret;
968 }
969
970 MEDCouplingStructuredMesh *MEDCouplingIMesh::buildStructuredSubPart(const std::vector< std::pair<int,int> >& cellPart) const
971 {
972   checkCoherency();
973   int dim(getSpaceDimension());
974   if(dim!=(int)cellPart.size())
975     {
976       std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : the space dimension is " << dim << " and cell part size is " << cellPart.size() << " !";
977       throw INTERP_KERNEL::Exception(oss.str().c_str());
978     }
979   double retOrigin[3]={0.,0.,0.};
980   int retStruct[3]={0,0,0};
981   MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(dynamic_cast<MEDCouplingIMesh *>(deepCpy()));
982   for(int i=0;i<dim;i++)
983     {
984       int startNode(cellPart[i].first),endNode(cellPart[i].second+1);
985       int myDelta(endNode-startNode);
986       if(startNode<0 || startNode>=_structure[i])
987         {
988           std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : At dimension #" << i << " the start node id is " << startNode << " it should be in [0," << _structure[i] << ") !";
989           throw INTERP_KERNEL::Exception(oss.str().c_str());
990         }
991       if(myDelta<0 || myDelta>_structure[i])
992         {
993           std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : Along dimension #" << i << " the number of nodes is " << _structure[i] << ", and you are requesting for " << myDelta << " nodes wide range !" << std::endl;
994           throw INTERP_KERNEL::Exception(oss.str().c_str());
995         }
996       retOrigin[i]=_origin[i]+startNode*_dxyz[i];
997       retStruct[i]=myDelta;
998     }
999   ret->setNodeStruct(retStruct,retStruct+dim);
1000   ret->setOrigin(retOrigin,retOrigin+dim);
1001   ret->checkCoherency();
1002   return ret.retn();
1003 }
1004
1005 /*!
1006  * Return the space dimension of \a this.
1007  */
1008 int MEDCouplingIMesh::getSpaceDimension() const
1009 {
1010   return _space_dim;
1011 }
1012
1013 void MEDCouplingIMesh::getCoordinatesOfNode(int nodeId, std::vector<double>& coo) const
1014 {
1015   int tmp[3];
1016   int spaceDim(getSpaceDimension());
1017   getSplitNodeValues(tmp);
1018   int tmp2[3];
1019   GetPosFromId(nodeId,spaceDim,tmp,tmp2);
1020   for(int j=0;j<spaceDim;j++)
1021     coo.push_back(_origin[j]+_dxyz[j]*tmp2[j]);
1022 }
1023
1024 std::string MEDCouplingIMesh::simpleRepr() const
1025 {
1026   std::ostringstream ret;
1027   ret << "Image grid with name : \"" << getName() << "\"\n";
1028   ret << "Description of mesh : \"" << getDescription() << "\"\n";
1029   int tmpp1,tmpp2;
1030   double tt(getTime(tmpp1,tmpp2));
1031   int spaceDim(_space_dim);
1032   ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
1033   ret << "Iteration : " << tmpp1  << " Order : " << tmpp2 << "\n";
1034   ret << "Space dimension : " << spaceDim << "\n";
1035   if(spaceDim<0 || spaceDim>3)
1036     return ret.str();
1037   ret << "The nodal structure is : "; std::copy(_structure,_structure+spaceDim,std::ostream_iterator<int>(ret," ")); ret << "\n";
1038   ret << "The origin position is [" << _axis_unit << "]: ";
1039   std::copy(_origin,_origin+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
1040   ret << "The intervals along axis are : ";
1041   std::copy(_dxyz,_dxyz+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
1042   return ret.str();
1043 }
1044
1045 std::string MEDCouplingIMesh::advancedRepr() const
1046 {
1047   return simpleRepr();
1048 }
1049
1050 void MEDCouplingIMesh::getBoundingBox(double *bbox) const
1051 {
1052   checkCoherency();
1053   int dim(getSpaceDimension());
1054   for(int idim=0; idim<dim; idim++)
1055     {
1056       bbox[2*idim]=_origin[idim];
1057       int coeff(_structure[idim]);
1058       if(_structure[idim]<0)
1059         {
1060           std::ostringstream oss; oss << "MEDCouplingIMesh::getBoundingBox : on axis #" << idim << " number of nodes in structure is < 0 !";
1061           throw INTERP_KERNEL::Exception(oss.str().c_str());
1062         }
1063       if(_structure[idim]>1)
1064         coeff=_structure[idim]-1;
1065       bbox[2*idim+1]=_origin[idim]+_dxyz[idim]*coeff;
1066     }
1067 }
1068
1069 /*!
1070  * Returns a new MEDCouplingFieldDouble containing volumes of cells constituting \a this
1071  * mesh.<br>
1072  * For 1D cells, the returned field contains lengths.<br>
1073  * For 2D cells, the returned field contains areas.<br>
1074  * For 3D cells, the returned field contains volumes.
1075  *  \param [in] isAbs - a not used parameter.
1076  *  \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble on cells
1077  *         and one time . The caller is to delete this field using decrRef() as it is no
1078  *         more needed.
1079  */
1080 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureField(bool isAbs) const
1081 {
1082   checkCoherency();
1083   std::string name="MeasureOfMesh_";
1084   name+=getName();
1085   int nbelem(getNumberOfCells());
1086   MEDCouplingFieldDouble *field(MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME));
1087   field->setName(name);
1088   DataArrayDouble* array(DataArrayDouble::New());
1089   array->alloc(nbelem,1);
1090   array->fillWithValue(getMeasureOfAnyCell());
1091   field->setArray(array) ;
1092   array->decrRef();
1093   field->setMesh(const_cast<MEDCouplingIMesh *>(this));
1094   field->synchronizeTimeWithMesh();
1095   return field;
1096 }
1097
1098 /*!
1099  * not implemented yet !
1100  */
1101 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureFieldOnNode(bool isAbs) const
1102 {
1103   throw INTERP_KERNEL::Exception("MEDCouplingIMesh::getMeasureFieldOnNode : not implemented yet !");
1104   //return 0;
1105 }
1106
1107 int MEDCouplingIMesh::getCellContainingPoint(const double *pos, double eps) const
1108 {
1109   int dim(getSpaceDimension()),ret(0),coeff(1);
1110   for(int i=0;i<dim;i++)
1111     {
1112       int nbOfCells(_structure[i]-1);
1113       double ref(pos[i]);
1114       int tmp((int)((ref-_origin[i])/_dxyz[i]));
1115       if(tmp>=0 && tmp<nbOfCells)
1116         {
1117           ret+=coeff*tmp;
1118           coeff*=nbOfCells;
1119         }
1120       else
1121         return -1;
1122     }
1123   return ret;
1124 }
1125
1126 void MEDCouplingIMesh::rotate(const double *center, const double *vector, double angle)
1127 {
1128   throw INTERP_KERNEL::Exception("No rotation available on IMesh : Traduce it to unstructured mesh to apply it !");
1129 }
1130
1131 /*!
1132  * Translates all nodes of \a this mesh by a given vector. Actually, it adds each
1133  * component of the \a vector to all node coordinates of a corresponding axis.
1134  *  \param [in] vector - the translation vector whose size must be not less than \a
1135  *         this->getSpaceDimension().
1136  */
1137 void MEDCouplingIMesh::translate(const double *vector)
1138 {
1139   checkSpaceDimension();
1140   int dim(getSpaceDimension());
1141   std::transform(_origin,_origin+dim,vector,_origin,std::plus<double>());
1142   declareAsNew();
1143 }
1144
1145 /*!
1146  * Applies scaling transformation to all nodes of \a this mesh.
1147  *  \param [in] point - coordinates of a scaling center. This array is to be of
1148  *         size \a this->getSpaceDimension() at least.
1149  *  \param [in] factor - a scale factor.
1150  */
1151 void MEDCouplingIMesh::scale(const double *point, double factor)
1152 {
1153   checkSpaceDimension();
1154   int dim(getSpaceDimension());
1155   std::transform(_origin,_origin+dim,point,_origin,std::minus<double>());
1156   std::transform(_origin,_origin+dim,_origin,std::bind2nd(std::multiplies<double>(),factor));
1157   std::transform(_dxyz,_dxyz+dim,_dxyz,std::bind2nd(std::multiplies<double>(),factor));
1158   std::transform(_origin,_origin+dim,point,_origin,std::plus<double>());
1159   declareAsNew();
1160 }
1161
1162 MEDCouplingMesh *MEDCouplingIMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
1163 {
1164   //not implemented yet !
1165   return 0;
1166 }
1167
1168 /*!
1169  * Returns a new DataArrayDouble holding coordinates of all nodes of \a this mesh.
1170  *  \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
1171  *          this->getNumberOfNodes() tuples per \a this->getSpaceDimension()
1172  *          components. The caller is to delete this array using decrRef() as it is
1173  *          no more needed.
1174  */
1175 DataArrayDouble *MEDCouplingIMesh::getCoordinatesAndOwner() const
1176 {
1177   checkCoherency();
1178   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
1179   int spaceDim(getSpaceDimension()),nbNodes(getNumberOfNodes());
1180   ret->alloc(nbNodes,spaceDim);
1181   double *pt(ret->getPointer());
1182   ret->setInfoOnComponents(buildInfoOnComponents());
1183   int tmp2[3],tmp[3];
1184   getSplitNodeValues(tmp);
1185   for(int i=0;i<nbNodes;i++)
1186     {
1187       GetPosFromId(i,spaceDim,tmp,tmp2);
1188       for(int j=0;j<spaceDim;j++)
1189         pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+_origin[j];
1190     }
1191   return ret.retn();
1192 }
1193
1194 /*!
1195  * Returns a new DataArrayDouble holding barycenters of all cells. The barycenter is
1196  * computed by averaging coordinates of cell nodes.
1197  *  \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
1198  *          this->getNumberOfCells() tuples per \a this->getSpaceDimension()
1199  *          components. The caller is to delete this array using decrRef() as it is
1200  *          no more needed.
1201  */
1202 DataArrayDouble *MEDCouplingIMesh::getBarycenterAndOwner() const
1203 {
1204   checkCoherency();
1205   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
1206   int spaceDim(getSpaceDimension()),nbCells(getNumberOfCells()),tmp[3],tmp2[3];
1207   ret->alloc(nbCells,spaceDim);
1208   double *pt(ret->getPointer()),shiftOrigin[3];
1209   std::transform(_dxyz,_dxyz+spaceDim,shiftOrigin,std::bind2nd(std::multiplies<double>(),0.5));
1210   std::transform(_origin,_origin+spaceDim,shiftOrigin,shiftOrigin,std::plus<double>());
1211   getSplitCellValues(tmp);
1212   ret->setInfoOnComponents(buildInfoOnComponents());
1213   for(int i=0;i<nbCells;i++)
1214     {
1215       GetPosFromId(i,spaceDim,tmp,tmp2);
1216       for(int j=0;j<spaceDim;j++)
1217         pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+shiftOrigin[j];
1218     }
1219   return ret.retn();
1220 }
1221
1222 DataArrayDouble *MEDCouplingIMesh::computeIsoBarycenterOfNodesPerCell() const
1223 {
1224   return MEDCouplingIMesh::getBarycenterAndOwner();
1225 }
1226
1227 void MEDCouplingIMesh::renumberCells(const int *old2NewBg, bool check)
1228 {
1229   throw INTERP_KERNEL::Exception("Functionnality of renumbering cell not available for IMesh !");
1230 }
1231
1232 void MEDCouplingIMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
1233 {
1234   int it,order;
1235   double time(getTime(it,order));
1236   tinyInfo.clear();
1237   tinyInfoD.clear();
1238   littleStrings.clear();
1239   littleStrings.push_back(getName());
1240   littleStrings.push_back(getDescription());
1241   littleStrings.push_back(getTimeUnit());
1242   littleStrings.push_back(getAxisUnit());
1243   tinyInfo.push_back(it);
1244   tinyInfo.push_back(order);
1245   tinyInfo.push_back(_space_dim);
1246   tinyInfo.insert(tinyInfo.end(),_structure,_structure+3);
1247   tinyInfoD.push_back(time);
1248   tinyInfoD.insert(tinyInfoD.end(),_dxyz,_dxyz+3);
1249   tinyInfoD.insert(tinyInfoD.end(),_origin,_origin+3);
1250 }
1251
1252 void MEDCouplingIMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
1253 {
1254   a1->alloc(0,1);
1255   a2->alloc(0,1);
1256 }
1257
1258 void MEDCouplingIMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
1259 {
1260   a1=DataArrayInt::New();
1261   a1->alloc(0,1);
1262   a2=DataArrayDouble::New();
1263   a2->alloc(0,1);
1264 }
1265
1266 void MEDCouplingIMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
1267                                        const std::vector<std::string>& littleStrings)
1268 {
1269   setName(littleStrings[0]);
1270   setDescription(littleStrings[1]);
1271   setTimeUnit(littleStrings[2]);
1272   setAxisUnit(littleStrings[3]);
1273   setTime(tinyInfoD[0],tinyInfo[0],tinyInfo[1]);
1274   _space_dim=tinyInfo[2];
1275   _structure[0]=tinyInfo[3]; _structure[1]=tinyInfo[4]; _structure[2]=tinyInfo[5];
1276   _dxyz[0]=tinyInfoD[1]; _dxyz[1]=tinyInfoD[2]; _dxyz[2]=tinyInfoD[3];
1277   _origin[0]=tinyInfoD[4]; _origin[1]=tinyInfoD[5]; _origin[2]=tinyInfoD[6];
1278   declareAsNew();
1279 }
1280
1281 void MEDCouplingIMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const
1282 {
1283   checkCoherency();
1284   std::ostringstream extent,origin,spacing;
1285   for(int i=0;i<3;i++)
1286     {
1287       if(i<_space_dim)
1288         { extent << "0 " <<  _structure[i]-1 << " "; origin << _origin[i] << " "; spacing << _dxyz[i] << " "; }
1289       else
1290         { extent << "0 0 "; origin << "0 "; spacing << "0 "; }
1291     }
1292   ofs << "  <" << getVTKDataSetType() << " WholeExtent=\"" << extent.str() << "\" Origin=\"" << origin.str() << "\" Spacing=\"" << spacing.str() << "\">\n";
1293   ofs << "    <Piece Extent=\"" << extent.str() << "\">\n";
1294   ofs << "      <PointData>\n" << pointData << std::endl;
1295   ofs << "      </PointData>\n";
1296   ofs << "      <CellData>\n" << cellData << std::endl;
1297   ofs << "      </CellData>\n";
1298   ofs << "      <Coordinates>\n";
1299   ofs << "      </Coordinates>\n";
1300   ofs << "    </Piece>\n";
1301   ofs << "  </" << getVTKDataSetType() << ">\n";
1302 }
1303
1304 void MEDCouplingIMesh::reprQuickOverview(std::ostream& stream) const
1305 {
1306   stream << "MEDCouplingIMesh C++ instance at " << this << ". Name : \"" << getName() << "\". Space dimension : " << _space_dim << ".";
1307   if(_space_dim<0 || _space_dim>3)
1308     return ;
1309   stream << "\n";
1310   std::ostringstream stream0,stream1;
1311   int nbNodes(1),nbCells(0);
1312   bool isPb(false);
1313   for(int i=0;i<_space_dim;i++)
1314     {
1315       char tmp('X'+i);
1316       int tmpNodes(_structure[i]);
1317       stream1 << "- Axis " << tmp << " : " << tmpNodes << " nodes (orig=" << _origin[i] << ", inter=" << _dxyz[i] << ").";
1318       if(i!=_space_dim-1)
1319         stream1 << std::endl;
1320       if(tmpNodes>=1)
1321         nbNodes*=tmpNodes;
1322       else
1323         isPb=true;
1324       if(tmpNodes>=2)
1325         nbCells=nbCells==0?tmpNodes-1:nbCells*(tmpNodes-1);
1326     }
1327   if(!isPb)
1328     {
1329       stream0 << "Number of cells : " << nbCells << ", Number of nodes : " << nbNodes;
1330       stream << stream0.str();
1331       if(_space_dim>0)
1332         stream << std::endl;
1333     }
1334   stream << stream1.str();
1335 }
1336
1337 std::string MEDCouplingIMesh::getVTKFileExtension() const
1338 {
1339   return std::string("vti");
1340 }
1341
1342 std::string MEDCouplingIMesh::getVTKDataSetType() const
1343 {
1344   return std::string("ImageData");
1345 }
1346
1347 std::vector<std::string> MEDCouplingIMesh::buildInfoOnComponents() const
1348 {
1349   checkSpaceDimension();
1350   int dim(getSpaceDimension());
1351   std::vector<std::string> ret(dim);
1352   for(int i=0;i<dim;i++)
1353     {
1354       std::ostringstream oss;
1355       char tmp('X'+i); oss << tmp;
1356       ret[i]=DataArray::BuildInfoFromVarAndUnit(oss.str(),_axis_unit);
1357     }
1358   return ret;
1359 }
1360
1361 void MEDCouplingIMesh::checkSpaceDimension() const
1362 {
1363   CheckSpaceDimension(_space_dim);
1364 }
1365
1366 void MEDCouplingIMesh::CheckSpaceDimension(int spaceDim)
1367 {
1368   if(spaceDim<0 || spaceDim>3)
1369     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CheckSpaceDimension : input spaceDim must be in [0,1,2,3] !");
1370 }
1371
1372 int MEDCouplingIMesh::FindIntRoot(int val, int order)
1373 {
1374   if(order==0)
1375     return 1;
1376   if(val<0)
1377     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : input val is < 0 ! Not possible to compute a root !");
1378   if(order==1)
1379     return val;
1380   if(order!=2 && order!=3)
1381     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the order available are 0,1,2 or 3 !");
1382   double valf((double)val);
1383   if(order==2)
1384     {
1385       double retf(sqrt(valf));
1386       int ret((int)retf);
1387       if(ret*ret!=val)
1388         throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect square root !");
1389       return ret;
1390     }
1391   else//order==3
1392     {
1393       double retf(std::pow(val,0.3333333333333333));
1394       int ret((int)retf),ret2(ret+1);
1395       if(ret*ret*ret!=val && ret2*ret2*ret2!=val)
1396         throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect cublic root !");
1397       if(ret*ret*ret==val)
1398         return ret;
1399       else
1400         return ret2;
1401     }
1402 }
1403
1404 void MEDCouplingIMesh::SpreadCoarseToFineGhost2D(const double *inPtr, double *outPtr, int nbCompo, const std::vector<int>& coarseSt, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, int ghostSize)
1405 {
1406   double *outPtrWork(outPtr);
1407   std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
1408   int nxwg(coarseSt[0]+2*ghostSize),fact0(facts[0]),fact1(facts[1]);
1409   int kk(fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].first+ghostSize-1));//kk is always >=0 thanks to the fact that ghostSize>=1 !
1410   for(int jg=0;jg<ghostSize;jg++)
1411     {
1412       for(int ig=0;ig<ghostSize;ig++)
1413         outPtrWork=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtrWork);
1414       int kk0(kk+1);
1415       for(int ig=0;ig<dims[0];ig++,kk0++)
1416         for(int ifact=0;ifact<fact0;ifact++)
1417           outPtrWork=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtrWork);
1418       for(int ik=0;ik<ghostSize;ik++)
1419         outPtrWork=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtrWork);
1420     }
1421   for(int j=0;j<dims[1];j++)
1422     {
1423       kk=fineLocInCoarse[0].first-1+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize+j);
1424       for(int jfact=0;jfact<fact1;jfact++)
1425         {
1426           for(int ig=0;ig<ghostSize;ig++)
1427             outPtrWork=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtrWork);
1428           int kk0(kk+1);//1 not ghost. We make the hypothesis that factors is >= ghostlev
1429           for(int i=0;i<dims[0];i++,kk0++)
1430             {
1431               const double *loc(inPtr+kk0*nbCompo);
1432               for(int ifact=0;ifact<fact0;ifact++)
1433                 outPtrWork=std::copy(loc,loc+nbCompo,outPtrWork);
1434             }
1435           for(int ig=0;ig<ghostSize;ig++)
1436             outPtrWork=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtrWork);
1437         }
1438     }
1439   kk=fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].second+ghostSize);
1440   for(int jg=0;jg<ghostSize;jg++)
1441     {
1442       for(int ig=0;ig<ghostSize;ig++)
1443         outPtrWork=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtrWork);
1444       int kk0(kk+1);
1445       for(int ig=0;ig<dims[0];ig++,kk0++)
1446         for(int ifact=0;ifact<fact0;ifact++)
1447           outPtrWork=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtrWork);
1448       for(int ik=0;ik<ghostSize;ik++)
1449         outPtrWork=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtrWork);
1450     }
1451 }
1452
1453 void MEDCouplingIMesh::SpreadCoarseToFineGhostZone2D(const double *inPtr, double *outPtr, int nbCompo, const std::vector<int>& coarseSt, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, int ghostSize)
1454 {
1455   double *outPtr2(outPtr);
1456   std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
1457   int nxwg(coarseSt[0]+2*ghostSize),fact0(facts[0]),fact1(facts[1]);
1458   int kk(fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].first+ghostSize-1));//kk is always >=0 thanks to the fact that ghostSize>=1 !
1459   for(int jg=0;jg<ghostSize;jg++)
1460     {
1461       for(int ig=0;ig<ghostSize;ig++)
1462         outPtr2=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr2);
1463       int kk0(kk+1);
1464       for(int ig=0;ig<dims[0];ig++,kk0++)
1465         for(int ifact=0;ifact<fact0;ifact++)
1466           outPtr2=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr2);
1467       for(int ik=0;ik<ghostSize;ik++)
1468         outPtr2=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr2);
1469     }
1470   for(int j=0;j<dims[1];j++)
1471     {
1472       kk=fineLocInCoarse[0].first-1+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize+j);
1473       for(int jfact=0;jfact<fact1;jfact++)
1474         {
1475           for(int ig=0;ig<ghostSize;ig++)
1476             outPtr2=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr2);
1477           int kk0(kk+1+dims[0]);//1 not ghost. We make the hypothesis that factors is >= ghostlev
1478           outPtr2+=fact0*nbCompo*dims[0];
1479           for(int ig=0;ig<ghostSize;ig++)
1480             outPtr2=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr2);
1481         }
1482     }
1483   kk=fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].second+ghostSize);
1484   for(int jg=0;jg<ghostSize;jg++)
1485     {
1486       for(int ig=0;ig<ghostSize;ig++)
1487         outPtr2=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr2);
1488       int kk0(kk+1);
1489       for(int ig=0;ig<dims[0];ig++,kk0++)
1490         for(int ifact=0;ifact<fact0;ifact++)
1491           outPtr2=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr2);
1492       for(int ik=0;ik<ghostSize;ik++)
1493         outPtr2=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr2);
1494     }
1495 }