src/MEDCoupling/MEDCouplingFieldDouble.cxx

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay (CEA/DEN)
  20
  21 #include "MEDCouplingFieldDouble.hxx"
  22 #include "MEDCouplingFieldTemplate.hxx"
  23 #include "MEDCouplingUMesh.hxx"
  24 #include "MEDCouplingTimeDiscretization.hxx"
  25 #include "MEDCouplingFieldDiscretization.hxx"
  26 #include "MEDCouplingAutoRefCountObjectPtr.hxx"
  27 #include "MEDCouplingNatureOfField.hxx"
  28
  29 #include "InterpKernelAutoPtr.hxx"
  30
  31 #include <sstream>
  32 #include <limits>
  33 #include <algorithm>
  34 #include <functional>
  35
  36 using namespace ParaMEDMEM;
  37
  38 /*!
  39  * Creates a new instance of MEDCouplingFieldDouble of given type. The caller is responsable for the returned field.
  40  *
  41  * \param [in] type type of spatial discretization of a created field (\ref ParaMEDMEM::ON_CELLS "ON_CELLS", \ref ParaMEDMEM::ON_NODES "ON_NODES", \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT", \ref ParaMEDMEM::ON_GAUSS_NE "ON_GAUSS_NE", \ref ParaMEDMEM::ON_NODES_KR "ON_NODES_KR").
  42  * \param [in] td type of time discretization of a created field (\ref ParaMEDMEM::NO_TIME "NO_TIME", \ref ParaMEDMEM::ONE_TIME "ONE_TIME", \ref ParaMEDMEM::LINEAR_TIME "LINEAR_TIME", \ref ParaMEDMEM::CONST_ON_TIME_INTERVAL "CONST_ON_TIME_INTERVAL").
  43
  44  * \return a newly allocated field the caller should deal with.
  45  */
  46 MEDCouplingFieldDouble *MEDCouplingFieldDouble::New(TypeOfField type, TypeOfTimeDiscretization td)
  47 {
  48   return new MEDCouplingFieldDouble(type,td);
  49 }
  50
  51 /*!
  52  * Creates a new instance of MEDCouplingFieldDouble of given type. The caller is responsable for the returned field.
  53  * ** WARINING : This method do not deeply copy neither mesh nor spatial discretization. Only a shallow copy (reference) is done for mesh and spatial discretization ! **
  54  *
  55  * \param [in] ft \ref MEDCouplingFieldTemplatesPage "field template" defining its spatial discretization and supporting mesh.
  56  * \param [in] td type of time discretization of a created field (\ref ParaMEDMEM::NO_TIME "NO_TIME", \ref ParaMEDMEM::ONE_TIME "ONE_TIME", \ref ParaMEDMEM::LINEAR_TIME "LINEAR_TIME", \ref ParaMEDMEM::CONST_ON_TIME_INTERVAL "CONST_ON_TIME_INTERVAL")
  57
  58  * \return a newly allocated field the caller should deal with.
  59  */
  60 MEDCouplingFieldDouble *MEDCouplingFieldDouble::New(const MEDCouplingFieldTemplate& ft, TypeOfTimeDiscretization td)
  61 {
  62   return new MEDCouplingFieldDouble(ft,td);
  63 }
  64
  65 /*!
  66  * Sets time \a unit of \a this field.
  67  *
  68  * \param [in] unit \a unit (string) in which time is measured.
  69  */
  70 void MEDCouplingFieldDouble::setTimeUnit(const char *unit)
  71 {
  72   _time_discr->setTimeUnit(unit);
  73 }
  74
  75 /*!
  76  * Returns a time unit of \a this field.
  77  *
  78  * \return a string describing units in which time is measured.
  79  */
  80 const char *MEDCouplingFieldDouble::getTimeUnit() const
  81 {
  82   return _time_discr->getTimeUnit();
  83 }
  84
  85 /*!
  86  * This method if possible the time information (time unit, time iteration, time unit and time value) with its support
  87  * that is to say its mesh.
  88  *
  89  * \throw  If \c this->_mesh is null an exception will be thrown. An exception will also be throw if the spatial discretization is
  90  *         NO_TIME.
  91  */
  92 void MEDCouplingFieldDouble::synchronizeTimeWithSupport() throw(INTERP_KERNEL::Exception)
  93 {
  94   _time_discr->synchronizeTimeWith(_mesh);
  95 }
  96
  97 /*!
  98  * This method performs a copy of \a this **without any copy of the underlying mesh** ( see warning section of this method).
  99  * The copy of arrays is deep if \b recDeepCpy equals to true, no copy of arrays is done if \b recDeepCpy equals to false.
 100  *
 101  * \c clone(false) is rather dedicated for advanced users that want to limit the amount of memory.
 102  *
 103  * It allows the user to perform methods
 104  * MEDCouplingFieldDouble::AddFields, MEDCouplingFieldDouble::MultiplyFields with \a this and the returned field.
 105  *
 106  * \warning The \b underlying \b mesh of the returned field is \b always the same (same pointer) than \a this **whatever the value** of \a recDeepCpy parameter.
 107  * If the user wants to duplicated deeply the underlying mesh he should call MEDCouplingFieldDouble::cloneWithMesh method or MEDCouplingFieldDouble::deepCpy instead.
 108  *
 109  * \param [in] recDeepCpy specifies if underlying arrays in \a this should be copied or only attached to the returned field.
 110  * \return a newly allocated MEDCouplingFieldDouble instance that the caller should deal with.
 111  * \sa ParaMEDMEM::MEDCouplingFieldDouble::cloneWithMesh(bool recDeepCpy) const
 112  */
 113 MEDCouplingFieldDouble *MEDCouplingFieldDouble::clone(bool recDeepCpy) const
 114 {
 115   return new MEDCouplingFieldDouble(*this,recDeepCpy);
 116 }
 117
 118 /*!
 119  * This method behaves exactly like MEDCouplingFieldDouble::clone method **except that here the underlying mesh is systematically **
 120  * (whatever the value of the input parameter \a recDeepCpy) **deeply duplicated**.
 121  *
 122  * The result of \c cloneWithMesh(true) is exactly the same than calling \ref MEDCouplingFieldDouble::deepCpy "deepCpy".
 123  *
 124  * So the resulting field of this call cannot be called with \a this with the following methods MEDCouplingFieldDouble::AddFields, MEDCouplingFieldDouble::MultiplyFields ...
 125  * To avoid to deep copy the underlying mesh the user should call MEDCouplingFieldDouble::clone method instead.
 126
 127  * \param [in] recDeepCpy specifies if underlying arrays in \a this should be copied or only attached to the returned field.
 128  * \return a newly allocated MEDCouplingFieldDouble instance that the caller should deal with.
 129  * \sa ParaMEDMEM::MEDCouplingFieldDouble::clone(bool recDeepCpy) const
 130  */
 131 MEDCouplingFieldDouble *MEDCouplingFieldDouble::cloneWithMesh(bool recDeepCpy) const
 132 {
 133   MEDCouplingFieldDouble *ret=clone(recDeepCpy);
 134   if(_mesh)
 135     {
 136       MEDCouplingMesh *mCpy=_mesh->deepCpy();
 137       ret->setMesh(mCpy);
 138       mCpy->decrRef();
 139     }
 140   return ret;
 141 }
 142
 143 /*!
 144  * This method performs a deepCpy of \a this (**mesh included**)!
 145  * So the resulting field of this call cannot be called with \a this with following methods MEDCouplingFieldDouble::AddFields, MEDCouplingFieldDouble::MultiplyFields ...
 146  * To avoid deep copying the underlying mesh the user should call MEDCouplingFieldDouble::clone method instead.
 147  * This method is exactly equivalent to MEDCouplingFieldDouble::cloneWithMesh called with parameter true.
 148  *
 149  * \return a newly allocated MEDCouplingFieldDouble instance that the caller should deal with.
 150  * \sa ParaMEDMEM::MEDCouplingFieldDouble::cloneWithMesh(bool recDeepCpy) const
 151  */
 152 MEDCouplingFieldDouble *MEDCouplingFieldDouble::deepCpy() const
 153 {
 154   return cloneWithMesh(true);
 155 }
 156
 157 /*!
 158  * TODOC
 159  *
 160  * \param [in] td type of time discretization of a created field (\ref ParaMEDMEM::NO_TIME "NO_TIME", \ref ParaMEDMEM::ONE_TIME "ONE_TIME", \ref ParaMEDMEM::LINEAR_TIME "LINEAR_TIME", \ref ParaMEDMEM::CONST_ON_TIME_INTERVAL "CONST_ON_TIME_INTERVAL").
 161  * \param [in] deepCopy specifies if underlying arrays in \a this should be copied or only attached to the returned field.
 162  * \return a newly allocated MEDCouplingFieldDouble instance that the caller should deal with.
 163  *
 164  * \ref cpp_mcfielddouble_buildnewtimereprfromthis "Here a C++ example."
 165
 166  * \ref py_mcfielddouble_buildnewtimereprfromthis "Here a Python example."
 167  * \sa ParaMEDMEM::MEDCouplingFieldDouble::clone(bool recDeepCpy) const
 168  */
 169 MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildNewTimeReprFromThis(TypeOfTimeDiscretization td, bool deepCopy) const
 170 {
 171   MEDCouplingTimeDiscretization *tdo=_time_discr->buildNewTimeReprFromThis(td,deepCopy);
 172   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),tdo,_type->clone());
 173   ret->setMesh(getMesh());
 174   ret->setName(getName());
 175   ret->setDescription(getDescription());
 176   return ret;
 177 }
 178
 179 /*!
 180  * Copy tiny info (component names, name, description) but warning the underlying mesh is not renamed (for safety reason).
 181  */
 182 void MEDCouplingFieldDouble::copyTinyStringsFrom(const MEDCouplingField *other) throw(INTERP_KERNEL::Exception)
 183 {
 184   MEDCouplingField::copyTinyStringsFrom(other);
 185   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
 186   if(otherC)
 187     {
 188       _time_discr->copyTinyStringsFrom(*otherC->_time_discr);
 189     }
 190 }
 191
 192 /*!
 193  * Copy only times, order, iteration from other. The underlying mesh is not impacted by this method.
 194  * Arrays are not impacted too.
 195  */
 196 void MEDCouplingFieldDouble::copyTinyAttrFrom(const MEDCouplingFieldDouble *other) throw(INTERP_KERNEL::Exception)
 197 {
 198   if(other)
 199     {
 200       _time_discr->copyTinyAttrFrom(*other->_time_discr);
 201     }
 202
 203 }
 204
 205 void MEDCouplingFieldDouble::copyAllTinyAttrFrom(const MEDCouplingFieldDouble *other) throw(INTERP_KERNEL::Exception)
 206 {
 207   copyTinyStringsFrom(other);
 208   copyTinyAttrFrom(other);
 209 }
 210
 211 std::string MEDCouplingFieldDouble::simpleRepr() const
 212 {
 213   std::ostringstream ret;
 214   ret << "FieldDouble with name : \"" << getName() << "\"\n";
 215   ret << "Description of field is : \"" << getDescription() << "\"\n";
 216   ret << "FieldDouble space discretization is : " << _type->getStringRepr() << "\n";
 217   ret << "FieldDouble time discretization is : " << _time_discr->getStringRepr() << "\n";
 218   ret << "FieldDouble nature of field is : " << MEDCouplingNatureOfField::GetRepr(_nature) << "\n";
 219   if(getArray())
 220     {
 221       int nbOfCompo=getArray()->getNumberOfComponents();
 222       ret << "FieldDouble default array has " << nbOfCompo << " components and " << getArray()->getNumberOfTuples() << " tuples.\n";
 223       ret << "FieldDouble default array has following info on components : ";
 224       for(int i=0;i<nbOfCompo;i++)
 225         ret << "\"" << getArray()->getInfoOnComponent(i) << "\" ";
 226       ret << "\n";
 227     }
 228   if(_mesh)
 229     ret << "Mesh support information :\n__________________________\n" << _mesh->simpleRepr();
 230   else
 231     ret << "Mesh support information : No mesh set !\n";
 232   return ret.str();
 233 }
 234
 235 std::string MEDCouplingFieldDouble::advancedRepr() const
 236 {
 237   std::ostringstream ret;
 238   ret << "FieldDouble with name : \"" << getName() << "\"\n";
 239   ret << "Description of field is : \"" << getDescription() << "\"\n";
 240   ret << "FieldDouble space discretization is : " << _type->getStringRepr() << "\n";
 241   ret << "FieldDouble time discretization is : " << _time_discr->getStringRepr() << "\n";
 242   if(getArray())
 243     ret << "FieldDouble default array has " << getArray()->getNumberOfComponents() << " components and " << getArray()->getNumberOfTuples() << " tuples.\n";
 244   if(_mesh)
 245     ret << "Mesh support information :\n__________________________\n" << _mesh->advancedRepr();
 246   else
 247     ret << "Mesh support information : No mesh set !\n";
 248   std::vector<DataArrayDouble *> arrays;
 249   _time_discr->getArrays(arrays);
 250   int arrayId=0;
 251   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++,arrayId++)
 252     {
 253       ret << "Array #" << arrayId << " :\n__________\n";
 254       if(*iter)
 255         (*iter)->reprWithoutNameStream(ret);
 256       else
 257         ret << "Array empty !";
 258       ret << "\n";
 259     }
 260   return ret.str();
 261 }
 262
 263 void MEDCouplingFieldDouble::writeVTK(const char *fileName) const throw(INTERP_KERNEL::Exception)
 264 {
 265   std::vector<const MEDCouplingFieldDouble *> fs(1,this);
 266   MEDCouplingFieldDouble::WriteVTK(fileName,fs);
 267 }
 268
 269 bool MEDCouplingFieldDouble::isEqualIfNotWhy(const MEDCouplingField *other, double meshPrec, double valsPrec, std::string& reason) const throw(INTERP_KERNEL::Exception)
 270 {
 271   if(!other)
 272     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::isEqualIfNotWhy : other instance is NULL !");
 273   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
 274   if(!otherC)
 275     {
 276       reason="field given in input is not castable in MEDCouplingFieldDouble !";
 277       return false;
 278     }
 279   if(!MEDCouplingField::isEqualIfNotWhy(other,meshPrec,valsPrec,reason))
 280     return false;
 281   if(!_time_discr->isEqualIfNotWhy(otherC->_time_discr,valsPrec,reason))
 282     {
 283       reason.insert(0,"In FieldDouble time discretizations differ :");
 284       return false;
 285     }
 286   return true;
 287 }
 288
 289 bool MEDCouplingFieldDouble::isEqualWithoutConsideringStr(const MEDCouplingField *other, double meshPrec, double valsPrec) const
 290 {
 291   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
 292   if(!otherC)
 293     return false;
 294   if(!MEDCouplingField::isEqualWithoutConsideringStr(other,meshPrec,valsPrec))
 295     return false;
 296   if(!_time_discr->isEqualWithoutConsideringStr(otherC->_time_discr,valsPrec))
 297     return false;
 298   return true;
 299 }
 300
 301 /*!
 302  * This method states if \a this and 'other' are compatibles each other before performing any treatment.
 303  * This method is good for methods like : mergeFields.
 304  * This method is not very demanding compared to areStrictlyCompatible that is better for operation on fields.
 305  */
 306 bool MEDCouplingFieldDouble::areCompatibleForMerge(const MEDCouplingField *other) const
 307 {
 308   if(!MEDCouplingField::areCompatibleForMerge(other))
 309     return false;
 310   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
 311   if(!otherC)
 312     return false;
 313   if(!_time_discr->areCompatible(otherC->_time_discr))
 314     return false;
 315   return true;
 316 }
 317
 318 /*!
 319  * This method is more strict than MEDCouplingField::areCompatibleForMerge method.
 320  * This method is used for operation on fields to operate a first check before attempting operation.
 321  */
 322 bool MEDCouplingFieldDouble::areStrictlyCompatible(const MEDCouplingField *other) const
 323 {
 324   std::string tmp;
 325   if(!MEDCouplingField::areStrictlyCompatible(other))
 326     return false;
 327   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
 328   if(!otherC)
 329     return false;
 330   if(!_time_discr->areStrictlyCompatible(otherC->_time_discr,tmp))
 331     return false;
 332   return true;
 333 }
 334
 335 /*!
 336  * Method with same principle than MEDCouplingFieldDouble::areStrictlyCompatible method except that
 337  * number of components between \a this and 'other' can be different here (for operator*).
 338  */
 339 bool MEDCouplingFieldDouble::areCompatibleForMul(const MEDCouplingField *other) const
 340 {
 341   if(!MEDCouplingField::areStrictlyCompatible(other))
 342     return false;
 343   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
 344   if(!otherC)
 345     return false;
 346   if(!_time_discr->areStrictlyCompatibleForMul(otherC->_time_discr))
 347     return false;
 348   return true;
 349 }
 350
 351 /*!
 352  * Method with same principle than MEDCouplingFieldDouble::areStrictlyCompatible method except that
 353  * number of components between \a this and 'other' can be different here (for operator/).
 354  */
 355 bool MEDCouplingFieldDouble::areCompatibleForDiv(const MEDCouplingField *other) const
 356 {
 357   if(!MEDCouplingField::areStrictlyCompatible(other))
 358     return false;
 359   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
 360   if(!otherC)
 361     return false;
 362   if(!_time_discr->areStrictlyCompatibleForDiv(otherC->_time_discr))
 363     return false;
 364   return true;
 365 }
 366
 367 /*!
 368  * This method is invocated before any attempt of melding. This method is very close to areStrictlyCompatible,
 369  * except that \a this and other can have different number of components.
 370  */
 371 bool MEDCouplingFieldDouble::areCompatibleForMeld(const MEDCouplingFieldDouble *other) const
 372 {
 373   if(!MEDCouplingField::areStrictlyCompatible(other))
 374     return false;
 375   if(!_time_discr->areCompatibleForMeld(other->_time_discr))
 376     return false;
 377   return true;
 378 }
 379
 380 /*!
 381  * This method performs a clone of mesh and a renumbering of underlying cells of it. The number of cells remains the same.
 382  * The values of field are impacted in consequence to have the same geometrical field.
 383  */
 384 void MEDCouplingFieldDouble::renumberCells(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
 385 {
 386   renumberCellsWithoutMesh(old2NewBg,check);
 387   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> m=_mesh->deepCpy();
 388   m->renumberCells(old2NewBg,check);
 389   setMesh(m);
 390   updateTime();
 391 }
 392
 393 /*!
 394  * \b WARNING : use this method with lot of care !
 395  * This method performs half job of MEDCouplingFieldDouble::renumberCells. That is to say no permutation of cells is done on underlying mesh.
 396  * That is to say, the field content is changed by this method. The reason of this method is only for multi-field instances lying on the same mesh to
 397  * avoid a systematic duplication and renumbering of _mesh attribute.
 398  */
 399 void MEDCouplingFieldDouble::renumberCellsWithoutMesh(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
 400 {
 401    if(!_mesh)
 402     throw INTERP_KERNEL::Exception("Expecting a defined mesh to be able to operate a renumbering !");
 403   //
 404   _type->renumberCells(old2NewBg,check);
 405   std::vector<DataArrayDouble *> arrays;
 406   _time_discr->getArrays(arrays);
 407   _type->renumberArraysForCell(_mesh,arrays,old2NewBg,check);
 408   //
 409   updateTime();
 410 }
 411
 412 /*!
 413  * This method performs a clone of mesh and a renumbering of underlying nodes of it. The number of nodes remains not compulsory the same as renumberCells method.
 414  * The values of field are impacted in consequence to have the same geometrical field.
 415  */
 416 void MEDCouplingFieldDouble::renumberNodes(const int *old2NewBg) throw(INTERP_KERNEL::Exception)
 417 {
 418   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
 419   if(!meshC)
 420     throw INTERP_KERNEL::Exception("Invalid mesh to apply renumberNodes on it !");
 421   int nbOfNodes=meshC->getNumberOfNodes();
 422   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
 423   int newNbOfNodes=*std::max_element(old2NewBg,old2NewBg+nbOfNodes)+1;
 424   renumberNodesWithoutMesh(old2NewBg,newNbOfNodes);
 425   meshC2->renumberNodes(old2NewBg,newNbOfNodes);
 426   setMesh(meshC2);
 427 }
 428
 429 /*!
 430  * \b WARNING : use this method with lot of care !
 431  * This method performs half job of MEDCouplingFieldDouble::renumberNodes. That is to say no permutation of cells is done on underlying mesh.
 432  * That is to say, the field content is changed by this method.
 433  */
 434 void MEDCouplingFieldDouble::renumberNodesWithoutMesh(const int *old2NewBg, int newNbOfNodes, double eps) throw(INTERP_KERNEL::Exception)
 435 {
 436   std::vector<DataArrayDouble *> arrays;
 437   _time_discr->getArrays(arrays);
 438   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
 439     if(*iter)
 440       _type->renumberValuesOnNodes(eps,old2NewBg,newNbOfNodes,*iter);
 441 }
 442
 443 /*!
 444  * This method makes the assumption that the default array is set. If not an exception will be thrown.
 445  * This method is usable only if the default array has exactly one component. If not an exception will be thrown too.
 446  * This method returns all tuples ids that fit the range [vmin,vmax].
 447  * The caller has the responsability of the returned DataArrayInt.
 448  */
 449 DataArrayInt *MEDCouplingFieldDouble::getIdsInRange(double vmin, double vmax) const throw(INTERP_KERNEL::Exception)
 450 {
 451   if(getArray()==0)
 452     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getIdsInRange : no default array set !");
 453   return getArray()->getIdsInRange(vmin,vmax);
 454 }
 455
 456 /*!
 457  * Builds a newly created field, that the caller will have the responsability to deal with (decrRef).
 458  * This method makes the assumption that the field is correctly defined when this method is called, no check of this will be done.
 459  * This method returns a restriction of \a this so that only tuples id specified in 'part' will be contained in returned field.
 460  * Parameter 'part' specifies \b cell \b ids \b whatever \b the \b spatial \b discretization of \a this (ON_CELLS, ON_NODES, ON_GAUSS_PT, ON_GAUSS_NE)
 461  *
 462  * If \a this is a field on cell lying on a mesh that have 10 cells. If part contains following cellIds [3,7,6].
 463  * In this case the returned field will lie on mesh having 3 cells and the returned field will contain 3 tuples.
 464  * Tuple#0 of return field will refer to the cell#0 of returned mesh. The cell #0 of returned mesh will be equal to the cell#3 of 'this->getMesh()'
 465  * Tuple#1 of return field will refer to the cell#1 of returned mesh. The cell #1 of returned mesh will be equal to the cell#7 of 'this->getMesh()'
 466  * Tuple#2 of return field will refer to the cell#2 of returned mesh. The cell #2 of returned mesh will be equal to the cell#6 of 'this->getMesh()'
 467  *
 468  * If \a this is field on node lying on a mesh that have 10 cells and 11 nodes for example. If part contains following cellIds [3,7,6].
 469  * \a this is currently contains 11 tuples. If the restriction of mesh to 3 cells leads to a mesh with 6 nodes, the returned field,
 470  * will contain 6 tuples and this field will lie on this restricted mesh.
 471  */
 472 MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPart(const DataArrayInt *part) const throw(INTERP_KERNEL::Exception)
 473 {
 474   if(part==0)
 475     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::buildSubPart : not empty array must be passed to this method !");
 476   const int *start=part->getConstPointer();
 477   const int *end=start+part->getNbOfElems();
 478   return buildSubPart(start,end);
 479 }
 480
 481 /*!
 482  * Builds a newly created field, that the caller will have the responsability to deal with.
 483  * \n This method makes the assumption that the field \a this is correctly defined when this method is called (\c this->checkCoherency() returns without any exception thrown), **no check of this will be done**.
 484  * \n This method returns a restriction of \a this so that only tuples id specified in [ \a partBg , \a partEnd ) will be contained in returned field.
 485  * \n Parameter [\a partBg, \a partEnd ) specifies \b cell \b ids \b whatever \b the \b spatial \b discretization of \a this
 486  * (\ref ParaMEDMEM::ON_CELLS "ON_CELLS", \ref ParaMEDMEM::ON_NODES "ON_CELLS", \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT", \ref ParaMEDMEM::ON_GAUSS_NE "ON_GAUSS_NE")
 487  *
 488  * If \a this is a field on cell lying on a mesh that have 10 cells. If part contains following cellIds [3,7,6].
 489  * In this case the returned field will lie on mesh having 3 cells and the returned field will contain 3 tuples.
 490  *
 491  *- Tuple#0 of return field will refer to the cell#0 of returned mesh. The cell #0 of returned mesh will be equal to the cell#3 of \c this->getMesh()
 492  *- Tuple#1 of return field will refer to the cell#1 of returned mesh. The cell #1 of returned mesh will be equal to the cell#7 of \c this->getMesh()
 493  *- Tuple#2 of return field will refer to the cell#2 of returned mesh. The cell #2 of returned mesh will be equal to the cell#6 of \c this->getMesh()
 494  *
 495  * If \a this is field on node lying on a mesh that have 10 cells and 11 nodes for example. So \a this is currently contains 11 tuples.
 496  * \n If part contains following cellIds [3,7,6].
 497  * \n If the restriction of mesh to 3 cells leads to a mesh with 6 nodes, the returned field,
 498  * will contain 6 tuples (and same number of components \c this->getArray()->getNumberOfComponents() ) and this field will lie on this restricted mesh.
 499  *
 500  * \param [in] partBg start (included) of input range cell ids to select [ \a partBg, \a partEnd )
 501  * \param [in] partEnd end (not included) of input range cell ids to select [ \a partBg, \a partEnd )
 502  * \return a newly allocated field the caller should deal with.
 503  *
 504  * \throw if there is presence of an invalid cell id in [ \a partBg, \a partEnd ) regarding the number of cells of \c this->getMesh()
 505  *
 506  * \ref cpp_mcfielddouble_subpart1 "Here a C++ example."
 507  *
 508  * \ref py_mcfielddouble_subpart1 "Here a Python example."
 509  * \sa ParaMEDMEM::MEDCouplingFieldDouble::buildSubPart(const DataArrayInt *) const
 510  */
 511 MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPart(const int *partBg, const int *partEnd) const throw(INTERP_KERNEL::Exception)
 512 {
 513   DataArrayInt *arrSelect;
 514   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> m=_type->buildSubMeshData(_mesh,partBg,partEnd,arrSelect);
 515   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arrSelect2(arrSelect);
 516   MEDCouplingFieldDouble *ret=clone(false);//quick shallow copy.
 517   const MEDCouplingFieldDiscretization *disc=getDiscretization();
 518   if(disc)
 519     ret->setDiscretization(MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDiscretization>(disc->clonePart(partBg,partEnd)));
 520   ret->setMesh(m);
 521   std::vector<DataArrayDouble *> arrays;
 522   _time_discr->getArrays(arrays);
 523   std::vector<DataArrayDouble *> arrs;
 524   const int *arrSelBg=arrSelect->begin();
 525   const int *arrSelEnd=arrSelect->end();
 526   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
 527     {
 528       DataArrayDouble *arr=0;
 529       if(*iter)
 530         arr=(*iter)->selectByTupleIdSafe(arrSelBg,arrSelEnd);
 531       arrs.push_back(arr);
 532     }
 533   ret->_time_discr->setArrays(arrs,0);
 534   for(std::vector<DataArrayDouble *>::const_iterator iter=arrs.begin();iter!=arrs.end();iter++)
 535     if(*iter)
 536       (*iter)->decrRef();
 537   return ret;
 538 }
 539
 540 TypeOfTimeDiscretization MEDCouplingFieldDouble::getTimeDiscretization() const
 541 {
 542   return _time_discr->getEnum();
 543 }
 544
 545 MEDCouplingFieldDouble::MEDCouplingFieldDouble(TypeOfField type, TypeOfTimeDiscretization td):MEDCouplingField(type),
 546                                                                                               _time_discr(MEDCouplingTimeDiscretization::New(td))
 547 {
 548 }
 549
 550 /*!
 551  * ** WARINING : This method do not deeply copy neither mesh nor spatial discretization. Only a shallow copy (reference) is done for mesh and spatial discretization ! **
 552  */
 553 MEDCouplingFieldDouble::MEDCouplingFieldDouble(const MEDCouplingFieldTemplate& ft, TypeOfTimeDiscretization td):MEDCouplingField(ft,false),
 554                                                                                                                 _time_discr(MEDCouplingTimeDiscretization::New(td))
 555 {
 556 }
 557
 558 MEDCouplingFieldDouble::MEDCouplingFieldDouble(const MEDCouplingFieldDouble& other, bool deepCopy):MEDCouplingField(other,deepCopy),
 559                                                                                                    _time_discr(other._time_discr->performCpy(deepCopy))
 560 {
 561 }
 562
 563 MEDCouplingFieldDouble::MEDCouplingFieldDouble(NatureOfField n, MEDCouplingTimeDiscretization *td, MEDCouplingFieldDiscretization *type):MEDCouplingField(type,n),_time_discr(td)
 564 {
 565 }
 566
 567 MEDCouplingFieldDouble::~MEDCouplingFieldDouble()
 568 {
 569   delete _time_discr;
 570 }
 571
 572 void MEDCouplingFieldDouble::checkCoherency() const throw(INTERP_KERNEL::Exception)
 573 {
 574   if(!_mesh)
 575     throw INTERP_KERNEL::Exception("Field invalid because no mesh specified !");
 576   _time_discr->checkCoherency();
 577   _type->checkCoherencyBetween(_mesh,getArray());
 578 }
 579
 580 /*!
 581  * Returns the accumulation (the sum) of comId_th component of each tuples of \b default and \b only \b default array.
 582  */
 583 double MEDCouplingFieldDouble::accumulate(int compId) const
 584 {
 585   if(getArray()==0)
 586     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::accumulate : no default array defined !");
 587   return getArray()->accumulate(compId);
 588 }
 589
 590 /*!
 591  * Returns the accumulation (the sum) of all tuples of \b default and \b only default array.
 592  * The res is expected to be of size getNumberOfComponents().
 593  */
 594 void MEDCouplingFieldDouble::accumulate(double *res) const
 595 {
 596   if(getArray()==0)
 597     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::accumulate : no default array defined !");
 598   getArray()->accumulate(res);
 599 }
 600
 601 /*!
 602  * This method returns the max value in \a this. \a this is expected to be a field with exactly \b one component. If not an exception will be thrown.
 603  * To getMaxValue on vector field applyFunc is needed before. This method looks only on all arrays stored in 'this->_time_discr'.
 604  * If no arrays exists, an exception will be thrown.
 605  */
 606 double MEDCouplingFieldDouble::getMaxValue() const throw(INTERP_KERNEL::Exception)
 607 {
 608   std::vector<DataArrayDouble *> arrays;
 609   _time_discr->getArrays(arrays);
 610   double ret=-std::numeric_limits<double>::max();
 611   bool isExistingArr=false;
 612   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
 613     {
 614       if(*iter)
 615         {
 616           isExistingArr=true;
 617           int loc;
 618           ret=std::max(ret,(*iter)->getMaxValue(loc));
 619         }
 620     }
 621   if(!isExistingArr)
 622     throw INTERP_KERNEL::Exception("getMaxValue : No arrays defined !");
 623   return ret;
 624 }
 625
 626 /*!
 627  * This method is an extension of ParaMEDMEM::MEDCouplingFieldDouble::getMaxValue method because the returned
 628  * value is the same but this method also returns to you a tupleIds object which the caller have the responsibility
 629  * to deal with. The main difference is that the returned tupleIds is those corresponding the first set array.
 630  * If you have more than one array set (in LINEAR_TIME instance for example) only the first not null array will be used
 631  * to compute tupleIds.
 632  */
 633 double MEDCouplingFieldDouble::getMaxValue2(DataArrayInt*& tupleIds) const throw(INTERP_KERNEL::Exception)
 634 {
 635   std::vector<DataArrayDouble *> arrays;
 636   _time_discr->getArrays(arrays);
 637   double ret=-std::numeric_limits<double>::max();
 638   bool isExistingArr=false;
 639   tupleIds=0;
 640   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
 641     {
 642       if(*iter)
 643         {
 644           isExistingArr=true;
 645           DataArrayInt *tmp;
 646           ret=std::max(ret,(*iter)->getMaxValue2(tmp));
 647           if(!tupleIds)
 648             tupleIds=tmp;
 649           else
 650             tmp->decrRef();
 651         }
 652     }
 653   if(!isExistingArr)
 654     throw INTERP_KERNEL::Exception("getMaxValue2 : No arrays defined !");
 655   return ret;
 656 }
 657
 658 /*!
 659  * This method returns the min value in \a this. \a this is expected to be a field with exactly \b one component. If not an exception will be thrown.
 660  * To getMinValue on vector field applyFunc is needed before. This method looks only on all arrays stored in 'this->_time_discr'.
 661  * If no arrays exists, an exception will be thrown.
 662  */
 663 double MEDCouplingFieldDouble::getMinValue() const throw(INTERP_KERNEL::Exception)
 664 {
 665   std::vector<DataArrayDouble *> arrays;
 666   _time_discr->getArrays(arrays);
 667   double ret=std::numeric_limits<double>::max();
 668   bool isExistingArr=false;
 669   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
 670     {
 671       if(*iter)
 672         {
 673           isExistingArr=true;
 674           int loc;
 675           ret=std::min(ret,(*iter)->getMinValue(loc));
 676         }
 677     }
 678   if(!isExistingArr)
 679     throw INTERP_KERNEL::Exception("getMinValue : No arrays defined !");
 680   return ret;
 681 }
 682
 683 /*!
 684  * This method is an extension of ParaMEDMEM::MEDCouplingFieldDouble::getMinValue method because the returned
 685  * value is the same but this method also returns to you a tupleIds object which the caller have the responsibility
 686  * to deal with. The main difference is that the returned tupleIds is those corresponding the first set array.
 687  * If you have more than one array set (in LINEAR_TIME instance for example) only the first not null array will be used
 688  * to compute tupleIds.
 689  */
 690 double MEDCouplingFieldDouble::getMinValue2(DataArrayInt*& tupleIds) const throw(INTERP_KERNEL::Exception)
 691 {
 692   std::vector<DataArrayDouble *> arrays;
 693   _time_discr->getArrays(arrays);
 694   double ret=-std::numeric_limits<double>::max();
 695   bool isExistingArr=false;
 696   tupleIds=0;
 697   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
 698     {
 699       if(*iter)
 700         {
 701           isExistingArr=true;
 702           DataArrayInt *tmp;
 703           ret=std::max(ret,(*iter)->getMinValue2(tmp));
 704           if(!tupleIds)
 705             tupleIds=tmp;
 706           else
 707             tmp->decrRef();
 708         }
 709     }
 710   if(!isExistingArr)
 711     throw INTERP_KERNEL::Exception("getMinValue2 : No arrays defined !");
 712   return ret;
 713 }
 714
 715 /*!
 716  * This method returns the average value in \a this. \a this is expected to be a field with exactly \b one component. If not an exception will be thrown.
 717  * To getAverageValue on vector field applyFunc is needed before. This method looks only \b default array \b and \b only \b default.
 718  * If default array does not exist, an exception will be thrown.
 719  */
 720 double MEDCouplingFieldDouble::getAverageValue() const throw(INTERP_KERNEL::Exception)
 721 {
 722   if(getArray()==0)
 723     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getAverageValue : no default array defined !");
 724   return getArray()->getAverageValue();
 725 }
 726
 727 /*!
 728  * This method returns the euclidean norm of \a this.
 729  * \f[
 730  * \sqrt{\sum_{0 \leq i < nbOfEntity}val[i]*val[i]}
 731  * \f]
 732  * If default array does not exist, an exception will be thrown.
 733  */
 734 double MEDCouplingFieldDouble::norm2() const throw(INTERP_KERNEL::Exception)
 735 {
 736   if(getArray()==0)
 737     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::norm2 : no default array defined !");
 738   return getArray()->norm2();
 739 }
 740
 741 /*!
 742  * This method returns the max norm of \a this.
 743  * \f[
 744  * \max_{0 \leq i < nbOfEntity}{abs(val[i])}
 745  * \f]
 746  * If default array does not exist, an exception will be thrown.
 747  */
 748 double MEDCouplingFieldDouble::normMax() const throw(INTERP_KERNEL::Exception)
 749 {
 750   if(getArray()==0)
 751     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::normMax : no default array defined !");
 752   return getArray()->normMax();
 753 }
 754
 755 /*!
 756  * This method returns the average value in \a this weighted by ParaMEDMEM::MEDCouplingField::buildMeasureField.
 757  * \a this is expected to be a field with exactly \b one component. If not an exception will be thrown.
 758  * To getAverageValue on vector field applyFunc is needed before. This method looks only \b default array \b and \b only \b default.
 759  * If default array does not exist, an exception will be thrown.
 760  *
 761  * \param [out] res the location where the result will be stored. \a res is expected to be a location with \c this->getNumberOfComponents() places available.
 762  * \param [in] isWAbs specifies if abs is applied on measure on underlying mesh before performing computation. For a user already sure that all cells of its underlying mesh
 763  *                    are all well oriented this parameter can be set to false to be 'faster'. By default this parameter is true.
 764  */
 765 void MEDCouplingFieldDouble::getWeightedAverageValue(double *res, bool isWAbs) const throw(INTERP_KERNEL::Exception)
 766 {
 767   if(getArray()==0)
 768     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getWeightedAverageValue : no default array defined !");
 769   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> w=buildMeasureField(isWAbs);
 770   double deno=w->getArray()->accumulate(0);
 771   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr=getArray()->deepCpy();
 772   arr->multiplyEqual(w->getArray());
 773   std::transform(arr->begin(),arr->end(),arr->getPointer(),std::bind2nd(std::multiplies<double>(),1./deno));
 774   arr->accumulate(res);
 775 }
 776
 777 /*!
 778  * This method returns the average value in \a this weighted by ParaMEDMEM::MEDCouplingField::buildMeasureField.
 779  * \a this is expected to be a field with exactly \b one component. If not an exception will be thrown.
 780  * To getAverageValue on vector field applyFunc is needed before. This method looks only \b default array \b and \b only \b default.
 781  * If default array does not exist, an exception will be thrown.
 782  *
 783  * \param [in] compId The component id that should be in [0, \c this->getNumberOfComponents() ). If not an INTERP_KERNEL::Exception will be thrown.
 784  * \param [in] isWAbs specifies if abs is applied on measure on underlying mesh before performing computation. For a user already sure that all cells of its underlying mesh
 785  *                    are all well oriented this parameter can be set to false to be 'faster'. By default this parameter is true in C++ not in python (overloading confusion).
 786  */
 787 double MEDCouplingFieldDouble::getWeightedAverageValue(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception)
 788 {
 789   int nbComps=getArray()->getNumberOfComponents();
 790   if(compId<0 || compId>=nbComps)
 791     {
 792       std::ostringstream oss; oss << "MEDCouplingFieldDouble::getWeightedAverageValue : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
 793       throw INTERP_KERNEL::Exception(oss.str().c_str());
 794     }
 795   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
 796   getWeightedAverageValue(res,isWAbs);
 797   return res[compId];
 798 }
 799
 800 /*!
 801  * Returns the normL1 of current field on compId component :
 802  * \f[
 803  * \frac{\sum_{0 \leq i < nbOfEntity}|val[i]*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}
 804  * \f]
 805  * If compId>=nbOfComponent an exception is thrown.
 806  */
 807 double MEDCouplingFieldDouble::normL1(int compId) const throw(INTERP_KERNEL::Exception)
 808 {
 809   if(!_mesh)
 810     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1");
 811   int nbComps=getArray()->getNumberOfComponents();
 812   if(compId<0 || compId>=nbComps)
 813     {
 814       std::ostringstream oss; oss << "MEDCouplingFieldDouble::normL1 : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
 815       throw INTERP_KERNEL::Exception(oss.str().c_str());
 816     }
 817   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
 818   _type->normL1(_mesh,getArray(),res);
 819   return res[compId];
 820 }
 821
 822 /*!
 823  * Returns the normL1 of current field on each components :
 824  * \f[
 825  * \frac{\sum_{0 \leq i < nbOfEntity}|val[i]*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}
 826  * \f]
 827  * The res is expected to be of size getNumberOfComponents().
 828  */
 829 void MEDCouplingFieldDouble::normL1(double *res) const throw(INTERP_KERNEL::Exception)
 830 {
 831   if(!_mesh)
 832     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1");
 833   _type->normL1(_mesh,getArray(),res);
 834 }
 835
 836 /*!
 837  * Returns the normL2 of current field on compId component :
 838  * \f[
 839  * \sqrt{\frac{\sum_{0 \leq i < nbOfEntity}|val[i]^{2}*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}}
 840  * \f]
 841  * If compId>=nbOfComponent an exception is thrown.
 842  */
 843 double MEDCouplingFieldDouble::normL2(int compId) const throw(INTERP_KERNEL::Exception)
 844 {
 845   if(!_mesh)
 846     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL2");
 847   int nbComps=getArray()->getNumberOfComponents();
 848   if(compId<0 || compId>=nbComps)
 849     {
 850       std::ostringstream oss; oss << "MEDCouplingFieldDouble::normL2 : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
 851       throw INTERP_KERNEL::Exception(oss.str().c_str());
 852     }
 853   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
 854   _type->normL2(_mesh,getArray(),res);
 855   return res[compId];
 856 }
 857
 858 /*!
 859  * Returns the normL2 of current field on each components :
 860  * \f[
 861  * \sqrt{\frac{\sum_{0 \leq i < nbOfEntity}|val[i]^{2}*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}}
 862  * \f]
 863  * The res is expected to be of size getNumberOfComponents().
 864  */
 865 void MEDCouplingFieldDouble::normL2(double *res) const throw(INTERP_KERNEL::Exception)
 866 {
 867   if(!_mesh)
 868     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL2");
 869   _type->normL2(_mesh,getArray(),res);
 870 }
 871
 872 /*!
 873  * Returns the accumulation (the sum) of comId_th component of each tuples weigthed by the field
 874  * returns by getWeightingField relative of the _type of field of default array.
 875  * This method is useful to check the conservativity of interpolation method.
 876  */
 877 double MEDCouplingFieldDouble::integral(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception)
 878 {
 879   if(!_mesh)
 880     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform integral");
 881   int nbComps=getArray()->getNumberOfComponents();
 882   if(compId<0 || compId>=nbComps)
 883     {
 884       std::ostringstream oss; oss << "MEDCouplingFieldDouble::integral : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
 885       throw INTERP_KERNEL::Exception(oss.str().c_str());
 886     }
 887   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
 888   _type->integral(_mesh,getArray(),isWAbs,res);
 889   return res[compId];
 890 }
 891
 892 /*!
 893  * Returns the accumulation (the sum) of each tuples weigthed by the field
 894  * returns by getWeightingField relative of the _type of field of default array.
 895  * This method is useful to check the conservativity of interpolation method.
 896  */
 897 void MEDCouplingFieldDouble::integral(bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception)
 898 {
 899   if(!_mesh)
 900     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform integral2");
 901   _type->integral(_mesh,getArray(),isWAbs,res);
 902 }
 903
 904 /*!
 905  * This method is reserved for field lying on structured mesh spatial support. It returns the value of cell localized by (i,j,k)
 906  * If spatial support is not structured mesh an exception will be thrown.
 907  * @param res out array expected to be equal to size getNumberOfComponents()
 908  */
 909 void MEDCouplingFieldDouble::getValueOnPos(int i, int j, int k, double *res) const throw(INTERP_KERNEL::Exception)
 910 {
 911   const DataArrayDouble *arr=_time_discr->getArray();
 912   if(!_mesh)
 913     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOnPos");
 914   _type->getValueOnPos(arr,_mesh,i,j,k,res);
 915 }
 916
 917 /*!
 918  * Returns value of \a this on default time of point 'spaceLoc' using spatial discretization.
 919  * If 'point' is outside the spatial discretization of this an exception will be thrown.
 920  */
 921 void MEDCouplingFieldDouble::getValueOn(const double *spaceLoc, double *res) const throw(INTERP_KERNEL::Exception)
 922 {
 923   const DataArrayDouble *arr=_time_discr->getArray();
 924   if(!_mesh)
 925     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOn");
 926   _type->getValueOn(arr,_mesh,spaceLoc,res);
 927 }
 928
 929 /*!
 930  * Returns a newly allocated array with 'nbOfPoints' tuples and nb of components equal to 'this->getNumberOfComponents()'.
 931  */
 932 DataArrayDouble *MEDCouplingFieldDouble::getValueOnMulti(const double *spaceLoc, int nbOfPoints) const throw(INTERP_KERNEL::Exception)
 933 {
 934   const DataArrayDouble *arr=_time_discr->getArray();
 935   if(!_mesh)
 936     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOnMulti");
 937   return _type->getValueOnMulti(arr,_mesh,spaceLoc,nbOfPoints);
 938 }
 939
 940 /*!
 941  * Returns value of \a this on time 'time' of point 'spaceLoc' using spatial discretization.
 942  * If 'time' is not covered by this->_time_discr an exception will be thrown.
 943  * If 'point' is outside the spatial discretization of this an exception will be thrown.
 944  */
 945 void MEDCouplingFieldDouble::getValueOn(const double *spaceLoc, double time, double *res) const throw(INTERP_KERNEL::Exception)
 946 {
 947   std::vector< const DataArrayDouble *> arrs=_time_discr->getArraysForTime(time);
 948   if(!_mesh)
 949     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOn");
 950   std::vector<double> res2;
 951   for(std::vector< const DataArrayDouble *>::const_iterator iter=arrs.begin();iter!=arrs.end();iter++)
 952     {
 953       int sz=(int)res2.size();
 954       res2.resize(sz+(*iter)->getNumberOfComponents());
 955       _type->getValueOn(*iter,_mesh,spaceLoc,&res2[sz]);
 956     }
 957   _time_discr->getValueForTime(time,res2,res);
 958 }
 959
 960 /*!
 961  * Applies a*x+b on 'compoId'th component of each cell.
 962  */
 963 void MEDCouplingFieldDouble::applyLin(double a, double b, int compoId)
 964 {
 965   _time_discr->applyLin(a,b,compoId);
 966 }
 967
 968 /*!
 969  * This method sets \a this to a uniform scalar field with one component.
 970  * All tuples will have the same value 'value'.
 971  * An exception is thrown if no underlying mesh is defined.
 972  */
 973 MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator=(double value) throw(INTERP_KERNEL::Exception)
 974 {
 975   if(!_mesh)
 976     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::operator= : no mesh defined !");
 977   int nbOfTuple=_type->getNumberOfTuples(_mesh);
 978   _time_discr->setUniformValue(nbOfTuple,1,value);
 979   return *this;
 980 }
 981
 982 /*!
 983  * This method is very similar to this one MEDCouplingMesh::fillFromAnalytic.
 984  * See MEDCouplingMesh::fillFromAnalytic method doc to have more details.
 985  * The main difference is that the field as been started to be constructed here.
 986  * An exception is thrown if no underlying mesh is set before the call of this method.
 987  */
 988 void MEDCouplingFieldDouble::fillFromAnalytic(int nbOfComp, FunctionToEvaluate func) throw(INTERP_KERNEL::Exception)
 989 {
 990   if(!_mesh)
 991     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic : no mesh defined !");
 992   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
 993   _time_discr->fillFromAnalytic(loc,nbOfComp,func);
 994 }
 995
 996 /*!
 997  * This method is very similar to this one MEDCouplingMesh::fillFromAnalytic.
 998  * See MEDCouplingMesh::fillFromAnalytic method doc to have more details.
 999  * The main difference is that the field as been started to be constructed here.
1000  * An exception is thrown if no underlying mesh is set before the call of this method.
1001  */
1002 void MEDCouplingFieldDouble::fillFromAnalytic(int nbOfComp, const char *func) throw(INTERP_KERNEL::Exception)
1003 {
1004   if(!_mesh)
1005     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic : no mesh defined !");
1006   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
1007   _time_discr->fillFromAnalytic(loc,nbOfComp,func);
1008 }
1009
1010 /*!
1011  * This method is very similar to this one MEDCouplingMesh::fillFromAnalytic2.
1012  * See MEDCouplingMesh::fillFromAnalytic method doc to have more details.
1013  * The main difference is that the field as been started to be constructed here.
1014  * An exception is throw if no underlying mesh is set before the call of this method.
1015  */
1016 void MEDCouplingFieldDouble::fillFromAnalytic2(int nbOfComp, const char *func) throw(INTERP_KERNEL::Exception)
1017 {
1018   if(!_mesh)
1019     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic2 : no mesh defined !");
1020   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
1021   _time_discr->fillFromAnalytic2(loc,nbOfComp,func);
1022 }
1023
1024 /*!
1025  * This method is very similar to this one MEDCouplingMesh::fillFromAnalytic3.
1026  * See MEDCouplingMesh::fillFromAnalytic method doc to have more details.
1027  * The main difference is that the field as been started to be constructed here.
1028  * An exception is thrown if no underlying mesh is set before the call of this method.
1029  */
1030 void MEDCouplingFieldDouble::fillFromAnalytic3(int nbOfComp, const std::vector<std::string>& varsOrder, const char *func) throw(INTERP_KERNEL::Exception)
1031 {
1032   if(!_mesh)
1033     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic2 : no mesh defined !");
1034   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
1035   _time_discr->fillFromAnalytic3(loc,nbOfComp,varsOrder,func);
1036 }
1037
1038 /*!
1039  * Applyies the function specified by pointer 'func' on each tuples on all arrays contained in _time_discr.
1040  * If '*func' returns false during one evaluation an exception will be thrown.
1041  */
1042 void MEDCouplingFieldDouble::applyFunc(int nbOfComp, FunctionToEvaluate func)
1043 {
1044   _time_discr->applyFunc(nbOfComp,func);
1045 }
1046
1047 /*!
1048  * This method is a specialization of other overloaded methods. When 'nbOfComp' equals 1 this method is equivalent to
1049  * ParaMEDMEM::MEDCouplingFieldDouble::operator=.
1050  */
1051 void MEDCouplingFieldDouble::applyFunc(int nbOfComp, double val)
1052 {
1053   if(!_mesh)
1054     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::applyFunc : no mesh defined !");
1055   int nbOfTuple=_type->getNumberOfTuples(_mesh);
1056   _time_discr->setUniformValue(nbOfTuple,nbOfComp,val);
1057 }
1058
1059 /*!
1060  * Applyies the function specified by the string repr 'func' on each tuples on all arrays contained in _time_discr.
1061  * If '*func' fails in evaluation during one evaluation an exception will be thrown.
1062  * The field will contain 'nbOfComp' components after the call.
1063  */
1064 void MEDCouplingFieldDouble::applyFunc(int nbOfComp, const char *func) throw(INTERP_KERNEL::Exception)
1065 {
1066   _time_discr->applyFunc(nbOfComp,func);
1067 }
1068
1069 /*!
1070  * This method is equivalent to MEDCouplingFieldDouble::applyFunc, except that here components info are used to determine variables position in 'func'.
1071  * If there is vars detected in 'func' that is not in an info on components an exception will be thrown.
1072  */
1073 void MEDCouplingFieldDouble::applyFunc2(int nbOfComp, const char *func) throw(INTERP_KERNEL::Exception)
1074 {
1075   _time_discr->applyFunc2(nbOfComp,func);
1076 }
1077
1078 /*!
1079  * This method is equivalent to MEDCouplingFieldDouble::applyFunc, except that here 'varsOrder' is used to determine variables position in 'func'.
1080  * If there is vars detected in 'func' that is not in 'varsOrder' an exception will be thrown.
1081  */
1082 void MEDCouplingFieldDouble::applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const char *func) throw(INTERP_KERNEL::Exception)
1083 {
1084   _time_discr->applyFunc3(nbOfComp,varsOrder,func);
1085 }
1086
1087 /*!
1088  * Applyies the function specified by the string repr 'func' on each tuples on all arrays contained in _time_discr.
1089  * If '*func' fails in evaluation during one evaluation an exception will be thrown.
1090  * The field will contain exactly the same number of components after the call.
1091  */
1092 void MEDCouplingFieldDouble::applyFunc(const char *func) throw(INTERP_KERNEL::Exception)
1093 {
1094   _time_discr->applyFunc(func);
1095 }
1096
1097 /*!
1098  * Applyies the function specified by the string repr 'func' on each tuples on all arrays contained in _time_discr.
1099  * The field will contain exactly the same number of components after the call.
1100  * Use is not warranted for the moment !
1101  */
1102 void MEDCouplingFieldDouble::applyFuncFast32(const char *func) throw(INTERP_KERNEL::Exception)
1103 {
1104   _time_discr->applyFuncFast32(func);
1105 }
1106
1107 /*!
1108  * Applyies the function specified by the string repr 'func' on each tuples on all arrays contained in _time_discr.
1109  * The field will contain exactly the same number of components after the call.
1110  * Use is not warranted for the moment !
1111  */
1112 void MEDCouplingFieldDouble::applyFuncFast64(const char *func) throw(INTERP_KERNEL::Exception)
1113 {
1114   _time_discr->applyFuncFast64(func);
1115 }
1116
1117 /*!
1118  * This method makes the assumption that the default array has been set before.
1119  * If not an exception will be sent.
1120  * If default array set, the number of components will be sent.
1121  */
1122 int MEDCouplingFieldDouble::getNumberOfComponents() const throw(INTERP_KERNEL::Exception)
1123 {
1124   if(getArray()==0)
1125     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getNumberOfComponents : No array specified !");
1126   return getArray()->getNumberOfComponents();
1127 }
1128
1129 /*!
1130  * This method makes the assumption that _mesh has be set before the call of this method and description of gauss
1131  * localizations in case of Gauss field. If not an exception will sent.
1132  * \b Contrary to MEDCouplingFieldDouble::getNumberOfComponents and MEDCouplingFieldDouble::getNumberOfValues is
1133  * \b not aware of the presence of the default array.
1134  * \b WARNING \b no coherency check is done here. MEDCouplingFieldDouble::checkCoherency method should be called to check that !
1135  */
1136 int MEDCouplingFieldDouble::getNumberOfTuples() const throw(INTERP_KERNEL::Exception)
1137 {
1138   if(!_mesh)
1139     throw INTERP_KERNEL::Exception("Impossible to retrieve number of tuples because no mesh specified !");
1140   return _type->getNumberOfTuples(_mesh);
1141 }
1142
1143 /*!
1144  * This method makes the assumption that the default array has been set before.
1145  * If not an exception will be sent.
1146  * If default array set, the number of values present in the default array will be sent.
1147  */
1148 int MEDCouplingFieldDouble::getNumberOfValues() const throw(INTERP_KERNEL::Exception)
1149 {
1150   if(getArray()==0)
1151     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getNumberOfValues : No array specified !");
1152   return getArray()->getNbOfElems();
1153 }
1154
1155 void MEDCouplingFieldDouble::updateTime() const
1156 {
1157   MEDCouplingField::updateTime();
1158   updateTimeWith(*_time_discr);
1159 }
1160
1161 std::size_t MEDCouplingFieldDouble::getHeapMemorySize() const
1162 {
1163   std::size_t ret=0;
1164   if(_time_discr)
1165     ret+=_time_discr->getHeapMemorySize();
1166   return MEDCouplingField::getHeapMemorySize()+ret;
1167 }
1168
1169 void MEDCouplingFieldDouble::setNature(NatureOfField nat) throw(INTERP_KERNEL::Exception)
1170 {
1171   MEDCouplingField::setNature(nat);
1172   _type->checkCompatibilityWithNature(nat);
1173 }
1174
1175 /*!
1176  * This method synchronizes time information (time, iteration, order, time unit) regarding the information in \c this->_mesh.
1177  * \throw If no mesh is set in this. Or if \a this is not compatible with time setting (typically NO_TIME)
1178  */
1179 void MEDCouplingFieldDouble::synchronizeTimeWithMesh() throw(INTERP_KERNEL::Exception)
1180 {
1181   if(!_mesh)
1182     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::synchronizeTimeWithMesh : no mesh set in this !");
1183   int it=-1,ordr=-1;
1184   double val=_mesh->getTime(it,ordr);
1185   std::string timeUnit(_mesh->getTimeUnit());
1186   setTime(val,it,ordr);
1187   setTimeUnit(timeUnit.c_str());
1188 }
1189
1190 double MEDCouplingFieldDouble::getIJK(int cellId, int nodeIdInCell, int compoId) const
1191 {
1192   return _type->getIJK(_mesh,getArray(),cellId,nodeIdInCell,compoId);
1193 }
1194
1195 void MEDCouplingFieldDouble::setArray(DataArrayDouble *array)
1196 {
1197   _time_discr->setArray(array,this);
1198 }
1199
1200 void MEDCouplingFieldDouble::setEndArray(DataArrayDouble *array)
1201 {
1202   _time_discr->setEndArray(array,this);
1203 }
1204
1205 void MEDCouplingFieldDouble::setArrays(const std::vector<DataArrayDouble *>& arrs) throw(INTERP_KERNEL::Exception)
1206 {
1207   _time_discr->setArrays(arrs,this);
1208 }
1209
1210 void MEDCouplingFieldDouble::getTinySerializationStrInformation(std::vector<std::string>& tinyInfo) const
1211 {
1212   tinyInfo.clear();
1213   _time_discr->getTinySerializationStrInformation(tinyInfo);
1214   tinyInfo.push_back(_name);
1215   tinyInfo.push_back(_desc);
1216   tinyInfo.push_back(getTimeUnit());
1217 }
1218
1219 /*!
1220  * This method retrieves some critical values to resize and prepare remote instance.
1221  * The first two elements returned in tinyInfo correspond to the parameters to give in constructor.
1222  * @param tinyInfo out parameter resized correctly after the call. The length of this vector is tiny.
1223  */
1224 void MEDCouplingFieldDouble::getTinySerializationIntInformation(std::vector<int>& tinyInfo) const
1225 {
1226   tinyInfo.clear();
1227   tinyInfo.push_back((int)_type->getEnum());
1228   tinyInfo.push_back((int)_time_discr->getEnum());
1229   tinyInfo.push_back((int)_nature);
1230   _time_discr->getTinySerializationIntInformation(tinyInfo);
1231   std::vector<int> tinyInfo2;
1232   _type->getTinySerializationIntInformation(tinyInfo2);
1233   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
1234   tinyInfo.push_back((int)tinyInfo2.size());
1235 }
1236
1237 /*!
1238  * This method retrieves some critical values to resize and prepare remote instance.
1239  * @param tinyInfo out parameter resized correctly after the call. The length of this vector is tiny.
1240  */
1241 void MEDCouplingFieldDouble::getTinySerializationDbleInformation(std::vector<double>& tinyInfo) const
1242 {
1243   tinyInfo.clear();
1244   _time_discr->getTinySerializationDbleInformation(tinyInfo);
1245   std::vector<double> tinyInfo2;
1246   _type->getTinySerializationDbleInformation(tinyInfo2);
1247   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
1248   tinyInfo.push_back((int)tinyInfo2.size());//very bad, lack of time to improve it
1249 }
1250
1251 /*!
1252  * This method has to be called to the new instance filled by CORBA, MPI, File...
1253  * @param tinyInfoI is the value retrieves from distant result of getTinySerializationIntInformation on source instance to be copied.
1254  * @param dataInt out parameter. If not null the pointer is already owned by \a this after the call of this method. In this case no decrRef must be applied.
1255  * @param arrays out parameter is a vector resized to the right size. The pointers in the vector is already owned by \a this after the call of this method.
1256  *               No decrRef must be applied to every instances in returned vector.
1257  */
1258 void MEDCouplingFieldDouble::resizeForUnserialization(const std::vector<int>& tinyInfoI, DataArrayInt *&dataInt, std::vector<DataArrayDouble *>& arrays)
1259 {
1260   dataInt=0;
1261   std::vector<int> tinyInfoITmp(tinyInfoI);
1262   int sz=tinyInfoITmp.back();
1263   tinyInfoITmp.pop_back();
1264   std::vector<int> tinyInfoITmp2(tinyInfoITmp.begin(),tinyInfoITmp.end()-sz);
1265   std::vector<int> tinyInfoI2(tinyInfoITmp2.begin()+3,tinyInfoITmp2.end());
1266   _time_discr->resizeForUnserialization(tinyInfoI2,arrays);
1267   std::vector<int> tinyInfoITmp3(tinyInfoITmp.end()-sz,tinyInfoITmp.end());
1268   _type->resizeForUnserialization(tinyInfoITmp3,dataInt);
1269 }
1270
1271 void MEDCouplingFieldDouble::finishUnserialization(const std::vector<int>& tinyInfoI, const std::vector<double>& tinyInfoD, const std::vector<std::string>& tinyInfoS)
1272 {
1273   std::vector<int> tinyInfoI2(tinyInfoI.begin()+3,tinyInfoI.end());
1274   //
1275   std::vector<double> tmp(tinyInfoD);
1276   int sz=(int)tinyInfoD.back();//very bad, lack of time to improve it
1277   tmp.pop_back();
1278   std::vector<double> tmp1(tmp.begin(),tmp.end()-sz);
1279   std::vector<double> tmp2(tmp.end()-sz,tmp.end());
1280   //
1281   _time_discr->finishUnserialization(tinyInfoI2,tmp1,tinyInfoS);
1282   _nature=(NatureOfField)tinyInfoI[2];
1283   _type->finishUnserialization(tmp2);
1284   int nbOfElemS=(int)tinyInfoS.size();
1285   _name=tinyInfoS[nbOfElemS-3];
1286   _desc=tinyInfoS[nbOfElemS-2];
1287   setTimeUnit(tinyInfoS[nbOfElemS-1].c_str());
1288 }
1289
1290 /*!
1291  * Contrary to MEDCouplingPointSet class the returned arrays are \b not the responsabilities of the caller.
1292  * The values returned must be consulted only in readonly mode.
1293  */
1294 void MEDCouplingFieldDouble::serialize(DataArrayInt *&dataInt, std::vector<DataArrayDouble *>& arrays) const
1295 {
1296   _time_discr->getArrays(arrays);
1297   _type->getSerializationIntArray(dataInt);
1298 }
1299
1300 /*!
1301  * This method tries to to change the mesh support of \a this following the parameter 'levOfCheck' and 'prec'.
1302  * Semantic of 'levOfCheck' is explained in MEDCouplingMesh::checkGeoEquivalWith method. This method is used to perform the job.
1303  * If this->_mesh is not defined or other an exeption will be throw.
1304  */
1305 void MEDCouplingFieldDouble::changeUnderlyingMesh(const MEDCouplingMesh *other, int levOfCheck, double prec) throw(INTERP_KERNEL::Exception)
1306 {
1307   if(_mesh==0 || other==0)
1308     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::changeUnderlyingMesh : is expected to operate on not null meshes !");
1309   DataArrayInt *cellCor=0,*nodeCor=0;
1310   other->checkGeoEquivalWith(_mesh,levOfCheck,prec,cellCor,nodeCor);
1311   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellCor2(cellCor),nodeCor2(nodeCor);
1312   if(cellCor)
1313     renumberCellsWithoutMesh(cellCor->getConstPointer(),false);
1314   if(nodeCor)
1315     renumberNodesWithoutMesh(nodeCor->getConstPointer(),_mesh->getNumberOfNodes());
1316   setMesh(const_cast<MEDCouplingMesh *>(other));
1317 }
1318
1319 /*!
1320  * This method is an extension of MEDCouplingFieldDouble::operator-=. It allows a user to operate a difference of 2 fields (\a this and 'f') even if they do not share same meshes.
1321  * No interpolation will be done here only an analyze of two underlying mesh will be done to see if the meshes are geometrically equivalent. If yes, the eventual renumbering will be done and operator-= applyed after.
1322  * This method requires that 'f' and \a this are coherent (check coherency) and that 'f' and \a this would be coherent for a merge.
1323  * Semantic of 'levOfCheck' is explained in MEDCouplingMesh::checkGeoEquivalWith method.
1324  */
1325 void MEDCouplingFieldDouble::substractInPlaceDM(const MEDCouplingFieldDouble *f, int levOfCheck, double prec) throw(INTERP_KERNEL::Exception)
1326 {
1327   checkCoherency();
1328   f->checkCoherency();
1329   if(!areCompatibleForMerge(f))
1330     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::diffWith : Fields are not compatible ; unable to apply mergeFields on them !");
1331   changeUnderlyingMesh(f->getMesh(),levOfCheck,prec);
1332   operator-=(*f);
1333 }
1334
1335 /*!
1336  * Merge nodes of underlying mesh. In case of some node will be merged the underlying mesh instance will change.
1337  * The first 'eps' stands for geometric approximation. The second 'epsOnVals' is for epsilon on values in case of node merging.
1338  * If 2 nodes distant from less than 'eps' and with value different with more than 'epsOnVals' an exception will be thrown.
1339  */
1340 bool MEDCouplingFieldDouble::mergeNodes(double eps, double epsOnVals) throw(INTERP_KERNEL::Exception)
1341 {
1342   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
1343   if(!meshC)
1344     throw INTERP_KERNEL::Exception("Invalid support mesh to apply mergeNodes on it : must be a MEDCouplingPointSet one !");
1345   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
1346   bool ret;
1347   int ret2;
1348   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->mergeNodes(eps,ret,ret2);
1349   if(!ret)//no nodes have been merged.
1350     return ret;
1351   std::vector<DataArrayDouble *> arrays;
1352   _time_discr->getArrays(arrays);
1353   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1354     if(*iter)
1355       _type->renumberValuesOnNodes(epsOnVals,arr->getConstPointer(),meshC2->getNumberOfNodes(),*iter);
1356   setMesh(meshC2);
1357   return true;
1358 }
1359
1360 /*!
1361  * Merge nodes with (barycenter computation) of underlying mesh. In case of some node will be merged the underlying mesh instance will change.
1362  * The first 'eps' stands for geometric approximation. The second 'epsOnVals' is for epsilon on values in case of node merging.
1363  * If 2 nodes distant from less than 'eps' and with value different with more than 'epsOnVals' an exception will be thrown.
1364  */
1365 bool MEDCouplingFieldDouble::mergeNodes2(double eps, double epsOnVals) throw(INTERP_KERNEL::Exception)
1366 {
1367   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
1368   if(!meshC)
1369     throw INTERP_KERNEL::Exception("Invalid support mesh to apply mergeNodes on it : must be a MEDCouplingPointSet one !");
1370   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
1371   bool ret;
1372   int ret2;
1373   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->mergeNodes2(eps,ret,ret2);
1374   if(!ret)//no nodes have been merged.
1375     return ret;
1376   std::vector<DataArrayDouble *> arrays;
1377   _time_discr->getArrays(arrays);
1378   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1379     if(*iter)
1380       _type->renumberValuesOnNodes(epsOnVals,arr->getConstPointer(),meshC2->getNumberOfNodes(),*iter);
1381   setMesh(meshC2);
1382   return true;
1383 }
1384
1385 /*!
1386  * This method applyies ParaMEDMEM::MEDCouplingPointSet::zipCoords method on 'this->_mesh' that should be set and of type ParaMEDMEM::MEDCouplingPointSet.
1387  * If some nodes have disappeared true is returned.
1388  * 'epsOnVals' stands for epsilon in case of merge of cells. This value is used as tolerance in case the corresponding values differ.
1389  */
1390 bool MEDCouplingFieldDouble::zipCoords(double epsOnVals) throw(INTERP_KERNEL::Exception)
1391 {
1392   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
1393   if(!meshC)
1394     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::zipCoords : Invalid support mesh to apply zipCoords on it : must be a MEDCouplingPointSet one !");
1395   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
1396   int oldNbOfNodes=meshC2->getNumberOfNodes();
1397   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->zipCoordsTraducer();
1398   if(meshC2->getNumberOfNodes()!=oldNbOfNodes)
1399     {
1400       std::vector<DataArrayDouble *> arrays;
1401       _time_discr->getArrays(arrays);
1402       for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1403         if(*iter)
1404           _type->renumberValuesOnNodes(epsOnVals,arr->getConstPointer(),meshC2->getNumberOfNodes(),*iter);
1405       setMesh(meshC2);
1406       return true;
1407     }
1408   return false;
1409 }
1410
1411 /*!
1412  * This method applyies ParaMEDMEM::MEDCouplingUMesh::zipConnectivityTraducer on 'this->_mesh' that should be set and of type ParaMEDMEM::MEDCouplingUMesh.
1413  * The semantic of 'compType' is given in ParaMEDMEM::MEDCouplingUMesh::zipConnectivityTraducer method.
1414  * 'epsOnVals' stands for epsilon in case of merge of cells. This value is used as tolerance in case the corresponding values differ.
1415  */
1416 bool MEDCouplingFieldDouble::zipConnectivity(int compType, double epsOnVals) throw(INTERP_KERNEL::Exception)
1417 {
1418   const MEDCouplingUMesh *meshC=dynamic_cast<const MEDCouplingUMesh *>(_mesh);
1419   if(!meshC)
1420     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::zipCoords : Invalid support mesh to apply zipCoords on it : must be a MEDCouplingPointSet one !");
1421   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> meshC2((MEDCouplingUMesh *)meshC->deepCpy());
1422   int oldNbOfCells=meshC2->getNumberOfCells();
1423   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->zipConnectivityTraducer(compType);
1424   if(meshC2->getNumberOfCells()!=oldNbOfCells)
1425     {
1426       std::vector<DataArrayDouble *> arrays;
1427       _time_discr->getArrays(arrays);
1428       for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1429         if(*iter)
1430           _type->renumberValuesOnCells(epsOnVals,meshC,arr->getConstPointer(),meshC2->getNumberOfCells(),*iter);
1431       setMesh(meshC2);
1432       return true;
1433     }
1434   return false;
1435 }
1436
1437 /*!
1438  * This method calls MEDCouplingUMesh::buildSlice3D method. So this method makes the assumption that underlying mesh exists.
1439  * For the moment, this method is implemented for fields on cells.
1440  *
1441  * \return a newly allocated field double containing the result that the user should deallocate.
1442  */
1443 MEDCouplingFieldDouble *MEDCouplingFieldDouble::extractSlice3D(const double *origin, const double *vec, double eps) const throw(INTERP_KERNEL::Exception)
1444 {
1445   const MEDCouplingMesh *mesh=getMesh();
1446   if(!mesh)
1447     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::extractSlice3D : underlying mesh is null !");
1448   if(getTypeOfField()!=ON_CELLS)
1449     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::extractSlice3D : only implemented for fields on cells !");
1450   const MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> umesh(mesh->buildUnstructured());
1451   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=clone(false);
1452   ret->setMesh(umesh);
1453   DataArrayInt *cellIds=0;
1454   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh2=umesh->buildSlice3D(origin,vec,eps,cellIds);
1455   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellIds2=cellIds;
1456   ret->setMesh(mesh2);
1457   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tupleIds=computeTupleIdsToSelectFromCellIds(cellIds->begin(),cellIds->end());
1458   std::vector<DataArrayDouble *> arrays;
1459   _time_discr->getArrays(arrays);
1460   int i=0;
1461   std::vector<DataArrayDouble *> newArr(arrays.size());
1462   std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > newArr2(arrays.size());
1463   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++,i++)
1464     {
1465       if(*iter)
1466         {
1467           newArr2[i]=(*iter)->selectByTupleIdSafe(cellIds->begin(),cellIds->end());
1468           newArr[i]=newArr2[i];
1469         }
1470     }
1471   ret->setArrays(newArr);
1472   return ret.retn();
1473 }
1474
1475 /*!
1476  * This method applyies ParaMEDMEM::MEDCouplingUMesh::simplexize on 'this->_mesh'.
1477  * The semantic of 'policy' is given in ParaMEDMEM::MEDCouplingUMesh::simplexize method.
1478  */
1479 bool MEDCouplingFieldDouble::simplexize(int policy) throw(INTERP_KERNEL::Exception)
1480 {
1481   int oldNbOfCells=_mesh->getNumberOfCells();
1482   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> meshC2(_mesh->deepCpy());
1483   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->simplexize(policy);
1484   int newNbOfCells=meshC2->getNumberOfCells();
1485   if(oldNbOfCells==newNbOfCells)
1486     return false;
1487   std::vector<DataArrayDouble *> arrays;
1488   _time_discr->getArrays(arrays);
1489   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1490     if(*iter)
1491       _type->renumberValuesOnCellsR(_mesh,arr->getConstPointer(),arr->getNbOfElems(),*iter);
1492   setMesh(meshC2);
1493   return true;
1494 }
1495
1496 MEDCouplingFieldDouble *MEDCouplingFieldDouble::doublyContractedProduct() const throw(INTERP_KERNEL::Exception)
1497 {
1498   MEDCouplingTimeDiscretization *td=_time_discr->doublyContractedProduct();
1499   td->copyTinyAttrFrom(*_time_discr);
1500   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1501   ret->setName("DoublyContractedProduct");
1502   ret->setMesh(getMesh());
1503   return ret;
1504 }
1505
1506 MEDCouplingFieldDouble *MEDCouplingFieldDouble::determinant() const throw(INTERP_KERNEL::Exception)
1507 {
1508   MEDCouplingTimeDiscretization *td=_time_discr->determinant();
1509   td->copyTinyAttrFrom(*_time_discr);
1510   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1511   ret->setName("Determinant");
1512   ret->setMesh(getMesh());
1513   return ret;
1514 }
1515
1516 MEDCouplingFieldDouble *MEDCouplingFieldDouble::eigenValues() const throw(INTERP_KERNEL::Exception)
1517 {
1518   MEDCouplingTimeDiscretization *td=_time_discr->eigenValues();
1519   td->copyTinyAttrFrom(*_time_discr);
1520   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1521   ret->setName("EigenValues");
1522   ret->setMesh(getMesh());
1523   return ret;
1524 }
1525
1526 MEDCouplingFieldDouble *MEDCouplingFieldDouble::eigenVectors() const throw(INTERP_KERNEL::Exception)
1527 {
1528   MEDCouplingTimeDiscretization *td=_time_discr->eigenVectors();
1529   td->copyTinyAttrFrom(*_time_discr);
1530   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1531   ret->setName("EigenVectors");
1532   ret->setMesh(getMesh());
1533   return ret;
1534 }
1535
1536 MEDCouplingFieldDouble *MEDCouplingFieldDouble::inverse() const throw(INTERP_KERNEL::Exception)
1537 {
1538   MEDCouplingTimeDiscretization *td=_time_discr->inverse();
1539   td->copyTinyAttrFrom(*_time_discr);
1540   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1541   ret->setName("Inversion");
1542   ret->setMesh(getMesh());
1543   return ret;
1544 }
1545
1546 MEDCouplingFieldDouble *MEDCouplingFieldDouble::trace() const throw(INTERP_KERNEL::Exception)
1547 {
1548   MEDCouplingTimeDiscretization *td=_time_discr->trace();
1549   td->copyTinyAttrFrom(*_time_discr);
1550   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1551   ret->setName("Trace");
1552   ret->setMesh(getMesh());
1553   return ret;
1554 }
1555
1556 MEDCouplingFieldDouble *MEDCouplingFieldDouble::deviator() const throw(INTERP_KERNEL::Exception)
1557 {
1558   MEDCouplingTimeDiscretization *td=_time_discr->deviator();
1559   td->copyTinyAttrFrom(*_time_discr);
1560   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1561   ret->setName("Trace");
1562   ret->setMesh(getMesh());
1563   return ret;
1564 }
1565
1566 MEDCouplingFieldDouble *MEDCouplingFieldDouble::magnitude() const throw(INTERP_KERNEL::Exception)
1567 {
1568   MEDCouplingTimeDiscretization *td=_time_discr->magnitude();
1569   td->copyTinyAttrFrom(*_time_discr);
1570   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1571   ret->setName("Magnitude");
1572   ret->setMesh(getMesh());
1573   return ret;
1574 }
1575
1576 MEDCouplingFieldDouble *MEDCouplingFieldDouble::maxPerTuple() const throw(INTERP_KERNEL::Exception)
1577 {
1578   MEDCouplingTimeDiscretization *td=_time_discr->maxPerTuple();
1579   td->copyTinyAttrFrom(*_time_discr);
1580   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1581   std::ostringstream oss;
1582   oss << "Max_" << getName();
1583   ret->setName(oss.str().c_str());
1584   ret->setMesh(getMesh());
1585   return ret;
1586 }
1587
1588 void MEDCouplingFieldDouble::changeNbOfComponents(int newNbOfComp, double dftValue) throw(INTERP_KERNEL::Exception)
1589 {
1590   _time_discr->changeNbOfComponents(newNbOfComp,dftValue);
1591 }
1592
1593 MEDCouplingFieldDouble *MEDCouplingFieldDouble::keepSelectedComponents(const std::vector<int>& compoIds) const throw(INTERP_KERNEL::Exception)
1594 {
1595   MEDCouplingTimeDiscretization *td=_time_discr->keepSelectedComponents(compoIds);
1596   td->copyTinyAttrFrom(*_time_discr);
1597   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
1598   ret->setName(getName());
1599   ret->setMesh(getMesh());
1600   return ret;
1601 }
1602
1603 void MEDCouplingFieldDouble::setSelectedComponents(const MEDCouplingFieldDouble *f, const std::vector<int>& compoIds) throw(INTERP_KERNEL::Exception)
1604 {
1605   _time_discr->setSelectedComponents(f->_time_discr,compoIds);
1606 }
1607
1608 void MEDCouplingFieldDouble::sortPerTuple(bool asc) throw(INTERP_KERNEL::Exception)
1609 {
1610   _time_discr->sortPerTuple(asc);
1611 }
1612
1613 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MergeFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1614 {
1615   if(!f1->areCompatibleForMerge(f2))
1616     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MergeFields on them !");
1617   const MEDCouplingMesh *m1=f1->getMesh();
1618   const MEDCouplingMesh *m2=f2->getMesh();
1619   MEDCouplingMesh *m=m1->mergeMyselfWith(m2);
1620   MEDCouplingTimeDiscretization *td=f1->_time_discr->aggregate(f2->_time_discr);
1621   td->copyTinyAttrFrom(*f1->_time_discr);
1622   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1623   ret->setMesh(m);
1624   m->decrRef();
1625   ret->setName(f1->getName());
1626   ret->setDescription(f1->getDescription());
1627   return ret;
1628 }
1629
1630 /*!
1631  * This method returns a newly created field that is the union of all fields in input array 'a'.
1632  * This method expects that 'a' is non empty. If not an exception will be thrown.
1633  * If there is only one field in 'a' a deepCopy (except time information of mesh and field) of the unique field instance in 'a' will be returned.
1634  * Generally speaking the first instance field in 'a' will be used to assign tiny attributes of returned field.
1635  */
1636 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MergeFields(const std::vector<const MEDCouplingFieldDouble *>& a) throw(INTERP_KERNEL::Exception)
1637 {
1638   if(a.size()<1)
1639     throw INTERP_KERNEL::Exception("FieldDouble::MergeFields : size of array must be >= 1 !");
1640   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> > ms(a.size());
1641   std::vector< const MEDCouplingUMesh *> ms2(a.size());
1642   std::vector< const MEDCouplingTimeDiscretization *> tds(a.size());
1643   std::vector<const MEDCouplingFieldDouble *>::const_iterator it=a.begin();
1644   const MEDCouplingFieldDouble *ref=(*it++);
1645   for(;it!=a.end();it++)
1646     if(!ref->areCompatibleForMerge(*it))
1647       throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MergeFields on them !");
1648   for(int i=0;i<(int)a.size();i++)
1649     {
1650       if(!a[i]->getMesh())
1651         throw INTERP_KERNEL::Exception("MergeFields : A field as no underlying mesh !");
1652       ms[i]=a[i]->getMesh()->buildUnstructured();
1653       ms2[i]=ms[i];
1654       tds[i]=a[i]->_time_discr;
1655     }
1656   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=MEDCouplingUMesh::MergeUMeshes(ms2);
1657   m->setName(ms2[0]->getName()); m->setDescription(ms2[0]->getDescription());
1658   MEDCouplingTimeDiscretization *td=tds[0]->aggregate(tds);
1659   td->copyTinyAttrFrom(*(a[0]->_time_discr));
1660   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(a[0]->getNature(),td,a[0]->_type->clone());
1661   ret->setMesh(m);
1662   ret->setName(a[0]->getName());
1663   ret->setDescription(a[0]->getDescription());
1664   return ret;
1665 }
1666
1667 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MeldFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1668 {
1669   if(!f1->areCompatibleForMeld(f2))
1670     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MeldFields on them !");
1671   MEDCouplingTimeDiscretization *td=f1->_time_discr->meld(f2->_time_discr);
1672   td->copyTinyAttrFrom(*f1->_time_discr);
1673   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1674   ret->setMesh(f1->getMesh());
1675   return ret;
1676 }
1677
1678 MEDCouplingFieldDouble *MEDCouplingFieldDouble::DotFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1679 {
1680   if(!f1->areStrictlyCompatible(f2))
1681     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply DotFields on them !");
1682   MEDCouplingTimeDiscretization *td=f1->_time_discr->dot(f2->_time_discr);
1683   td->copyTinyAttrFrom(*f1->_time_discr);
1684   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1685   ret->setMesh(f1->getMesh());
1686   return ret;
1687 }
1688
1689 MEDCouplingFieldDouble *MEDCouplingFieldDouble::CrossProductFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1690 {
1691   if(!f1->areStrictlyCompatible(f2))
1692     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply CrossProductFields on them !");
1693   MEDCouplingTimeDiscretization *td=f1->_time_discr->crossProduct(f2->_time_discr);
1694   td->copyTinyAttrFrom(*f1->_time_discr);
1695   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1696   ret->setMesh(f1->getMesh());
1697   return ret;
1698 }
1699
1700 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MaxFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1701 {
1702   if(!f1->areStrictlyCompatible(f2))
1703     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MaxFields on them !");
1704   MEDCouplingTimeDiscretization *td=f1->_time_discr->max(f2->_time_discr);
1705   td->copyTinyAttrFrom(*f1->_time_discr);
1706   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1707   ret->setMesh(f1->getMesh());
1708   return ret;
1709 }
1710
1711 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MinFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1712 {
1713   if(!f1->areStrictlyCompatible(f2))
1714     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MinFields on them !");
1715   MEDCouplingTimeDiscretization *td=f1->_time_discr->min(f2->_time_discr);
1716   td->copyTinyAttrFrom(*f1->_time_discr);
1717   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1718   ret->setMesh(f1->getMesh());
1719   return ret;
1720 }
1721
1722 MEDCouplingFieldDouble *MEDCouplingFieldDouble::AddFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1723 {
1724   if(!f1->areStrictlyCompatible(f2))
1725     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply AddFields on them !");
1726   MEDCouplingTimeDiscretization *td=f1->_time_discr->add(f2->_time_discr);
1727   td->copyTinyAttrFrom(*f1->_time_discr);
1728   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1729   ret->setMesh(f1->getMesh());
1730   return ret;
1731 }
1732
1733 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator+=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
1734 {
1735   if(!areStrictlyCompatible(&other))
1736     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply += on them !");
1737   _time_discr->addEqual(other._time_discr);
1738   return *this;
1739 }
1740
1741 MEDCouplingFieldDouble *MEDCouplingFieldDouble::SubstractFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1742 {
1743   if(!f1->areStrictlyCompatible(f2))
1744     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply SubstractFields on them !");
1745   MEDCouplingTimeDiscretization *td=f1->_time_discr->substract(f2->_time_discr);
1746   td->copyTinyAttrFrom(*f1->_time_discr);
1747   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1748   ret->setMesh(f1->getMesh());
1749   return ret;
1750 }
1751
1752 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator-=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
1753 {
1754   if(!areStrictlyCompatible(&other))
1755     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply -= on them !");
1756   _time_discr->substractEqual(other._time_discr);
1757   return *this;
1758 }
1759
1760 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MultiplyFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1761 {
1762   if(!f1->areCompatibleForMul(f2))
1763     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MultiplyFields on them !");
1764   MEDCouplingTimeDiscretization *td=f1->_time_discr->multiply(f2->_time_discr);
1765   td->copyTinyAttrFrom(*f1->_time_discr);
1766   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1767   ret->setMesh(f1->getMesh());
1768   return ret;
1769 }
1770
1771 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator*=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
1772 {
1773   if(!areCompatibleForMul(&other))
1774     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply *= on them !");
1775   _time_discr->multiplyEqual(other._time_discr);
1776   return *this;
1777 }
1778
1779 MEDCouplingFieldDouble *MEDCouplingFieldDouble::DivideFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception)
1780 {
1781   if(!f1->areCompatibleForDiv(f2))
1782     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply DivideFields on them !");
1783   MEDCouplingTimeDiscretization *td=f1->_time_discr->divide(f2->_time_discr);
1784   td->copyTinyAttrFrom(*f1->_time_discr);
1785   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
1786   ret->setMesh(f1->getMesh());
1787   return ret;
1788 }
1789
1790 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator/=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
1791 {
1792   if(!areCompatibleForDiv(&other))
1793     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply /= on them !");
1794   _time_discr->divideEqual(other._time_discr);
1795   return *this;
1796 }
1797
1798 /*!
1799  * This method writes the field series 'fs' in the VTK file 'fileName'.
1800  * If 'fs' is empty no file is written. If fields lies on more than one mesh an exception will be thrown and no file will be written too.
1801  * If the single mesh is empty an exception will be thrown.
1802  * Finally there is a field in 'fs' with no name an exception will be thrown too.
1803  */
1804 void MEDCouplingFieldDouble::WriteVTK(const char *fileName, const std::vector<const MEDCouplingFieldDouble *>& fs) throw(INTERP_KERNEL::Exception)
1805 {
1806   if(fs.empty())
1807     return;
1808   std::size_t nfs=fs.size();
1809   const MEDCouplingMesh *m=fs[0]->getMesh();
1810   for(std::size_t i=1;i<nfs;i++)
1811     if(fs[i]->getMesh()!=m)
1812       throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::WriteVTK : Fields are not lying on a same mesh ! Expected by VTK ! MEDCouplingFieldDouble::setMesh or MEDCouplingFieldDouble::changeUnderlyingMesh can help to that.");
1813   if(!m)
1814     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::WriteVTK : Fields are lying on a same mesh but it is empty !");
1815   std::ostringstream coss,noss;
1816   for(std::size_t i=0;i<nfs;i++)
1817     {
1818       const MEDCouplingFieldDouble *cur=fs[i];
1819       std::string name(cur->getName());
1820       if(name.empty())
1821         {
1822           std::ostringstream oss; oss << "MEDCouplingFieldDouble::WriteVTK : Field in pos #" << i << " has no name !";
1823           throw INTERP_KERNEL::Exception(oss.str().c_str());
1824         }
1825       TypeOfField typ=cur->getTypeOfField();
1826       if(typ==ON_CELLS)
1827         cur->getArray()->writeVTK(coss,8,cur->getName());
1828       else if(typ==ON_NODES)
1829         cur->getArray()->writeVTK(noss,8,cur->getName());
1830     }
1831   m->writeVTKAdvanced(fileName,coss.str(),noss.str());
1832 }