src/MEDCoupling/MEDCouplingIMesh.cxx

   1 // Copyright (C) 2007-2014  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 #include <cmath>
  31
  32 using namespace ParaMEDMEM;
  33
  34 MEDCouplingIMesh::MEDCouplingIMesh():_space_dim(-1)
  35 {
  36   _origin[0]=0.; _origin[1]=0.; _origin[2]=0.;
  37   _dxyz[0]=0.; _dxyz[1]=0.; _dxyz[2]=0.;
  38   _structure[0]=0; _structure[1]=0; _structure[2]=0;
  39 }
  40
  41 MEDCouplingIMesh::MEDCouplingIMesh(const MEDCouplingIMesh& other, bool deepCopy):MEDCouplingStructuredMesh(other,deepCopy),_space_dim(other._space_dim),_axis_unit(other._axis_unit)
  42 {
  43   _origin[0]=other._origin[0]; _origin[1]=other._origin[1]; _origin[2]=other._origin[2];
  44   _dxyz[0]=other._dxyz[0]; _dxyz[1]=other._dxyz[1]; _dxyz[2]=other._dxyz[2];
  45   _structure[0]=other._structure[0]; _structure[1]=other._structure[1]; _structure[2]=other._structure[2];
  46 }
  47
  48 MEDCouplingIMesh::~MEDCouplingIMesh()
  49 {
  50 }
  51
  52 MEDCouplingIMesh *MEDCouplingIMesh::New()
  53 {
  54   return new MEDCouplingIMesh;
  55 }
  56
  57 MEDCouplingIMesh *MEDCouplingIMesh::New(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
  58                                         const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
  59 {
  60   MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(new MEDCouplingIMesh);
  61   ret->setName(meshName);
  62   ret->setSpaceDimension(spaceDim);
  63   ret->setNodeStruct(nodeStrctStart,nodeStrctStop);
  64   ret->setOrigin(originStart,originStop);
  65   ret->setDXYZ(dxyzStart,dxyzStop);
  66   return ret.retn();
  67 }
  68
  69 MEDCouplingMesh *MEDCouplingIMesh::deepCpy() const
  70 {
  71   return clone(true);
  72 }
  73
  74 MEDCouplingIMesh *MEDCouplingIMesh::clone(bool recDeepCpy) const
  75 {
  76   return new MEDCouplingIMesh(*this,recDeepCpy);
  77 }
  78
  79 void MEDCouplingIMesh::setNodeStruct(const int *nodeStrctStart, const int *nodeStrctStop)
  80 {
  81   checkSpaceDimension();
  82   int sz((int)std::distance(nodeStrctStart,nodeStrctStop));
  83   if(sz!=_space_dim)
  84     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setNodeStruct : input vector of node structure has not the right size ! Or change space dimension before calling it !");
  85   std::copy(nodeStrctStart,nodeStrctStop,_structure);
  86   declareAsNew();
  87 }
  88
  89 std::vector<int> MEDCouplingIMesh::getNodeStruct() const
  90 {
  91   checkSpaceDimension();
  92   return std::vector<int>(_structure,_structure+_space_dim);
  93 }
  94
  95 void MEDCouplingIMesh::setOrigin(const double *originStart, const double *originStop)
  96 {
  97   checkSpaceDimension();
  98   int sz((int)std::distance(originStart,originStop));
  99   if(sz!=_space_dim)
 100     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setOrigin : input vector of origin vector has not the right size ! Or change space dimension before calling it !");
 101   std::copy(originStart,originStop,_origin);
 102   declareAsNew();
 103 }
 104
 105 std::vector<double> MEDCouplingIMesh::getOrigin() const
 106 {
 107   checkSpaceDimension();
 108   return std::vector<double>(_origin,_origin+_space_dim);
 109 }
 110
 111 void MEDCouplingIMesh::setDXYZ(const double *dxyzStart, const double *dxyzStop)
 112 {
 113   checkSpaceDimension();
 114   int sz((int)std::distance(dxyzStart,dxyzStop));
 115   if(sz!=_space_dim)
 116     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setDXYZ : input vector of dxyz vector has not the right size ! Or change space dimension before calling it !");
 117   std::copy(dxyzStart,dxyzStop,_dxyz);
 118   declareAsNew();
 119 }
 120
 121 std::vector<double> MEDCouplingIMesh::getDXYZ() const
 122 {
 123   checkSpaceDimension();
 124   return std::vector<double>(_dxyz,_dxyz+_space_dim);
 125 }
 126
 127 void MEDCouplingIMesh::setAxisUnit(const std::string& unitName)
 128 {
 129   _axis_unit=unitName;
 130   declareAsNew();
 131 }
 132
 133 std::string MEDCouplingIMesh::getAxisUnit() const
 134 {
 135   return _axis_unit;
 136 }
 137
 138 /*!
 139  * This method returns the measure of any cell in \a this.
 140  * This specific method of image grid mesh utilizes the fact that any cell in \a this have the same measure.
 141  * The value returned by this method is those used to feed the returned field in the MEDCouplingIMesh::getMeasureField.
 142  *
 143  * \sa getMeasureField
 144  */
 145 double MEDCouplingIMesh::getMeasureOfAnyCell() const
 146 {
 147   checkCoherency();
 148   int dim(getSpaceDimension());
 149   double ret(1.);
 150   for(int i=0;i<dim;i++)
 151     ret*=fabs(_dxyz[i]);
 152   return ret;
 153 }
 154
 155 /*!
 156  * This method is allows to convert \a this into MEDCouplingCMesh instance.
 157  * This method is the middle level between MEDCouplingIMesh and the most general MEDCouplingUMesh.
 158  * This method is useful for MED writers that do not have still the image grid support.
 159  *
 160  * \sa MEDCouplingMesh::buildUnstructured
 161  */
 162 MEDCouplingCMesh *MEDCouplingIMesh::convertToCartesian() const
 163 {
 164   checkCoherency();
 165   MEDCouplingAutoRefCountObjectPtr<MEDCouplingCMesh> ret(MEDCouplingCMesh::New());
 166   try
 167   { ret->copyTinyInfoFrom(this); }
 168   catch(INTERP_KERNEL::Exception& ) { }
 169   int spaceDim(getSpaceDimension());
 170   std::vector<std::string> infos(buildInfoOnComponents());
 171   for(int i=0;i<spaceDim;i++)
 172     {
 173       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr(DataArrayDouble::New()); arr->alloc(_structure[i],1); arr->setInfoOnComponent(0,infos[i]);
 174       arr->iota(); arr->applyLin(_dxyz[i],_origin[i]);
 175       ret->setCoordsAt(i,arr);
 176     }
 177   return ret.retn();
 178 }
 179
 180 /*!
 181  * This method refines \a this uniformaly along all of its dimensions. In case of success the space covered by \a this will remain
 182  * the same before the invocation except that the number of cells will be multiplied by \a factor ^ this->getMeshDimension().
 183  * The origin of \a this will be not touched only spacing and node structure will be changed.
 184  * This method can be useful for AMR users.
 185  */
 186 void MEDCouplingIMesh::refineWithFactor(int factor)
 187 {
 188   if(factor==0)
 189     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::refineWithFactor : refinement factor must be != 0 !");
 190   checkCoherency();
 191   int factAbs(std::abs(factor));
 192   double fact2(1./(double)factor);
 193   std::transform(_structure,_structure+_space_dim,_structure,std::bind2nd(std::plus<int>(),-1));
 194   std::transform(_structure,_structure+_space_dim,_structure,std::bind2nd(std::multiplies<int>(),factAbs));
 195   std::transform(_structure,_structure+_space_dim,_structure,std::bind2nd(std::plus<int>(),1));
 196   std::transform(_dxyz,_dxyz+_space_dim,_dxyz,std::bind2nd(std::multiplies<double>(),fact2));
 197   declareAsNew();
 198 }
 199
 200 /*!
 201  * This static method is useful to condense field on cells of a MEDCouplingIMesh instance coming from a refinement ( MEDCouplingIMesh::refineWithFactor for example)
 202  * to a coarse MEDCouplingIMesh instance. So this method can be seen as a specialization in P0P0 conservative interpolation non overlaping from fine image mesh
 203  * 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.
 204  *
 205  * \param [in,out] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
 206  * \param [in] coarseSt The cell structure of coarse mesh.
 207  * \param [in] fineDA The DataArray containing the cell field on uniformly refined mesh
 208  * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
 209  */
 210 void MEDCouplingIMesh::CondenseFineToCoarse(DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, const DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse)
 211 {
 212   if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
 213     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the parameters 1 or 3 are NULL or not allocated !");
 214   int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
 215   int nbCompo(fineDA->getNumberOfComponents());
 216   if(coarseDA->getNumberOfComponents()!=nbCompo)
 217     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the number of components of fine DA and coarse one mismatches !");
 218   if(meshDim!=(int)fineLocInCoarse.size())
 219     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the size of fineLocInCoarse (4th param) must be equal to the sier of coarseSt (2nd param) !");
 220   if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
 221     {
 222       std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Expecting " << nbOfTuplesInCoarseExp << " having " << coarseDA->getNumberOfTuples() << " !";
 223       throw INTERP_KERNEL::Exception(oss.str().c_str());
 224     }
 225   int nbTuplesFine(fineDA->getNumberOfTuples());
 226   if(nbTuplesFine%nbOfTuplesInCoarseExp!=0)
 227     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : Invalid nb of tuples in fine DataArray regarding its structure !");
 228   int factN(nbTuplesFine/nbOfTuplesInFineExp);
 229   int fact(FindIntRoot(factN,meshDim));
 230   // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
 231   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ids(BuildExplicitIdsFrom(coarseSt,fineLocInCoarse));
 232   const int *idsPtr(ids->begin());
 233   double *outPtr(coarseDA->getPointer());
 234   const double *inPtr(fineDA->begin());
 235   coarseDA->setPartOfValuesSimple3(0.,ids->begin(),ids->end(),0,nbCompo,1);
 236   //
 237   switch(meshDim)
 238   {
 239     case 2:
 240       {
 241         int kk(0);
 242         std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
 243         for(int it=0;it<dims[1];it++)
 244           {
 245             for(int i=0;i<fact;i++)
 246               {
 247                 for(int j=0;j<dims[0];j++,inPtr+=fact)
 248                   {
 249                     double *loc(outPtr+idsPtr[kk+j]*nbCompo);
 250                     std::transform(inPtr,inPtr+fact,loc,loc,std::plus<double>());
 251                   }
 252               }
 253             kk+=it;
 254           }
 255         break;
 256       }
 257     default:
 258       throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : only dimensions 2 supported !");
 259   }
 260 }
 261
 262 void MEDCouplingIMesh::setSpaceDimension(int spaceDim)
 263 {
 264   if(spaceDim==_space_dim)
 265     return ;
 266   CheckSpaceDimension(spaceDim);
 267   _space_dim=spaceDim;
 268   declareAsNew();
 269 }
 270
 271 void MEDCouplingIMesh::updateTime() const
 272 {
 273 }
 274
 275 std::size_t MEDCouplingIMesh::getHeapMemorySizeWithoutChildren() const
 276 {
 277   return MEDCouplingStructuredMesh::getHeapMemorySizeWithoutChildren();
 278 }
 279
 280 std::vector<const BigMemoryObject *> MEDCouplingIMesh::getDirectChildren() const
 281 {
 282   return std::vector<const BigMemoryObject *>();
 283 }
 284
 285 /*!
 286  * This method copyies all tiny strings from other (name and components name).
 287  * @throw if other and this have not same mesh type.
 288  */
 289 void MEDCouplingIMesh::copyTinyStringsFrom(const MEDCouplingMesh *other)
 290 {
 291   const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
 292   if(!otherC)
 293     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::copyTinyStringsFrom : meshes have not same type !");
 294   MEDCouplingStructuredMesh::copyTinyStringsFrom(other);
 295   declareAsNew();
 296 }
 297
 298 bool MEDCouplingIMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
 299 {
 300   if(!other)
 301     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::isEqualIfNotWhy : input other pointer is null !");
 302   const MEDCouplingIMesh *otherC(dynamic_cast<const MEDCouplingIMesh *>(other));
 303   if(!otherC)
 304     {
 305       reason="mesh given in input is not castable in MEDCouplingIMesh !";
 306       return false;
 307     }
 308   if(!MEDCouplingStructuredMesh::isEqualIfNotWhy(other,prec,reason))
 309     return false;
 310   if(!isEqualWithoutConsideringStrInternal(otherC,prec,reason))
 311     return false;
 312   if(_axis_unit!=otherC->_axis_unit)
 313     {
 314       reason="The units of axis are not the same !";
 315       return false;
 316     }
 317   return true;
 318 }
 319
 320 bool MEDCouplingIMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
 321 {
 322   const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
 323   if(!otherC)
 324     return false;
 325   std::string tmp;
 326   return isEqualWithoutConsideringStrInternal(other,prec,tmp);
 327 }
 328
 329 bool MEDCouplingIMesh::isEqualWithoutConsideringStrInternal(const MEDCouplingMesh *other, double prec, std::string& reason) const
 330 {
 331   const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
 332   if(!otherC)
 333     return false;
 334   if(_space_dim!=otherC->_space_dim)
 335     {
 336       std::ostringstream oss;
 337       oss << "The spaceDimension of this (" << _space_dim << ") is not equal to those of other (" << otherC->_space_dim << ") !";
 338       return false;
 339     }
 340   checkSpaceDimension();
 341   for(int i=0;i<_space_dim;i++)
 342     {
 343       if(fabs(_origin[i]-otherC->_origin[i])>prec)
 344         {
 345           std::ostringstream oss;
 346           oss << "The origin of this and other differs at " << i << " !";
 347           reason=oss.str();
 348           return false;
 349         }
 350     }
 351   for(int i=0;i<_space_dim;i++)
 352     {
 353       if(fabs(_dxyz[i]-otherC->_dxyz[i])>prec)
 354         {
 355           std::ostringstream oss;
 356           oss << "The delta of this and other differs at " << i << " !";
 357           reason=oss.str();
 358           return false;
 359         }
 360     }
 361   for(int i=0;i<_space_dim;i++)
 362     {
 363       if(_structure[i]!=otherC->_structure[i])
 364         {
 365           std::ostringstream oss;
 366           oss << "The structure of this and other differs at " << i << " !";
 367           reason=oss.str();
 368           return false;
 369         }
 370     }
 371   return true;
 372 }
 373
 374 void MEDCouplingIMesh::checkDeepEquivalWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
 375                                             DataArrayInt *&cellCor, DataArrayInt *&nodeCor) const
 376 {
 377   if(!isEqualWithoutConsideringStr(other,prec))
 378     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalWith : Meshes are not the same !");
 379 }
 380
 381 /*!
 382  * Nothing is done here (except to check that the other is a ParaMEDMEM::MEDCouplingIMesh instance too).
 383  * The user intend that the nodes are the same, so by construction of ParaMEDMEM::MEDCouplingIMesh, \a this and \a other are the same !
 384  */
 385 void MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
 386                                                        DataArrayInt *&cellCor) const
 387 {
 388   if(!isEqualWithoutConsideringStr(other,prec))
 389     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith : Meshes are not the same !");
 390 }
 391
 392 void MEDCouplingIMesh::checkCoherency() const
 393 {
 394   checkSpaceDimension();
 395   for(int i=0;i<_space_dim;i++)
 396     if(_structure[i]<1)
 397       {
 398         std::ostringstream oss; oss << "MEDCouplingIMesh::checkCoherency : On axis " << i << "/" << _space_dim << ", number of nodes is equal to " << _structure[i] << " ! must be >=1 !";
 399         throw INTERP_KERNEL::Exception(oss.str().c_str());
 400       }
 401 }
 402
 403 void MEDCouplingIMesh::checkCoherency1(double eps) const
 404 {
 405   checkCoherency();
 406 }
 407
 408 void MEDCouplingIMesh::checkCoherency2(double eps) const
 409 {
 410   checkCoherency1(eps);
 411 }
 412
 413 void MEDCouplingIMesh::getNodeGridStructure(int *res) const
 414 {
 415   checkSpaceDimension();
 416   std::copy(_structure,_structure+_space_dim,res);
 417 }
 418
 419 std::vector<int> MEDCouplingIMesh::getNodeGridStructure() const
 420 {
 421   checkSpaceDimension();
 422   std::vector<int> ret(_structure,_structure+_space_dim);
 423   return ret;
 424 }
 425
 426 MEDCouplingStructuredMesh *MEDCouplingIMesh::buildStructuredSubPart(const std::vector< std::pair<int,int> >& cellPart) const
 427 {
 428   checkCoherency();
 429   int dim(getSpaceDimension());
 430   if(dim!=(int)cellPart.size())
 431     {
 432       std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : the space dimension is " << dim << " and cell part size is " << cellPart.size() << " !";
 433       throw INTERP_KERNEL::Exception(oss.str().c_str());
 434     }
 435   double retOrigin[3]={0.,0.,0.};
 436   int retStruct[3]={0,0,0};
 437   MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(dynamic_cast<MEDCouplingIMesh *>(deepCpy()));
 438   for(int i=0;i<dim;i++)
 439     {
 440       int startNode(cellPart[i].first),endNode(cellPart[i].second+1);
 441       int myDelta(endNode-startNode);
 442       if(startNode<0 || startNode>=_structure[i])
 443         {
 444           std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : At dimension #" << i << " the start node id is " << startNode << " it should be in [0," << _structure[i] << ") !";
 445           throw INTERP_KERNEL::Exception(oss.str().c_str());
 446         }
 447       if(myDelta<0 || myDelta>_structure[i])
 448         {
 449           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;
 450           throw INTERP_KERNEL::Exception(oss.str().c_str());
 451         }
 452       retOrigin[i]=_origin[i]+startNode*_dxyz[i];
 453       retStruct[i]=myDelta;
 454     }
 455   ret->setNodeStruct(retStruct,retStruct+dim);
 456   ret->setOrigin(retOrigin,retOrigin+dim);
 457   ret->checkCoherency();
 458   return ret.retn();
 459 }
 460
 461 /*!
 462  * Return the space dimension of \a this.
 463  */
 464 int MEDCouplingIMesh::getSpaceDimension() const
 465 {
 466   return _space_dim;
 467 }
 468
 469 void MEDCouplingIMesh::getCoordinatesOfNode(int nodeId, std::vector<double>& coo) const
 470 {
 471   int tmp[3];
 472   int spaceDim(getSpaceDimension());
 473   getSplitNodeValues(tmp);
 474   int tmp2[3];
 475   GetPosFromId(nodeId,spaceDim,tmp,tmp2);
 476   for(int j=0;j<spaceDim;j++)
 477     coo.push_back(_origin[j]+_dxyz[j]*tmp2[j]);
 478 }
 479
 480 std::string MEDCouplingIMesh::simpleRepr() const
 481 {
 482   std::ostringstream ret;
 483   ret << "Image grid with name : \"" << getName() << "\"\n";
 484   ret << "Description of mesh : \"" << getDescription() << "\"\n";
 485   int tmpp1,tmpp2;
 486   double tt(getTime(tmpp1,tmpp2));
 487   int spaceDim(_space_dim);
 488   ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
 489   ret << "Iteration : " << tmpp1  << " Order : " << tmpp2 << "\n";
 490   ret << "Space dimension : " << spaceDim << "\n";
 491   if(spaceDim<0 || spaceDim>3)
 492     return ret.str();
 493   ret << "The nodal structure is : "; std::copy(_structure,_structure+spaceDim,std::ostream_iterator<int>(ret," ")); ret << "\n";
 494   ret << "The origin position is [" << _axis_unit << "]: ";
 495   std::copy(_origin,_origin+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
 496   ret << "The intervals along axis are : ";
 497   std::copy(_dxyz,_dxyz+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
 498   return ret.str();
 499 }
 500
 501 std::string MEDCouplingIMesh::advancedRepr() const
 502 {
 503   return simpleRepr();
 504 }
 505
 506 void MEDCouplingIMesh::getBoundingBox(double *bbox) const
 507 {
 508   checkCoherency();
 509   int dim(getSpaceDimension());
 510   for(int idim=0; idim<dim; idim++)
 511     {
 512       bbox[2*idim]=_origin[idim];
 513       bbox[2*idim+1]=_origin[idim]+_dxyz[idim]*_structure[idim];
 514     }
 515 }
 516
 517 /*!
 518  * Returns a new MEDCouplingFieldDouble containing volumes of cells constituting \a this
 519  * mesh.<br>
 520  * For 1D cells, the returned field contains lengths.<br>
 521  * For 2D cells, the returned field contains areas.<br>
 522  * For 3D cells, the returned field contains volumes.
 523  *  \param [in] isAbs - a not used parameter.
 524  *  \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble on cells
 525  *         and one time . The caller is to delete this field using decrRef() as it is no
 526  *         more needed.
 527  */
 528 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureField(bool isAbs) const
 529 {
 530   checkCoherency();
 531   std::string name="MeasureOfMesh_";
 532   name+=getName();
 533   int nbelem(getNumberOfCells());
 534   MEDCouplingFieldDouble *field(MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME));
 535   field->setName(name);
 536   DataArrayDouble* array(DataArrayDouble::New());
 537   array->alloc(nbelem,1);
 538   array->fillWithValue(getMeasureOfAnyCell());
 539   field->setArray(array) ;
 540   array->decrRef();
 541   field->setMesh(const_cast<MEDCouplingIMesh *>(this));
 542   field->synchronizeTimeWithMesh();
 543   return field;
 544 }
 545
 546 /*!
 547  * not implemented yet !
 548  */
 549 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureFieldOnNode(bool isAbs) const
 550 {
 551   throw INTERP_KERNEL::Exception("MEDCouplingIMesh::getMeasureFieldOnNode : not implemented yet !");
 552   //return 0;
 553 }
 554
 555 int MEDCouplingIMesh::getCellContainingPoint(const double *pos, double eps) const
 556 {
 557   int dim(getSpaceDimension()),ret(0),coeff(1);
 558   for(int i=0;i<dim;i++)
 559     {
 560       int nbOfCells(_structure[i]-1);
 561       double ref(pos[i]);
 562       int tmp((int)((ref-_origin[i])/_dxyz[i]));
 563       if(tmp>=0 && tmp<nbOfCells)
 564         {
 565           ret+=coeff*tmp;
 566           coeff*=nbOfCells;
 567         }
 568       else
 569         return -1;
 570     }
 571   return ret;
 572 }
 573
 574 void MEDCouplingIMesh::rotate(const double *center, const double *vector, double angle)
 575 {
 576   throw INTERP_KERNEL::Exception("No rotation available on IMesh : Traduce it to unstructured mesh to apply it !");
 577 }
 578
 579 /*!
 580  * Translates all nodes of \a this mesh by a given vector. Actually, it adds each
 581  * component of the \a vector to all node coordinates of a corresponding axis.
 582  *  \param [in] vector - the translation vector whose size must be not less than \a
 583  *         this->getSpaceDimension().
 584  */
 585 void MEDCouplingIMesh::translate(const double *vector)
 586 {
 587   checkSpaceDimension();
 588   int dim(getSpaceDimension());
 589   std::transform(_origin,_origin+dim,vector,_origin,std::plus<double>());
 590   declareAsNew();
 591 }
 592
 593 /*!
 594  * Applies scaling transformation to all nodes of \a this mesh.
 595  *  \param [in] point - coordinates of a scaling center. This array is to be of
 596  *         size \a this->getSpaceDimension() at least.
 597  *  \param [in] factor - a scale factor.
 598  */
 599 void MEDCouplingIMesh::scale(const double *point, double factor)
 600 {
 601   checkSpaceDimension();
 602   int dim(getSpaceDimension());
 603   std::transform(_origin,_origin+dim,point,_origin,std::minus<double>());
 604   std::transform(_origin,_origin+dim,_origin,std::bind2nd(std::multiplies<double>(),factor));
 605   std::transform(_dxyz,_dxyz+dim,_dxyz,std::bind2nd(std::multiplies<double>(),factor));
 606   std::transform(_origin,_origin+dim,point,_origin,std::plus<double>());
 607   declareAsNew();
 608 }
 609
 610 MEDCouplingMesh *MEDCouplingIMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
 611 {
 612   //not implemented yet !
 613   return 0;
 614 }
 615
 616 /*!
 617  * Returns a new DataArrayDouble holding coordinates of all nodes of \a this mesh.
 618  *  \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
 619  *          this->getNumberOfNodes() tuples per \a this->getSpaceDimension()
 620  *          components. The caller is to delete this array using decrRef() as it is
 621  *          no more needed.
 622  */
 623 DataArrayDouble *MEDCouplingIMesh::getCoordinatesAndOwner() const
 624 {
 625   checkCoherency();
 626   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
 627   int spaceDim(getSpaceDimension()),nbNodes(getNumberOfNodes());
 628   ret->alloc(nbNodes,spaceDim);
 629   double *pt(ret->getPointer());
 630   ret->setInfoOnComponents(buildInfoOnComponents());
 631   int tmp2[3],tmp[3];
 632   getSplitNodeValues(tmp);
 633   for(int i=0;i<nbNodes;i++)
 634     {
 635       GetPosFromId(i,spaceDim,tmp,tmp2);
 636       for(int j=0;j<spaceDim;j++)
 637         pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+_origin[j];
 638     }
 639   return ret.retn();
 640 }
 641
 642 /*!
 643  * Returns a new DataArrayDouble holding barycenters of all cells. The barycenter is
 644  * computed by averaging coordinates of cell nodes.
 645  *  \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
 646  *          this->getNumberOfCells() tuples per \a this->getSpaceDimension()
 647  *          components. The caller is to delete this array using decrRef() as it is
 648  *          no more needed.
 649  */
 650 DataArrayDouble *MEDCouplingIMesh::getBarycenterAndOwner() const
 651 {
 652   checkCoherency();
 653   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
 654   int spaceDim(getSpaceDimension()),nbCells(getNumberOfCells()),tmp[3],tmp2[3];
 655   ret->alloc(nbCells,spaceDim);
 656   double *pt(ret->getPointer()),shiftOrigin[3];
 657   std::transform(_dxyz,_dxyz+spaceDim,shiftOrigin,std::bind2nd(std::multiplies<double>(),0.5));
 658   std::transform(_origin,_origin+spaceDim,shiftOrigin,shiftOrigin,std::plus<double>());
 659   getSplitCellValues(tmp);
 660   ret->setInfoOnComponents(buildInfoOnComponents());
 661   for(int i=0;i<nbCells;i++)
 662     {
 663       GetPosFromId(i,spaceDim,tmp,tmp2);
 664       for(int j=0;j<spaceDim;j++)
 665         pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+shiftOrigin[j];
 666     }
 667   return ret.retn();
 668 }
 669
 670 DataArrayDouble *MEDCouplingIMesh::computeIsoBarycenterOfNodesPerCell() const
 671 {
 672   return MEDCouplingIMesh::getBarycenterAndOwner();
 673 }
 674
 675 void MEDCouplingIMesh::renumberCells(const int *old2NewBg, bool check)
 676 {
 677   throw INTERP_KERNEL::Exception("Functionnality of renumbering cell not available for IMesh !");
 678 }
 679
 680 void MEDCouplingIMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
 681 {
 682   int it,order;
 683   double time(getTime(it,order));
 684   tinyInfo.clear();
 685   tinyInfoD.clear();
 686   littleStrings.clear();
 687   littleStrings.push_back(getName());
 688   littleStrings.push_back(getDescription());
 689   littleStrings.push_back(getTimeUnit());
 690   littleStrings.push_back(getAxisUnit());
 691   tinyInfo.push_back(it);
 692   tinyInfo.push_back(order);
 693   tinyInfo.push_back(_space_dim);
 694   tinyInfo.insert(tinyInfo.end(),_structure,_structure+3);
 695   tinyInfoD.push_back(time);
 696   tinyInfoD.insert(tinyInfoD.end(),_dxyz,_dxyz+3);
 697   tinyInfoD.insert(tinyInfoD.end(),_origin,_origin+3);
 698 }
 699
 700 void MEDCouplingIMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
 701 {
 702   a1->alloc(0,1);
 703   a2->alloc(0,1);
 704 }
 705
 706 void MEDCouplingIMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
 707 {
 708   a1=DataArrayInt::New();
 709   a1->alloc(0,1);
 710   a2=DataArrayDouble::New();
 711   a2->alloc(0,1);
 712 }
 713
 714 void MEDCouplingIMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
 715                                        const std::vector<std::string>& littleStrings)
 716 {
 717   setName(littleStrings[0]);
 718   setDescription(littleStrings[1]);
 719   setTimeUnit(littleStrings[2]);
 720   setAxisUnit(littleStrings[3]);
 721   setTime(tinyInfoD[0],tinyInfo[0],tinyInfo[1]);
 722   _space_dim=tinyInfo[2];
 723   _structure[0]=tinyInfo[3]; _structure[1]=tinyInfo[4]; _structure[2]=tinyInfo[5];
 724   _dxyz[0]=tinyInfoD[1]; _dxyz[1]=tinyInfoD[2]; _dxyz[2]=tinyInfoD[3];
 725   _origin[0]=tinyInfoD[4]; _origin[1]=tinyInfoD[5]; _origin[2]=tinyInfoD[6];
 726   declareAsNew();
 727 }
 728
 729 void MEDCouplingIMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const
 730 {
 731   checkCoherency();
 732   std::ostringstream extent,origin,spacing;
 733   for(int i=0;i<3;i++)
 734     {
 735       if(i<_space_dim)
 736         { extent << "0 " <<  _structure[i]-1 << " "; origin << _origin[i] << " "; spacing << _dxyz[i] << " "; }
 737       else
 738         { extent << "0 0 "; origin << "0 "; spacing << "0 "; }
 739     }
 740   ofs << "  <" << getVTKDataSetType() << " WholeExtent=\"" << extent.str() << "\" Origin=\"" << origin.str() << "\" Spacing=\"" << spacing.str() << "\">\n";
 741   ofs << "    <Piece Extent=\"" << extent.str() << "\">\n";
 742   ofs << "      <PointData>\n" << pointData << std::endl;
 743   ofs << "      </PointData>\n";
 744   ofs << "      <CellData>\n" << cellData << std::endl;
 745   ofs << "      </CellData>\n";
 746   ofs << "      <Coordinates>\n";
 747   ofs << "      </Coordinates>\n";
 748   ofs << "    </Piece>\n";
 749   ofs << "  </" << getVTKDataSetType() << ">\n";
 750 }
 751
 752 void MEDCouplingIMesh::reprQuickOverview(std::ostream& stream) const
 753 {
 754   stream << "MEDCouplingIMesh C++ instance at " << this << ". Name : \"" << getName() << "\". Space dimension : " << _space_dim << ".";
 755   if(_space_dim<0 || _space_dim>3)
 756     return ;
 757   stream << "\n";
 758   std::ostringstream stream0,stream1;
 759   int nbNodes(1),nbCells(0);
 760   bool isPb(false);
 761   for(int i=0;i<_space_dim;i++)
 762     {
 763       char tmp('X'+i);
 764       int tmpNodes(_structure[i]);
 765       stream1 << "- Axis " << tmp << " : " << tmpNodes << " nodes (orig=" << _origin[i] << ", inter=" << _dxyz[i] << ").";
 766       if(i!=_space_dim-1)
 767         stream1 << std::endl;
 768       if(tmpNodes>=1)
 769         nbNodes*=tmpNodes;
 770       else
 771         isPb=true;
 772       if(tmpNodes>=2)
 773         nbCells=nbCells==0?tmpNodes-1:nbCells*(tmpNodes-1);
 774     }
 775   if(!isPb)
 776     {
 777       stream0 << "Number of cells : " << nbCells << ", Number of nodes : " << nbNodes;
 778       stream << stream0.str();
 779       if(_space_dim>0)
 780         stream << std::endl;
 781     }
 782   stream << stream1.str();
 783 }
 784
 785 std::string MEDCouplingIMesh::getVTKDataSetType() const
 786 {
 787   return std::string("ImageData");
 788 }
 789
 790 std::vector<std::string> MEDCouplingIMesh::buildInfoOnComponents() const
 791 {
 792   checkSpaceDimension();
 793   int dim(getSpaceDimension());
 794   std::vector<std::string> ret(dim);
 795   for(int i=0;i<dim;i++)
 796     {
 797       std::ostringstream oss;
 798       char tmp('X'+i); oss << tmp;
 799       ret[i]=DataArray::BuildInfoFromVarAndUnit(oss.str(),_axis_unit);
 800     }
 801   return ret;
 802 }
 803
 804 void MEDCouplingIMesh::checkSpaceDimension() const
 805 {
 806   CheckSpaceDimension(_space_dim);
 807 }
 808
 809 void MEDCouplingIMesh::CheckSpaceDimension(int spaceDim)
 810 {
 811   if(spaceDim<0 || spaceDim>3)
 812     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CheckSpaceDimension : input spaceDim must be in [0,1,2,3] !");
 813 }
 814
 815 int MEDCouplingIMesh::FindIntRoot(int val, int order)
 816 {
 817   if(order==0)
 818     return 1;
 819   if(val<0)
 820     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : input val is < 0 ! Not possible to compute a root !");
 821   if(order==1)
 822     return val;
 823   if(order!=2 && order!=3)
 824     throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the order available are 0,1,2 or 3 !");
 825   double valf((double)val);
 826   if(order==2)
 827     {
 828       double retf(sqrt(valf));
 829       int ret((int)retf);
 830       if(ret*ret!=val)
 831         throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect square root !");
 832       return ret;
 833     }
 834   else//order==3
 835     {
 836       double retf(std::pow(val,0.3333333333333333));
 837       int ret((int)round(retf));
 838       if(ret*ret*ret!=val)
 839         throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect cublic root !");
 840       return ret;
 841     }
 842 }