src/MEDCoupling/MEDCouplingIMesh.cxx

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