Salome HOME
Merge branch 'origin/agy/V731rc2'
[tools/medcoupling.git] / 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 /*!
40  * Creates a new MEDCouplingFieldDouble, of given spatial type and time discretization.
41  * For more info, see \ref MEDCouplingFirstSteps3.
42  * \param [in] type - the type of spatial discretization of the created field, one of
43  *        (\ref ParaMEDMEM::ON_CELLS "ON_CELLS", 
44  *         \ref ParaMEDMEM::ON_NODES "ON_NODES",
45  *         \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT", 
46  *         \ref ParaMEDMEM::ON_GAUSS_NE "ON_GAUSS_NE",
47  *         \ref ParaMEDMEM::ON_NODES_KR "ON_NODES_KR").
48  * \param [in] td - the type of time discretization of the created field, one of
49  *        (\ref ParaMEDMEM::NO_TIME "NO_TIME", 
50  *         \ref ParaMEDMEM::ONE_TIME "ONE_TIME", 
51  *         \ref ParaMEDMEM::LINEAR_TIME "LINEAR_TIME", 
52  *         \ref ParaMEDMEM::CONST_ON_TIME_INTERVAL "CONST_ON_TIME_INTERVAL").
53  * \return MEDCouplingFieldDouble* - a new instance of MEDCouplingFieldDouble. The
54  *         caller is to delete this field using decrRef() as it is no more needed. 
55  */
56 MEDCouplingFieldDouble* MEDCouplingFieldDouble::New(TypeOfField type, TypeOfTimeDiscretization td)
57 {
58   return new MEDCouplingFieldDouble(type,td);
59 }
60
61 /*!
62  * Creates a new MEDCouplingFieldDouble, of a given time discretization and with a
63  * spatial type and supporting mesh copied from a given 
64  * \ref MEDCouplingFieldTemplatesPage "field template".
65  * For more info, see \ref MEDCouplingFirstSteps3.
66  * \warning This method does not deeply copy neither the mesh nor the spatial
67  * discretization. Only a shallow copy (reference) is done for the mesh and the spatial
68  * discretization!
69  * \param [in] ft - the \ref MEDCouplingFieldTemplatesPage "field template" defining
70  *        the spatial discretization and the supporting mesh.
71  * \param [in] td - the type of time discretization of the created field, one of
72  *        (\ref ParaMEDMEM::NO_TIME "NO_TIME", 
73  *         \ref ParaMEDMEM::ONE_TIME "ONE_TIME", 
74  *         \ref ParaMEDMEM::LINEAR_TIME "LINEAR_TIME", 
75  *         \ref ParaMEDMEM::CONST_ON_TIME_INTERVAL "CONST_ON_TIME_INTERVAL").
76  * \return MEDCouplingFieldDouble* - a new instance of MEDCouplingFieldDouble. The
77  *         caller is to delete this field using decrRef() as it is no more needed. 
78  */
79 MEDCouplingFieldDouble *MEDCouplingFieldDouble::New(const MEDCouplingFieldTemplate& ft, TypeOfTimeDiscretization td)
80 {
81   return new MEDCouplingFieldDouble(ft,td);
82 }
83
84 /*!
85  * Sets a time \a unit of \a this field. For more info, see \ref MEDCouplingFirstSteps3.
86  * \param [in] unit \a unit (string) in which time is measured.
87  */
88 void MEDCouplingFieldDouble::setTimeUnit(const std::string& unit)
89 {
90   _time_discr->setTimeUnit(unit);
91 }
92
93 /*!
94  * Returns a time unit of \a this field.
95  * \return a string describing units in which time is measured.
96  */
97 std::string MEDCouplingFieldDouble::getTimeUnit() const
98 {
99   return _time_discr->getTimeUnit();
100 }
101
102 /*!
103  * This method if possible the time information (time unit, time iteration, time unit and time value) with its support
104  * that is to say its mesh.
105  * 
106  * \throw  If \c this->_mesh is null an exception will be thrown. An exception will also be throw if the spatial discretization is
107  *         NO_TIME.
108  */
109 void MEDCouplingFieldDouble::synchronizeTimeWithSupport()
110 {
111   _time_discr->synchronizeTimeWith(_mesh);
112 }
113
114 /*!
115  * Returns a new MEDCouplingFieldDouble which is a copy of \a this one. The data
116  * of \a this field is copied either deep or shallow depending on \a recDeepCpy
117  * parameter. But the underlying mesh is always shallow copied.
118  * Data that can be copied either deeply or shallow are:
119  * - \ref MEDCouplingTemporalDisc "temporal discretization" data that holds array(s)
120  * of field values,
121  * - \ref MEDCouplingSpatialDisc "a spatial discretization".
122  * 
123  * \c clone(false) is rather dedicated for advanced users that want to limit the amount
124  * of memory. It allows the user to perform methods like operator+(), operator*()
125  * etc. with \a this and the returned field. If the user wants to duplicate deeply the
126  * underlying mesh he should call cloneWithMesh() method or deepCpy() instead. 
127  * \warning The underlying \b mesh of the returned field is **always the same**
128  *         (pointer) as \a this one **whatever the value** of \a recDeepCpy parameter.
129  *  \param [in] recDeepCpy - if \c true, the copy of the underlying data arrays is
130  *         deep, else all data arrays of \a this field are shared by the new field.
131  *  \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble. The
132  *         caller is to delete this field using decrRef() as it is no more needed.
133  * \sa cloneWithMesh()
134  */
135 MEDCouplingFieldDouble *MEDCouplingFieldDouble::clone(bool recDeepCpy) const
136 {
137   return new MEDCouplingFieldDouble(*this,recDeepCpy);
138 }
139
140 /*!
141  * Returns a new MEDCouplingFieldDouble which is a copy of \a this one. The data
142  * of \a this field is copied either deep or shallow depending on \a recDeepCpy
143  * parameter. But the underlying mesh is always deep copied.
144  * Data that can be copied either deeply or shallow are:
145  * - \ref MEDCouplingTemporalDisc "temporal discretization" data that holds array(s)
146  * of field values,
147  * - \ref MEDCouplingSpatialDisc "a spatial discretization".
148  * 
149  * This method behaves exactly like clone() except that here the underlying **mesh is
150  * always deeply duplicated**, whatever the value \a recDeepCpy parameter.
151  * The result of \c cloneWithMesh(true) is exactly the same as that of deepCpy().
152  * So the resulting field can not be used together with \a this one in the methods
153  * like operator+(), operator*() etc. To avoid deep copying the underlying mesh,
154  * the user can call clone().
155  *  \param [in] recDeepCpy - if \c true, the copy of the underlying data arrays is
156  *         deep, else all data arrays of \a this field are shared by the new field.
157  *  \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble. The
158  *         caller is to delete this field using decrRef() as it is no more needed.
159  * \sa clone()
160  */
161 MEDCouplingFieldDouble *MEDCouplingFieldDouble::cloneWithMesh(bool recDeepCpy) const
162 {
163   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=clone(recDeepCpy);
164   if(_mesh)
165     {
166       MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> mCpy=_mesh->deepCpy();
167       ret->setMesh(mCpy);
168     }
169   return ret.retn();
170 }
171
172 /*!
173  * Returns a new MEDCouplingFieldDouble which is a deep copy of \a this one **including
174  * the mesh**.
175  * The result of this method is exactly the same as that of \c cloneWithMesh(true).
176  * So the resulting field can not be used together with \a this one in the methods
177  * like operator+(), operator*() etc. To avoid deep copying the underlying mesh,
178  * the user can call clone().
179  *  \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble. The
180  *         caller is to delete this field using decrRef() as it is no more needed.
181  * \sa cloneWithMesh()
182  */
183 MEDCouplingFieldDouble *MEDCouplingFieldDouble::deepCpy() const
184 {
185   return cloneWithMesh(true);
186 }
187
188 /*!
189  * Creates a new MEDCouplingFieldDouble of given
190  * \ref MEDCouplingTemporalDisc "temporal discretization". The result field either
191  * shares the data array(s) with \a this field, or holds a deep copy of it, depending on
192  * \a deepCopy parameter. But the underlying \b mesh is always **shallow copied**.
193  * \param [in] td - the type of time discretization of the created field, one of
194  *        (\ref ParaMEDMEM::NO_TIME "NO_TIME", 
195  *         \ref ParaMEDMEM::ONE_TIME "ONE_TIME", 
196  *         \ref ParaMEDMEM::LINEAR_TIME "LINEAR_TIME", 
197  *         \ref ParaMEDMEM::CONST_ON_TIME_INTERVAL "CONST_ON_TIME_INTERVAL").
198  * \param [in] deepCopy - if \c true, the copy of the underlying data arrays is
199  *         deep, else all data arrays of \a this field are shared by the new field.
200  * \return MEDCouplingFieldDouble* - a new instance of MEDCouplingFieldDouble. The
201  *         caller is to delete this field using decrRef() as it is no more needed. 
202  * 
203  * \if ENABLE_EXAMPLES
204  * \ref cpp_mcfielddouble_buildNewTimeReprFromThis "Here is a C++ example."<br>
205  * \ref py_mcfielddouble_buildNewTimeReprFromThis "Here is a Python example."
206  * \endif
207  * \sa clone()
208  */
209 MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildNewTimeReprFromThis(TypeOfTimeDiscretization td, bool deepCopy) const
210 {
211   MEDCouplingTimeDiscretization *tdo=_time_discr->buildNewTimeReprFromThis(td,deepCopy);
212   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDiscretization> disc;
213   if(_type)
214     disc=_type->clone();
215   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),tdo,disc.retn());
216   ret->setMesh(getMesh());
217   ret->setName(getName());
218   ret->setDescription(getDescription());
219   return ret.retn();
220 }
221
222 /*!
223  * This method converts a field on nodes (\a this) to a cell field (returned field). The convertion is a \b non \b conservative remapping !
224  * This method is useful only for users that need a fast convertion from node to cell spatial discretization. The algorithm applied is simply to attach
225  * to each cell the average of values on nodes constituting this cell.
226  *
227  * \return MEDCouplingFieldDouble* - a new instance of MEDCouplingFieldDouble. The
228  *         caller is to delete this field using decrRef() as it is no more needed. The returned field will share the same mesh object object than those in \a this.
229  * \throw If \a this spatial discretization is empty or not ON_NODES.
230  * \throw If \a this is not coherent (see MEDCouplingFieldDouble::checkCoherency).
231  * 
232  * \warning This method is a \b non \b conservative method of remapping from node spatial discretization to cell spatial discretization.
233  * If a conservative method of interpolation is required ParaMEDMEM::MEDCouplingRemapper class should be used instead with "P1P0" method.
234  */
235 MEDCouplingFieldDouble *MEDCouplingFieldDouble::nodeToCellDiscretization() const
236 {
237   checkCoherency();
238   TypeOfField tf(getTypeOfField());
239   if(tf!=ON_NODES)
240     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::nodeToCellDiscretization : this field is expected to be on ON_NODES !");
241   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret(clone(false));
242   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDiscretizationP0> nsp(new MEDCouplingFieldDiscretizationP0);
243   ret->setDiscretization(nsp);
244   const MEDCouplingMesh *m(getMesh());//m is non empty thanks to checkCoherency call
245   int nbCells(m->getNumberOfCells());
246   std::vector<DataArrayDouble *> arrs(getArrays());
247   std::size_t sz(arrs.size());
248   std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > outArrsSafe(sz); std::vector<DataArrayDouble *> outArrs(sz);
249   for(std::size_t j=0;j<sz;j++)
250     {
251       int nbCompo(arrs[j]->getNumberOfComponents());
252       outArrsSafe[j]=DataArrayDouble::New(); outArrsSafe[j]->alloc(nbCells,nbCompo);
253       outArrsSafe[j]->copyStringInfoFrom(*arrs[j]);
254       outArrs[j]=outArrsSafe[j];
255       double *pt(outArrsSafe[j]->getPointer());
256       const double *srcPt(arrs[j]->begin());
257       for(int i=0;i<nbCells;i++,pt+=nbCompo)
258         {
259           std::vector<int> nodeIds;
260           m->getNodeIdsOfCell(i,nodeIds);
261           std::fill(pt,pt+nbCompo,0.);
262           std::size_t nbNodesInCell(nodeIds.size());
263           for(std::size_t k=0;k<nbNodesInCell;k++)
264             std::transform(srcPt+nodeIds[k]*nbCompo,srcPt+(nodeIds[k]+1)*nbCompo,pt,pt,std::plus<double>());
265           if(nbNodesInCell!=0)
266             std::transform(pt,pt+nbCompo,pt,std::bind2nd(std::multiplies<double>(),1./((double)nbNodesInCell)));
267           else
268             {
269               std::ostringstream oss; oss << "MEDCouplingFieldDouble::nodeToCellDiscretization : Cell id #" << i << " has been detected to have no nodes !";
270               throw INTERP_KERNEL::Exception(oss.str().c_str());
271             }
272         }
273     }
274   ret->setArrays(outArrs);
275   return ret.retn();
276 }
277
278 /*!
279  * This method converts a field on cell (\a this) to a node field (returned field). The convertion is a \b non \b conservative remapping !
280  * This method is useful only for users that need a fast convertion from cell to node spatial discretization. The algorithm applied is simply to attach
281  * to each node the average of values on cell sharing this node. If \a this lies on a mesh having orphan nodes the values applied on them will be NaN (division by 0.).
282  *
283  * \return MEDCouplingFieldDouble* - a new instance of MEDCouplingFieldDouble. The
284  *         caller is to delete this field using decrRef() as it is no more needed. The returned field will share the same mesh object object than those in \a this.
285  * \throw If \a this spatial discretization is empty or not ON_CELLS.
286  * \throw If \a this is not coherent (see MEDCouplingFieldDouble::checkCoherency).
287  *
288  * \warning This method is a \b non \b conservative method of remapping from cell spatial discretization to node spatial discretization.
289  * If a conservative method of interpolation is required ParaMEDMEM::MEDCouplingRemapper class should be used instead with "P0P1" method.
290  */
291 MEDCouplingFieldDouble *MEDCouplingFieldDouble::cellToNodeDiscretization() const
292 {
293   checkCoherency();
294   TypeOfField tf(getTypeOfField());
295   if(tf!=ON_CELLS)
296     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::cellToNodeDiscretization : this field is expected to be on ON_CELLS !");
297   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret(clone(false));
298   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDiscretizationP1> nsp(new MEDCouplingFieldDiscretizationP1);
299   ret->setDiscretization(nsp);
300   const MEDCouplingMesh *m(getMesh());//m is non empty thanks to checkCoherency call
301   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> rn(DataArrayInt::New()),rni(DataArrayInt::New());
302   m->getReverseNodalConnectivity(rn,rni);
303   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> rni2(rni->deltaShiftIndex());
304   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> rni3(rni2->convertToDblArr()); rni2=0;
305   std::vector<DataArrayDouble *> arrs(getArrays());
306   std::size_t sz(arrs.size());
307   std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > outArrsSafe(sz); std::vector<DataArrayDouble *> outArrs(sz);
308   for(std::size_t j=0;j<sz;j++)
309     {
310       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> tmp(arrs[j]->selectByTupleIdSafe(rn->begin(),rn->end()));
311       outArrsSafe[j]=(tmp->accumulatePerChunck(rni->begin(),rni->end())); tmp=0;
312       outArrsSafe[j]->divideEqual(rni3);
313       outArrsSafe[j]->copyStringInfoFrom(*arrs[j]);
314       outArrs[j]=outArrsSafe[j];
315     }
316   ret->setArrays(outArrs);
317   return ret.retn();
318 }
319
320 /*!
321  * Copies tiny info (component names, name and description) from an \a other field to
322  * \a this one.
323  * \warning The underlying mesh is not renamed (for safety reason).
324  *  \param [in] other - the field to copy the tiny info from.
325  *  \throw If \a this->getNumberOfComponents() != \a other->getNumberOfComponents()
326  */
327 void MEDCouplingFieldDouble::copyTinyStringsFrom(const MEDCouplingField *other)
328 {
329   MEDCouplingField::copyTinyStringsFrom(other);
330   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
331   if(otherC)
332     {
333       _time_discr->copyTinyStringsFrom(*otherC->_time_discr);
334     }
335 }
336
337 /*!
338  * Copies only times, order and iteration from an \a other field to
339  * \a this one. The underlying mesh is not impacted by this method.
340  * Arrays are not impacted neither.
341  *  \param [in] other - the field to tiny attributes from.
342  *  \throw If \a this->getNumberOfComponents() != \a other->getNumberOfComponents()
343  */
344 void MEDCouplingFieldDouble::copyTinyAttrFrom(const MEDCouplingFieldDouble *other)
345 {
346   if(other)
347     {
348       _time_discr->copyTinyAttrFrom(*other->_time_discr);
349     }
350   
351 }
352
353 void MEDCouplingFieldDouble::copyAllTinyAttrFrom(const MEDCouplingFieldDouble *other)
354 {
355   copyTinyStringsFrom(other);
356   copyTinyAttrFrom(other);
357 }
358
359 /*!
360  * Returns a string describing \a this field. This string is outputted by \c print
361  * Python command. The string includes info on
362  * - name,
363  * - description,
364  * - \ref MEDCouplingSpatialDisc "spatial discretization",
365  * - \ref MEDCouplingTemporalDisc "time discretization",
366  * - \ref NatureOfField,
367  * - components,
368  * - mesh.
369  *
370  *  \return std::string - the string describing \a this field.
371  */
372 std::string MEDCouplingFieldDouble::simpleRepr() const
373 {
374   std::ostringstream ret;
375   ret << "FieldDouble with name : \"" << getName() << "\"\n";
376   ret << "Description of field is : \"" << getDescription() << "\"\n";
377   if(_type)
378     { ret << "FieldDouble space discretization is : " << _type->getStringRepr() << "\n"; }
379   else
380     { ret << "FieldDouble has no spatial discretization !\n"; }
381   if(_time_discr)
382     { ret << "FieldDouble time discretization is : " << _time_discr->getStringRepr() << "\n"; }
383   else
384     { ret << "FieldDouble has no time discretization !\n"; }
385   ret << "FieldDouble nature of field is : \"" << MEDCouplingNatureOfField::GetReprNoThrow(_nature) << "\"\n";
386   if(getArray())
387     {
388       if(getArray()->isAllocated())
389         {
390           int nbOfCompo=getArray()->getNumberOfComponents();
391           ret << "FieldDouble default array has " << nbOfCompo << " components and " << getArray()->getNumberOfTuples() << " tuples.\n";
392           ret << "FieldDouble default array has following info on components : ";
393           for(int i=0;i<nbOfCompo;i++)
394             ret << "\"" << getArray()->getInfoOnComponent(i) << "\" ";
395           ret << "\n";
396         }
397       else
398         {
399           ret << "Array set but not allocated !\n";
400         }
401     }
402   if(_mesh)
403     ret << "Mesh support information :\n__________________________\n" << _mesh->simpleRepr();
404   else
405     ret << "Mesh support information : No mesh set !\n";
406   return ret.str();
407 }
408
409 /*!
410  * Returns a string describing \a this field. The string includes info on
411  * - name,
412  * - description,
413  * - \ref MEDCouplingSpatialDisc "spatial discretization",
414  * - \ref MEDCouplingTemporalDisc "time discretization",
415  * - components,
416  * - mesh,
417  * - contents of data arrays.
418  *
419  *  \return std::string - the string describing \a this field.
420  */
421 std::string MEDCouplingFieldDouble::advancedRepr() const
422 {
423   std::ostringstream ret;
424   ret << "FieldDouble with name : \"" << getName() << "\"\n";
425   ret << "Description of field is : \"" << getDescription() << "\"\n";
426   if(_type)
427     { ret << "FieldDouble space discretization is : " << _type->getStringRepr() << "\n"; }
428   else
429     { ret << "FieldDouble has no space discretization set !\n"; }
430   if(_time_discr)
431     { ret << "FieldDouble time discretization is : " << _time_discr->getStringRepr() << "\n"; }
432   else
433     { ret << "FieldDouble has no time discretization set !\n"; }
434   if(getArray())
435     ret << "FieldDouble default array has " << getArray()->getNumberOfComponents() << " components and " << getArray()->getNumberOfTuples() << " tuples.\n";
436   if(_mesh)
437     ret << "Mesh support information :\n__________________________\n" << _mesh->advancedRepr();
438   else
439     ret << "Mesh support information : No mesh set !\n";
440   std::vector<DataArrayDouble *> arrays;
441   _time_discr->getArrays(arrays);
442   int arrayId=0;
443   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++,arrayId++)
444     {
445       ret << "Array #" << arrayId << " :\n__________\n";
446       if(*iter)
447         (*iter)->reprWithoutNameStream(ret);
448       else
449         ret << "Array empty !";
450       ret << "\n";
451     }
452   return ret.str();
453 }
454
455 void MEDCouplingFieldDouble::writeVTK(const std::string& fileName, bool isBinary) const
456 {
457   std::vector<const MEDCouplingFieldDouble *> fs(1,this);
458   MEDCouplingFieldDouble::WriteVTK(fileName,fs,isBinary);
459 }
460
461 bool MEDCouplingFieldDouble::isEqualIfNotWhy(const MEDCouplingField *other, double meshPrec, double valsPrec, std::string& reason) const
462 {
463   if(!other)
464     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::isEqualIfNotWhy : other instance is NULL !");
465   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
466   if(!otherC)
467     {
468       reason="field given in input is not castable in MEDCouplingFieldDouble !";
469       return false;
470     }
471   if(!MEDCouplingField::isEqualIfNotWhy(other,meshPrec,valsPrec,reason))
472     return false;
473   if(!_time_discr->isEqualIfNotWhy(otherC->_time_discr,valsPrec,reason))
474     {
475       reason.insert(0,"In FieldDouble time discretizations differ :");
476       return false;
477     }
478   return true;
479 }
480
481 /*!
482  * Checks equality of \a this and \a other field. Only numeric data is considered,
483  * i.e. names, description etc are not compared.
484  *  \param [in] other - the field to compare with.
485  *  \param [in] meshPrec - a precision used to compare node coordinates of meshes.
486  *  \param [in] valsPrec - a precision used to compare data arrays of the two fields.
487  *  \return bool - \c true if the two fields are equal, \c false else.
488  *  \throw If \a other == NULL.
489  *  \throw If the spatial discretization of \a this field is NULL.
490  */
491 bool MEDCouplingFieldDouble::isEqualWithoutConsideringStr(const MEDCouplingField *other, double meshPrec, double valsPrec) const
492 {
493   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
494   if(!otherC)
495     return false;
496   if(!MEDCouplingField::isEqualWithoutConsideringStr(other,meshPrec,valsPrec))
497     return false;
498   if(!_time_discr->isEqualWithoutConsideringStr(otherC->_time_discr,valsPrec))
499     return false;
500   return true;
501 }
502
503 /*!
504  * This method states if \a this and 'other' are compatibles each other before performing any treatment.
505  * This method is good for methods like : mergeFields.
506  * This method is not very demanding compared to areStrictlyCompatible that is better for operation on fields.
507  */
508 bool MEDCouplingFieldDouble::areCompatibleForMerge(const MEDCouplingField *other) const
509 {
510   if(!MEDCouplingField::areCompatibleForMerge(other))
511     return false;
512   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
513   if(!otherC)
514     return false;
515   if(!_time_discr->areCompatible(otherC->_time_discr))
516     return false;
517   return true;
518 }
519
520 /*!
521  * This method is more strict than MEDCouplingField::areCompatibleForMerge method.
522  * This method is used for operation on fields to operate a first check before attempting operation.
523  */
524 bool MEDCouplingFieldDouble::areStrictlyCompatible(const MEDCouplingField *other) const
525 {
526   std::string tmp;
527   if(!MEDCouplingField::areStrictlyCompatible(other))
528     return false;
529   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
530   if(!otherC)
531     return false;
532   if(!_time_discr->areStrictlyCompatible(otherC->_time_discr,tmp))
533     return false;
534   return true;
535 }
536
537 /*!
538  * Method with same principle than MEDCouplingFieldDouble::areStrictlyCompatible method except that
539  * number of components between \a this and 'other' can be different here (for operator*).
540  */
541 bool MEDCouplingFieldDouble::areCompatibleForMul(const MEDCouplingField *other) const
542 {
543   if(!MEDCouplingField::areStrictlyCompatible(other))
544     return false;
545   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
546   if(!otherC)
547     return false;
548   if(!_time_discr->areStrictlyCompatibleForMul(otherC->_time_discr))
549     return false;
550   return true;
551 }
552
553 /*!
554  * Method with same principle than MEDCouplingFieldDouble::areStrictlyCompatible method except that
555  * number of components between \a this and 'other' can be different here (for operator/).
556  */
557 bool MEDCouplingFieldDouble::areCompatibleForDiv(const MEDCouplingField *other) const
558 {
559   if(!MEDCouplingField::areStrictlyCompatible(other))
560     return false;
561   const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
562   if(!otherC)
563     return false;
564   if(!_time_discr->areStrictlyCompatibleForDiv(otherC->_time_discr))
565     return false;
566   return true;
567 }
568
569 /*!
570  * This method is invocated before any attempt of melding. This method is very close to areStrictlyCompatible,
571  * except that \a this and other can have different number of components.
572  */
573 bool MEDCouplingFieldDouble::areCompatibleForMeld(const MEDCouplingFieldDouble *other) const
574 {
575   if(!MEDCouplingField::areStrictlyCompatible(other))
576     return false;
577   if(!_time_discr->areCompatibleForMeld(other->_time_discr))
578     return false;
579   return true;
580 }
581
582 /*!
583  * Permutes values of \a this field according to a given permutation array for cells
584  * renumbering. The underlying mesh is deeply copied and its cells are also permuted. 
585  * The number of cells remains the same; for that the permutation array \a old2NewBg
586  * should not contain equal ids.
587  * ** Warning, this method modifies the mesh aggreagated by \a this (by performing a deep copy ) **.
588  *
589  *  \param [in] old2NewBg - the permutation array in "Old to New" mode. Its length is
590  *         to be equal to \a this->getMesh()->getNumberOfCells().
591  *  \param [in] check - if \c true, \a old2NewBg is transformed to a new permutation
592  *         array, so that its maximal cell id to correspond to (be less than) the number
593  *         of cells in mesh. This new array is then used for the renumbering. If \a 
594  *         check == \c false, \a old2NewBg is used as is, that is less secure as validity 
595  *         of ids in \a old2NewBg is not checked.
596  *  \throw If the mesh is not set.
597  *  \throw If the spatial discretization of \a this field is NULL.
598  *  \throw If \a check == \c true and \a old2NewBg contains equal ids.
599  *  \throw If mesh nature does not allow renumbering (e.g. structured mesh).
600  * 
601  *  \if ENABLE_EXAMPLES
602  *  \ref cpp_mcfielddouble_renumberCells "Here is a C++ example".<br>
603  *  \ref  py_mcfielddouble_renumberCells "Here is a Python example".
604  *  \endif
605  */
606 void MEDCouplingFieldDouble::renumberCells(const int *old2NewBg, bool check)
607 {
608   renumberCellsWithoutMesh(old2NewBg,check);
609   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> m=_mesh->deepCpy();
610   m->renumberCells(old2NewBg,check);
611   setMesh(m);
612   updateTime();
613 }
614
615 /*!
616  * Permutes values of \a this field according to a given permutation array for cells
617  * renumbering. The underlying mesh is \b not permuted. 
618  * The number of cells remains the same; for that the permutation array \a old2NewBg
619  * should not contain equal ids.
620  * This method performs a part of job of renumberCells(). The reasonable use of this
621  * method is only for multi-field instances lying on the same mesh to avoid a
622  * systematic duplication and renumbering of _mesh attribute. 
623  * \warning Use this method with a lot of care!
624  *  \param [in] old2NewBg - the permutation array in "Old to New" mode. Its length is
625  *         to be equal to \a this->getMesh()->getNumberOfCells().
626  *  \param [in] check - if \c true, \a old2NewBg is transformed to a new permutation
627  *         array, so that its maximal cell id to correspond to (be less than) the number
628  *         of cells in mesh. This new array is then used for the renumbering. If \a 
629  *         check == \c false, \a old2NewBg is used as is, that is less secure as validity 
630  *         of ids in \a old2NewBg is not checked.
631  *  \throw If the mesh is not set.
632  *  \throw If the spatial discretization of \a this field is NULL.
633  *  \throw If \a check == \c true and \a old2NewBg contains equal ids.
634  *  \throw If mesh nature does not allow renumbering (e.g. structured mesh).
635  */
636 void MEDCouplingFieldDouble::renumberCellsWithoutMesh(const int *old2NewBg, bool check)
637 {
638    if(!_mesh)
639      throw INTERP_KERNEL::Exception("Expecting a defined mesh to be able to operate a renumbering !");
640    if(!((const MEDCouplingFieldDiscretization *)_type))
641      throw INTERP_KERNEL::Exception("Expecting a spatial discretization to be able to operate a renumbering !");
642   //
643   _type->renumberCells(old2NewBg,check);
644   std::vector<DataArrayDouble *> arrays;
645   _time_discr->getArrays(arrays);
646   std::vector<DataArray *> arrays2(arrays.size()); std::copy(arrays.begin(),arrays.end(),arrays2.begin());
647   _type->renumberArraysForCell(_mesh,arrays2,old2NewBg,check);
648   //
649   updateTime();
650 }
651
652 /*!
653  * Permutes values of \a this field according to a given permutation array for node
654  * renumbering. The underlying mesh is deeply copied and its nodes are also permuted. 
655  * The number of nodes can change, contrary to renumberCells().
656  *  \param [in] old2NewBg - the permutation array in "Old to New" mode. Its length is
657  *         to be equal to \a this->getMesh()->getNumberOfNodes().
658  *  \param [in] eps - a precision used to compare field values at merged nodes. If
659  *         the values differ more than \a eps, an exception is thrown.
660  *  \throw If the mesh is not set.
661  *  \throw If the spatial discretization of \a this field is NULL.
662  *  \throw If \a check == \c true and \a old2NewBg contains equal ids.
663  *  \throw If mesh nature does not allow renumbering (e.g. structured mesh).
664  *  \throw If values at merged nodes deffer more than \a eps.
665  * 
666  *  \if ENABLE_EXAMPLES
667  *  \ref cpp_mcfielddouble_renumberNodes "Here is a C++ example".<br>
668  *  \ref  py_mcfielddouble_renumberNodes "Here is a Python example".
669  *  \endif
670  */
671 void MEDCouplingFieldDouble::renumberNodes(const int *old2NewBg, double eps)
672 {
673   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
674   if(!meshC)
675     throw INTERP_KERNEL::Exception("Invalid mesh to apply renumberNodes on it !");
676   int nbOfNodes=meshC->getNumberOfNodes();
677   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
678   int newNbOfNodes=*std::max_element(old2NewBg,old2NewBg+nbOfNodes)+1;
679   renumberNodesWithoutMesh(old2NewBg,newNbOfNodes,eps);
680   meshC2->renumberNodes(old2NewBg,newNbOfNodes);
681   setMesh(meshC2);
682 }
683
684 /*!
685  * Permutes values of \a this field according to a given permutation array for nodes
686  * renumbering. The underlying mesh is \b not permuted. 
687  * The number of nodes can change, contrary to renumberCells().
688  * A given epsilon specifies a threshold of error in case of two nodes are merged but
689  * the difference of values on these nodes are higher than \a eps.
690  * This method performs a part of job of renumberNodes(), excluding node renumbering
691  * in mesh. The reasonable use of this
692  * method is only for multi-field instances lying on the same mesh to avoid a
693  * systematic duplication and renumbering of _mesh attribute. 
694  * \warning Use this method with a lot of care!
695  * \warning In case of an exception thrown, the contents of the data array can be
696  *         partially modified until the exception occurs. 
697  *  \param [in] old2NewBg - the permutation array in "Old to New" mode. Its length is
698  *         to be equal to \a this->getMesh()->getNumberOfNodes().
699  *  \param [in] newNbOfNodes - a number of nodes in the mesh after renumbering.
700  *  \param [in] eps - a precision used to compare field values at merged nodes. If
701  *         the values differ more than \a eps, an exception is thrown.
702  *  \throw If the mesh is not set.
703  *  \throw If the spatial discretization of \a this field is NULL.
704  *  \throw If values at merged nodes deffer more than \a eps.
705  */
706 void MEDCouplingFieldDouble::renumberNodesWithoutMesh(const int *old2NewBg, int newNbOfNodes, double eps)
707 {
708   if(!((const MEDCouplingFieldDiscretization *)_type))
709     throw INTERP_KERNEL::Exception("Expecting a spatial discretization to be able to operate a renumbering !");
710   std::vector<DataArrayDouble *> arrays;
711   _time_discr->getArrays(arrays);
712   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
713     if(*iter)
714       _type->renumberValuesOnNodes(eps,old2NewBg,newNbOfNodes,*iter);
715 }
716
717 /*!
718  * Returns all tuple ids of \a this scalar field that fit the range [\a vmin,
719  * \a vmax]. This method calls DataArrayDouble::getIdsInRange().
720  *  \param [in] vmin - a lower boundary of the range. Tuples with values less than \a
721  *         vmin are not included in the result array.
722  *  \param [in] vmax - an upper boundary of the range. Tuples with values more than \a
723  *         vmax are not included in the result array.
724  *  \return DataArrayInt * - a new instance of DataArrayInt holding ids of selected
725  *          tuples. The caller is to delete this array using decrRef() as it is no
726  *          more needed.
727  *  \throw If the data array is not set.
728  *  \throw If \a this->getNumberOfComponents() != 1.
729  */
730 DataArrayInt *MEDCouplingFieldDouble::getIdsInRange(double vmin, double vmax) const
731 {
732   if(getArray()==0)
733     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getIdsInRange : no default array set !");
734   return getArray()->getIdsInRange(vmin,vmax);
735 }
736
737 /*!
738  * Builds a newly created field, that the caller will have the responsability to deal with (decrRef()).
739  * This method makes the assumption that the field is correctly defined when this method is called, no check of this will be done.
740  * This method returns a restriction of \a this so that only tuples with ids specified in \a part will be contained in the returned field.
741  * Parameter \a part specifies **cell ids whatever the spatial discretization of this** (
742  * \ref ParaMEDMEM::ON_CELLS "ON_CELLS", 
743  * \ref ParaMEDMEM::ON_NODES "ON_NODES",
744  * \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT", 
745  * \ref ParaMEDMEM::ON_GAUSS_NE "ON_GAUSS_NE",
746  * \ref ParaMEDMEM::ON_NODES_KR "ON_NODES_KR").
747  *
748  * For example, \a this is a field on cells lying on a mesh that have 10 cells, \a part contains following cell ids [3,7,6].
749  * Then the returned field will lie on mesh having 3 cells and the returned field will contain 3 tuples.<br>
750  * Tuple #0 of the result field will refer to the cell #0 of returned mesh. The cell #0 of returned mesh will be equal to the cell #3 of \a this->getMesh().<br>
751  * Tuple #1 of the result field will refer to the cell #1 of returned mesh. The cell #1 of returned mesh will be equal to the cell #7 of \a this->getMesh().<br>
752  * Tuple #2 of the result field will refer to the cell #2 of returned mesh. The cell #2 of returned mesh will be equal to the cell #6 of \a this->getMesh().
753  *
754  * Let, for example, \a this be a field on nodes lying on a mesh that have 10 cells and 11 nodes, and \a part contains following cellIds [3,7,6].
755  * Thus \a this currently contains 11 tuples. If the restriction of mesh to 3 cells leads to a mesh with 6 nodes, then the returned field
756  * will contain 6 tuples and \a this field will lie on this restricted mesh. 
757  *
758  *  \param [in] part - an array of cell ids to include to the result field.
759  *  \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble. The caller is to delete this field using decrRef() as it is no more needed.
760  *
761  *  \if ENABLE_EXAMPLES
762  *  \ref cpp_mcfielddouble_subpart1 "Here is a C++ example".<br>
763  *  \ref  py_mcfielddouble_subpart1 "Here is a Python example".
764  *  \endif
765  *  \sa MEDCouplingFieldDouble::buildSubPartRange
766  */
767
768 MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPart(const DataArrayInt *part) const
769 {
770   if(part==0)
771     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::buildSubPart : not empty array must be passed to this method !");
772   return buildSubPart(part->begin(),part->end());
773 }
774
775 /*!
776  * Builds a newly created field, that the caller will have the responsability to deal with.
777  * \n This method makes the assumption that \a this field is correctly defined when this method is called (\a this->checkCoherency() returns without any exception thrown), **no check of this will be done**.
778  * \n This method returns a restriction of \a this so that only tuple ids specified in [ \a partBg , \a partEnd ) will be contained in the returned field.
779  * \n Parameter [\a partBg, \a partEnd ) specifies **cell ids whatever the spatial discretization** of \a this (
780  * \ref ParaMEDMEM::ON_CELLS "ON_CELLS", 
781  * \ref ParaMEDMEM::ON_NODES "ON_NODES",
782  * \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT", 
783  * \ref ParaMEDMEM::ON_GAUSS_NE "ON_GAUSS_NE",
784  * \ref ParaMEDMEM::ON_NODES_KR "ON_NODES_KR").
785  *
786  * For example, \a this is a field on cells lying on a mesh that have 10 cells, \a partBg contains the following cell ids [3,7,6].
787  * Then the returned field will lie on mesh having 3 cells and will contain 3 tuples.
788  *- Tuple #0 of the result field will refer to the cell #0 of returned mesh. The cell #0 of returned mesh will be equal to the cell #3 of \a this->getMesh().
789  *- Tuple #1 of the result field will refer to the cell #1 of returned mesh. The cell #1 of returned mesh will be equal to the cell #7 of \a this->getMesh().
790  *- Tuple #2 of the result field will refer to the cell #2 of returned mesh. The cell #2 of returned mesh will be equal to the cell #6 of \a this->getMesh().
791  *
792  * Let, for example, \a this be a field on nodes lying on a mesh that have 10 cells and 11 nodes, and \a partBg contains following cellIds [3,7,6].
793  * Thus \a this currently contains 11 tuples. If the restriction of mesh to 3 cells leads to a mesh with 6 nodes, then the returned field
794  * will contain 6 tuples and \a this field will lie on this restricted mesh. 
795  *
796  * \param [in] partBg - start (included) of input range of cell ids to select [ \a partBg, \a partEnd )
797  * \param [in] partEnd - end (not included) of input range of cell ids to select [ \a partBg, \a partEnd )
798  * \return a newly allocated field the caller should deal with.
799  * 
800  * \throw if there is presence of an invalid cell id in [ \a partBg, \a partEnd ) regarding the number of cells of \a this->getMesh().
801  *
802  * \if ENABLE_EXAMPLES
803  * \ref cpp_mcfielddouble_subpart1 "Here a C++ example."<br>
804  * \ref py_mcfielddouble_subpart1 "Here a Python example."
805  * \endif
806  * \sa ParaMEDMEM::MEDCouplingFieldDouble::buildSubPart(const DataArrayInt *) const, MEDCouplingFieldDouble::buildSubPartRange
807  */
808 MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPart(const int *partBg, const int *partEnd) const
809 {
810   if(!((const MEDCouplingFieldDiscretization *)_type))
811     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::buildSubPart : Expecting a not NULL spatial discretization !");
812   DataArrayInt *arrSelect;
813   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> m=_type->buildSubMeshData(_mesh,partBg,partEnd,arrSelect);
814   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arrSelect2(arrSelect);
815   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=clone(false);//quick shallow copy.
816   const MEDCouplingFieldDiscretization *disc=getDiscretization();
817   if(disc)
818     ret->setDiscretization(MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDiscretization>(disc->clonePart(partBg,partEnd)));
819   ret->setMesh(m);
820   std::vector<DataArrayDouble *> arrays;
821   _time_discr->getArrays(arrays);
822   std::vector<DataArrayDouble *> arrs;
823   std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > arrsSafe;
824   const int *arrSelBg=arrSelect->begin();
825   const int *arrSelEnd=arrSelect->end();
826   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
827     {
828       DataArrayDouble *arr=0;
829       if(*iter)
830         arr=(*iter)->selectByTupleIdSafe(arrSelBg,arrSelEnd);
831       arrs.push_back(arr); arrsSafe.push_back(arr);
832     }
833   ret->_time_discr->setArrays(arrs,0);
834   return ret.retn();
835 }
836
837 /*!
838  * This method is equivalent to MEDCouplingFieldDouble::buildSubPart, the only difference is that the input range of cell ids is
839  * given using a range given \a begin, \a end and \a step to optimize the part computation.
840  * 
841  * \sa MEDCouplingFieldDouble::buildSubPart
842  */
843 MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPartRange(int begin, int end, int step) const
844 {
845   if(!((const MEDCouplingFieldDiscretization *)_type))
846     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::buildSubPart : Expecting a not NULL spatial discretization !");
847   DataArrayInt *arrSelect;
848   int beginOut,endOut,stepOut;
849   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> m=_type->buildSubMeshDataRange(_mesh,begin,end,step,beginOut,endOut,stepOut,arrSelect);
850   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arrSelect2(arrSelect);
851   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=clone(false);//quick shallow copy.
852   const MEDCouplingFieldDiscretization *disc=getDiscretization();
853   if(disc)
854     ret->setDiscretization(MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDiscretization>(disc->clonePartRange(begin,end,step)));
855   ret->setMesh(m);
856   std::vector<DataArrayDouble *> arrays;
857   _time_discr->getArrays(arrays);
858   std::vector<DataArrayDouble *> arrs;
859   std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > arrsSafe;
860   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
861     {
862       DataArrayDouble *arr=0;
863       if(*iter)
864         {
865           if(arrSelect)
866             {
867               const int *arrSelBg=arrSelect->begin();
868               const int *arrSelEnd=arrSelect->end();
869               arr=(*iter)->selectByTupleIdSafe(arrSelBg,arrSelEnd);
870             }
871           else
872             arr=(*iter)->selectByTupleId2(beginOut,endOut,stepOut);
873         }
874       arrs.push_back(arr); arrsSafe.push_back(arr);
875     }
876   ret->_time_discr->setArrays(arrs,0);
877   return ret.retn();
878 }
879
880 /*!
881  * Returns a type of \ref MEDCouplingTemporalDisc "time discretization" of \a this field.
882  *  \return ParaMEDMEM::TypeOfTimeDiscretization - an enum item describing the time
883  *          discretization type.
884  */
885 TypeOfTimeDiscretization MEDCouplingFieldDouble::getTimeDiscretization() const
886 {
887   return _time_discr->getEnum();
888 }
889
890 MEDCouplingFieldDouble::MEDCouplingFieldDouble(TypeOfField type, TypeOfTimeDiscretization td):MEDCouplingField(type),
891                                                                                               _time_discr(MEDCouplingTimeDiscretization::New(td))
892 {
893 }
894
895 /*!
896  * ** WARINING : This method do not deeply copy neither mesh nor spatial discretization. Only a shallow copy (reference) is done for mesh and spatial discretization ! **
897  */
898 MEDCouplingFieldDouble::MEDCouplingFieldDouble(const MEDCouplingFieldTemplate& ft, TypeOfTimeDiscretization td):MEDCouplingField(ft,false),
899                                                                                                                 _time_discr(MEDCouplingTimeDiscretization::New(td))
900 {
901 }
902
903 MEDCouplingFieldDouble::MEDCouplingFieldDouble(const MEDCouplingFieldDouble& other, bool deepCopy):MEDCouplingField(other,deepCopy),
904                                                                                                    _time_discr(other._time_discr->performCpy(deepCopy))
905 {
906 }
907
908 MEDCouplingFieldDouble::MEDCouplingFieldDouble(NatureOfField n, MEDCouplingTimeDiscretization *td, MEDCouplingFieldDiscretization *type):MEDCouplingField(type,n),_time_discr(td)
909 {
910 }
911
912 MEDCouplingFieldDouble::~MEDCouplingFieldDouble()
913 {
914   delete _time_discr;
915 }
916
917 /*!
918  * Checks if \a this field is correctly defined, else an exception is thrown.
919  *  \throw If the mesh is not set.
920  *  \throw If the data array is not set.
921  *  \throw If the spatial discretization of \a this field is NULL.
922  *  \throw If \a this->getTimeTolerance() < 0.
923  *  \throw If the temporal discretization data is incorrect.
924  *  \throw If mesh data does not correspond to field data.
925  */
926 void MEDCouplingFieldDouble::checkCoherency() const
927 {
928   if(!_mesh)
929     throw INTERP_KERNEL::Exception("Field invalid because no mesh specified !");
930   if(!((const MEDCouplingFieldDiscretization *)_type))
931     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::checkCoherency : no spatial discretization !");
932   _time_discr->checkCoherency();
933   _type->checkCoherencyBetween(_mesh,getArray());
934 }
935
936 /*!
937  * Accumulate values of a given component of \a this field.
938  *  \param [in] compId - the index of the component of interest.
939  *  \return double - a sum value of *compId*-th component.
940  *  \throw If the data array is not set.
941  *  \throw If \a the condition ( 0 <= \a compId < \a this->getNumberOfComponents() ) is
942  *         not respected.
943  */
944 double MEDCouplingFieldDouble::accumulate(int compId) const
945 {
946   if(getArray()==0)
947     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::accumulate : no default array defined !");
948   return getArray()->accumulate(compId);
949 }
950
951 /*!
952  * Accumulates values of each component of \a this array.
953  *  \param [out] res - an array of length \a this->getNumberOfComponents(), allocated 
954  *         by the caller, that is filled by this method with sum value for each
955  *         component.
956  *  \throw If the data array is not set.
957  */
958 void MEDCouplingFieldDouble::accumulate(double *res) const
959 {
960   if(getArray()==0)
961     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::accumulate : no default array defined !");
962   getArray()->accumulate(res);
963 }
964
965 /*!
966  * Returns the maximal value within \a this scalar field. Values of all arrays stored
967  * in \a this->_time_discr are checked.
968  *  \return double - the maximal value among all values of \a this field.
969  *  \throw If \a this->getNumberOfComponents() != 1
970  *  \throw If the data array is not set.
971  *  \throw If there is an empty data array in \a this field.
972  */
973 double MEDCouplingFieldDouble::getMaxValue() const
974 {
975   std::vector<DataArrayDouble *> arrays;
976   _time_discr->getArrays(arrays);
977   double ret=-std::numeric_limits<double>::max();
978   bool isExistingArr=false;
979   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
980     {
981       if(*iter)
982         {
983           isExistingArr=true;
984           int loc;
985           ret=std::max(ret,(*iter)->getMaxValue(loc));
986         }
987     }
988   if(!isExistingArr)
989     throw INTERP_KERNEL::Exception("getMaxValue : No arrays defined !");
990   return ret;
991 }
992
993 /*!
994  * Returns the maximal value and all its locations within \a this scalar field.
995  * Only the first of available data arrays is checked.
996  *  \param [out] tupleIds - a new instance of DataArrayInt containg indices of
997  *               tuples holding the maximal value. The caller is to delete it using
998  *               decrRef() as it is no more needed.
999  *  \return double - the maximal value among all values of the first array of \a this filed.
1000  *  \throw If \a this->getNumberOfComponents() != 1.
1001  *  \throw If there is an empty data array in \a this field.
1002  */
1003 double MEDCouplingFieldDouble::getMaxValue2(DataArrayInt*& tupleIds) const
1004 {
1005   std::vector<DataArrayDouble *> arrays;
1006   _time_discr->getArrays(arrays);
1007   double ret=-std::numeric_limits<double>::max();
1008   bool isExistingArr=false;
1009   tupleIds=0;
1010   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret1;
1011   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1012     {
1013       if(*iter)
1014         {
1015           isExistingArr=true;
1016           DataArrayInt *tmp;
1017           ret=std::max(ret,(*iter)->getMaxValue2(tmp));
1018           MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmpSafe(tmp);
1019           if(!((const DataArrayInt *)ret1))
1020             ret1=tmpSafe;
1021         }
1022     }
1023   if(!isExistingArr)
1024     throw INTERP_KERNEL::Exception("getMaxValue2 : No arrays defined !");
1025   tupleIds=ret1.retn();
1026   return ret;
1027 }
1028
1029 /*!
1030  * Returns the minimal value within \a this scalar field. Values of all arrays stored
1031  * in \a this->_time_discr are checked.
1032  *  \return double - the minimal value among all values of \a this field.
1033  *  \throw If \a this->getNumberOfComponents() != 1
1034  *  \throw If the data array is not set.
1035  *  \throw If there is an empty data array in \a this field.
1036  */
1037 double MEDCouplingFieldDouble::getMinValue() const
1038 {
1039   std::vector<DataArrayDouble *> arrays;
1040   _time_discr->getArrays(arrays);
1041   double ret=std::numeric_limits<double>::max();
1042   bool isExistingArr=false;
1043   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1044     {
1045       if(*iter)
1046         {
1047           isExistingArr=true;
1048           int loc;
1049           ret=std::min(ret,(*iter)->getMinValue(loc));
1050         }
1051     }
1052   if(!isExistingArr)
1053     throw INTERP_KERNEL::Exception("getMinValue : No arrays defined !");
1054   return ret;
1055 }
1056
1057 /*!
1058  * Returns the minimal value and all its locations within \a this scalar field.
1059  * Only the first of available data arrays is checked.
1060  *  \param [out] tupleIds - a new instance of DataArrayInt containg indices of
1061  *               tuples holding the minimal value. The caller is to delete it using
1062  *               decrRef() as it is no more needed.
1063  *  \return double - the minimal value among all values of the first array of \a this filed.
1064  *  \throw If \a this->getNumberOfComponents() != 1.
1065  *  \throw If there is an empty data array in \a this field.
1066  */
1067 double MEDCouplingFieldDouble::getMinValue2(DataArrayInt*& tupleIds) const
1068 {
1069   std::vector<DataArrayDouble *> arrays;
1070   _time_discr->getArrays(arrays);
1071   double ret=-std::numeric_limits<double>::max();
1072   bool isExistingArr=false;
1073   tupleIds=0;
1074   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret1;
1075   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
1076     {
1077       if(*iter)
1078         {
1079           isExistingArr=true;
1080           DataArrayInt *tmp;
1081           ret=std::max(ret,(*iter)->getMinValue2(tmp));
1082           MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmpSafe(tmp);
1083           if(!((const DataArrayInt *)ret1))
1084             ret1=tmpSafe;
1085         }
1086     }
1087   if(!isExistingArr)
1088     throw INTERP_KERNEL::Exception("getMinValue2 : No arrays defined !");
1089   tupleIds=ret1.retn();
1090   return ret;
1091 }
1092
1093 /*!
1094  * Returns the average value of \a this scalar field.
1095  *  \return double - the average value over all values of the data array.
1096  *  \throw If \a this->getNumberOfComponents() != 1
1097  *  \throw If the data array is not set or it is empty.
1098  */
1099 double MEDCouplingFieldDouble::getAverageValue() const
1100 {
1101   if(getArray()==0)
1102     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getAverageValue : no default array defined !");
1103   return getArray()->getAverageValue();
1104 }
1105
1106 /*!
1107  * This method returns the euclidean norm of \a this field.
1108  * \f[
1109  * \sqrt{\sum_{0 \leq i < nbOfEntity}val[i]*val[i]}
1110  * \f]
1111  *  \throw If the data array is not set.
1112  */
1113 double MEDCouplingFieldDouble::norm2() const
1114 {
1115   if(getArray()==0)
1116     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::norm2 : no default array defined !");
1117   return getArray()->norm2();
1118 }
1119
1120 /*!
1121  * This method returns the max norm of \a this field.
1122  * \f[
1123  * \max_{0 \leq i < nbOfEntity}{abs(val[i])}
1124  * \f]
1125  *  \throw If the data array is not set.
1126  */
1127 double MEDCouplingFieldDouble::normMax() const
1128 {
1129   if(getArray()==0)
1130     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::normMax : no default array defined !");
1131   return getArray()->normMax();
1132 }
1133
1134 /*!
1135  * Computes sums of values of each component of \a this field wighted with
1136  * values returned by buildMeasureField().  
1137  *  \param [out] res - pointer to an array of result sum values, of size at least \a
1138  *         this->getNumberOfComponents(), that is to be allocated by the caller.
1139  *  \param [in] isWAbs - if \c true (default), \c abs() is applied to the weighs computed by
1140  *         buildMeasureField() that makes this method slower. If a user is sure that all
1141  *         cells of the underlying mesh have correct orientation, he can put \a isWAbs ==
1142  *         \c false that speeds up this method.
1143  *  \throw If the mesh is not set.
1144  *  \throw If the data array is not set.
1145  */
1146 void MEDCouplingFieldDouble::getWeightedAverageValue(double *res, bool isWAbs) const
1147 {
1148   if(getArray()==0)
1149     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getWeightedAverageValue : no default array defined !");
1150   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> w=buildMeasureField(isWAbs);
1151   double deno=w->getArray()->accumulate(0);
1152   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr=getArray()->deepCpy();
1153   arr->multiplyEqual(w->getArray());
1154   std::transform(arr->begin(),arr->end(),arr->getPointer(),std::bind2nd(std::multiplies<double>(),1./deno));
1155   arr->accumulate(res);
1156 }
1157
1158 /*!
1159  * Computes a sum of values of a given component of \a this field wighted with
1160  * values returned by buildMeasureField().
1161  *  \param [in] compId - an index of the component of interest.
1162  *  \param [in] isWAbs - if \c true (default), \c abs() is applied to the weighs computed by
1163  *         buildMeasureField() that makes this method slower. If a user is sure that all
1164  *         cells of the underlying mesh have correct orientation, he can put \a isWAbs ==
1165  *         \c false that speeds up this method.
1166  *  \throw If the mesh is not set.
1167  *  \throw If the data array is not set.
1168  *  \throw If \a compId is not valid.
1169            A valid range is ( 0 <= \a compId < \a this->getNumberOfComponents() ).
1170  */
1171 double MEDCouplingFieldDouble::getWeightedAverageValue(int compId, bool isWAbs) const
1172 {
1173   int nbComps=getArray()->getNumberOfComponents();
1174   if(compId<0 || compId>=nbComps)
1175     {
1176       std::ostringstream oss; oss << "MEDCouplingFieldDouble::getWeightedAverageValue : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
1177       throw INTERP_KERNEL::Exception(oss.str().c_str());
1178     }
1179   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
1180   getWeightedAverageValue(res,isWAbs);
1181   return res[compId];
1182 }
1183
1184 /*!
1185  * Returns the \c normL1 of values of a given component of \a this field:
1186  * \f[
1187  * \frac{\sum_{0 \leq i < nbOfEntity}|val[i]*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}
1188  * \f]
1189  *  \param [in] compId - an index of the component of interest.
1190  *  \throw If the mesh is not set.
1191  *  \throw If the spatial discretization of \a this field is NULL.
1192  *  \throw If \a compId is not valid.
1193            A valid range is ( 0 <= \a compId < \a this->getNumberOfComponents() ).
1194  */
1195 double MEDCouplingFieldDouble::normL1(int compId) const
1196 {
1197   if(!_mesh)
1198     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1 !");
1199   if(!((const MEDCouplingFieldDiscretization *)_type))
1200     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform normL1 !");
1201   int nbComps=getArray()->getNumberOfComponents();
1202   if(compId<0 || compId>=nbComps)
1203     {
1204       std::ostringstream oss; oss << "MEDCouplingFieldDouble::normL1 : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
1205       throw INTERP_KERNEL::Exception(oss.str().c_str());
1206     }
1207   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
1208   _type->normL1(_mesh,getArray(),res);
1209   return res[compId];
1210 }
1211
1212 /*!
1213  * Returns the \c normL1 of values of each component of \a this field:
1214  * \f[
1215  * \frac{\sum_{0 \leq i < nbOfEntity}|val[i]*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}
1216  * \f]
1217  *  \param [out] res - pointer to an array of result values, of size at least \a
1218  *         this->getNumberOfComponents(), that is to be allocated by the caller.
1219  *  \throw If the mesh is not set.
1220  *  \throw If the spatial discretization of \a this field is NULL.
1221  */
1222 void MEDCouplingFieldDouble::normL1(double *res) const
1223 {
1224   if(!_mesh)
1225     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1");
1226   if(!((const MEDCouplingFieldDiscretization *)_type))
1227     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform normL1 !");
1228   _type->normL1(_mesh,getArray(),res);
1229 }
1230
1231 /*!
1232  * Returns the \c normL2 of values of a given component of \a this field:
1233  * \f[
1234  * \sqrt{\frac{\sum_{0 \leq i < nbOfEntity}|val[i]^{2}*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}}
1235  * \f]
1236  *  \param [in] compId - an index of the component of interest.
1237  *  \throw If the mesh is not set.
1238  *  \throw If the spatial discretization of \a this field is NULL.
1239  *  \throw If \a compId is not valid.
1240            A valid range is ( 0 <= \a compId < \a this->getNumberOfComponents() ).
1241  */
1242 double MEDCouplingFieldDouble::normL2(int compId) const
1243 {
1244   if(!_mesh)
1245     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL2");
1246   if(!((const MEDCouplingFieldDiscretization *)_type))
1247     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform normL2 !");
1248   int nbComps=getArray()->getNumberOfComponents();
1249   if(compId<0 || compId>=nbComps)
1250     {
1251       std::ostringstream oss; oss << "MEDCouplingFieldDouble::normL2 : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
1252       throw INTERP_KERNEL::Exception(oss.str().c_str());
1253     }
1254   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
1255   _type->normL2(_mesh,getArray(),res);
1256   return res[compId];
1257 }
1258
1259 /*!
1260  * Returns the \c normL2 of values of each component of \a this field:
1261  * \f[
1262  * \sqrt{\frac{\sum_{0 \leq i < nbOfEntity}|val[i]^{2}*Vol[i]|}{\sum_{0 \leq i < nbOfEntity}|Vol[i]|}}
1263  * \f]
1264  *  \param [out] res - pointer to an array of result values, of size at least \a
1265  *         this->getNumberOfComponents(), that is to be allocated by the caller.
1266  *  \throw If the mesh is not set.
1267  *  \throw If the spatial discretization of \a this field is NULL.
1268  */
1269 void MEDCouplingFieldDouble::normL2(double *res) const
1270 {
1271   if(!_mesh)
1272     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL2");
1273   if(!((const MEDCouplingFieldDiscretization *)_type))
1274     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform normL2 !");
1275   _type->normL2(_mesh,getArray(),res);
1276 }
1277
1278 /*!
1279  * Computes a sum of values of a given component of \a this field multiplied by
1280  * values returned by buildMeasureField().
1281  * This method is useful to check the conservativity of interpolation method.
1282  *  \param [in] compId - an index of the component of interest.
1283  *  \param [in] isWAbs - if \c true (default), \c abs() is applied to the weighs computed by
1284  *         buildMeasureField() that makes this method slower. If a user is sure that all
1285  *         cells of the underlying mesh have correct orientation, he can put \a isWAbs ==
1286  *         \c false that speeds up this method.
1287  *  \throw If the mesh is not set.
1288  *  \throw If the data array is not set.
1289  *  \throw If \a compId is not valid.
1290            A valid range is ( 0 <= \a compId < \a this->getNumberOfComponents() ).
1291  */
1292 double MEDCouplingFieldDouble::integral(int compId, bool isWAbs) const
1293 {
1294   if(!_mesh)
1295     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform integral");
1296   if(!((const MEDCouplingFieldDiscretization *)_type))
1297     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform integral !");
1298   int nbComps=getArray()->getNumberOfComponents();
1299   if(compId<0 || compId>=nbComps)
1300     {
1301       std::ostringstream oss; oss << "MEDCouplingFieldDouble::integral : Invalid compId specified : No such nb of components ! Should be in [0," << nbComps << ") !";
1302       throw INTERP_KERNEL::Exception(oss.str().c_str());
1303     }
1304   INTERP_KERNEL::AutoPtr<double> res=new double[nbComps];
1305   _type->integral(_mesh,getArray(),isWAbs,res);
1306   return res[compId];
1307 }
1308
1309 /*!
1310  * Computes a sum of values of each component of \a this field multiplied by
1311  * values returned by buildMeasureField().
1312  * This method is useful to check the conservativity of interpolation method.
1313  *  \param [in] isWAbs - if \c true (default), \c abs() is applied to the weighs computed by
1314  *         buildMeasureField() that makes this method slower. If a user is sure that all
1315  *         cells of the underlying mesh have correct orientation, he can put \a isWAbs ==
1316  *         \c false that speeds up this method.
1317  *  \param [out] res - pointer to an array of result sum values, of size at least \a
1318  *         this->getNumberOfComponents(), that is to be allocated by the caller.
1319  *  \throw If the mesh is not set.
1320  *  \throw If the data array is not set.
1321  *  \throw If the spatial discretization of \a this field is NULL.
1322  */
1323 void MEDCouplingFieldDouble::integral(bool isWAbs, double *res) const
1324 {
1325   if(!_mesh)
1326     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform integral2");
1327   if(!((const MEDCouplingFieldDiscretization *)_type))
1328     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform integral2 !");
1329   _type->integral(_mesh,getArray(),isWAbs,res);
1330 }
1331
1332 /*!
1333  * Returns a value at a given cell of a structured mesh. The cell is specified by its
1334  * (i,j,k) index.
1335  *  \param [in] i - a index of node coordinates array along X axis. The cell is
1336  *         located between the i-th and ( i + 1 )-th nodes along X axis.
1337  *  \param [in] j - a index of node coordinates array along Y axis. The cell is
1338  *         located between the j-th and ( j + 1 )-th nodes along Y axis.
1339  *  \param [in] k - a index of node coordinates array along Z axis. The cell is
1340  *         located between the k-th and ( k + 1 )-th nodes along Z axis.
1341  *  \param [out] res - pointer to an array returning a feild value, of size at least
1342  *         \a this->getNumberOfComponents(), that is to be allocated by the caller.
1343  *  \throw If the spatial discretization of \a this field is NULL.
1344  *  \throw If the mesh is not set.
1345  *  \throw If the mesh is not a structured one.
1346  *
1347  *  \if ENABLE_EXAMPLES
1348  *  \ref cpp_mcfielddouble_getValueOnPos "Here is a C++ example".<br>
1349  *  \ref  py_mcfielddouble_getValueOnPos "Here is a Python example".
1350  *  \endif
1351  */
1352 void MEDCouplingFieldDouble::getValueOnPos(int i, int j, int k, double *res) const
1353 {
1354   const DataArrayDouble *arr=_time_discr->getArray();
1355   if(!_mesh)
1356     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOnPos");
1357   if(!((const MEDCouplingFieldDiscretization *)_type))
1358     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getValueOnPos !");
1359   _type->getValueOnPos(arr,_mesh,i,j,k,res);
1360 }
1361
1362 /*!
1363  * Returns a value of \a this at a given point using spatial discretization.
1364  *  \param [in] spaceLoc - the point of interest.
1365  *  \param [out] res - pointer to an array returning a feild value, of size at least
1366  *         \a this->getNumberOfComponents(), that is to be allocated by the caller.
1367  *  \throw If the spatial discretization of \a this field is NULL.
1368  *  \throw If the mesh is not set.
1369  *  \throw If \a spaceLoc is out of the spatial discretization.
1370  *
1371  *  \if ENABLE_EXAMPLES
1372  *  \ref cpp_mcfielddouble_getValueOn "Here is a C++ example".<br>
1373  *  \ref  py_mcfielddouble_getValueOn "Here is a Python example".
1374  *  \endif
1375  */
1376 void MEDCouplingFieldDouble::getValueOn(const double *spaceLoc, double *res) const
1377 {
1378   const DataArrayDouble *arr=_time_discr->getArray();
1379   if(!_mesh)
1380     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOn");
1381   if(!((const MEDCouplingFieldDiscretization *)_type))
1382     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getValueOnPos !");
1383   _type->getValueOn(arr,_mesh,spaceLoc,res);
1384 }
1385
1386 /*!
1387  * Returns values of \a this at given points using spatial discretization.
1388  *  \param [in] spaceLoc - coordinates of points of interest in full-interlace
1389  *          mode. This array is to be of size ( \a nbOfPoints * \a this->getNumberOfComponents() ).
1390  *  \param [in] nbOfPoints - number of points of interest.
1391  *  \return DataArrayDouble * - a new instance of DataArrayDouble holding field
1392  *         values relating to the input points. This array is of size \a nbOfPoints
1393  *         tuples per \a this->getNumberOfComponents() components. The caller is to 
1394  *         delete this array using decrRef() as it is no more needed.
1395  *  \throw If the spatial discretization of \a this field is NULL.
1396  *  \throw If the mesh is not set.
1397  *  \throw If any point in \a spaceLoc is out of the spatial discretization.
1398  *
1399  *  \if ENABLE_EXAMPLES
1400  *  \ref cpp_mcfielddouble_getValueOnMulti "Here is a C++ example".<br>
1401  *  \ref  py_mcfielddouble_getValueOnMulti "Here is a Python example".
1402  *  \endif
1403  */
1404 DataArrayDouble *MEDCouplingFieldDouble::getValueOnMulti(const double *spaceLoc, int nbOfPoints) const
1405 {
1406   const DataArrayDouble *arr=_time_discr->getArray();
1407   if(!_mesh)
1408     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOnMulti");
1409   if(!((const MEDCouplingFieldDiscretization *)_type))
1410     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getValueOnMulti !");
1411   return _type->getValueOnMulti(arr,_mesh,spaceLoc,nbOfPoints);
1412 }
1413
1414 /*!
1415  * Returns a value of \a this field at a given point at a given time using spatial discretization.
1416  * If the time is not covered by \a this->_time_discr, an exception is thrown.
1417  *  \param [in] spaceLoc - the point of interest.
1418  *  \param [in] time - the time of interest.
1419  *  \param [out] res - pointer to an array returning a feild value, of size at least
1420  *         \a this->getNumberOfComponents(), that is to be allocated by the caller.
1421  *  \throw If the spatial discretization of \a this field is NULL.
1422  *  \throw If the mesh is not set.
1423  *  \throw If \a spaceLoc is out of the spatial discretization.
1424  *  \throw If \a time is not covered by \a this->_time_discr.
1425  *
1426  *  \if ENABLE_EXAMPLES
1427  *  \ref cpp_mcfielddouble_getValueOn_time "Here is a C++ example".<br>
1428  *  \ref  py_mcfielddouble_getValueOn_time "Here is a Python example".
1429  *  \endif
1430  */
1431 void MEDCouplingFieldDouble::getValueOn(const double *spaceLoc, double time, double *res) const
1432 {
1433   std::vector< const DataArrayDouble *> arrs=_time_discr->getArraysForTime(time);
1434   if(!_mesh)
1435     throw INTERP_KERNEL::Exception("No mesh underlying this field to perform getValueOn");
1436   if(!((const MEDCouplingFieldDiscretization *)_type))
1437     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getValueOn !");
1438   std::vector<double> res2;
1439   for(std::vector< const DataArrayDouble *>::const_iterator iter=arrs.begin();iter!=arrs.end();iter++)
1440     {
1441       int sz=(int)res2.size();
1442       res2.resize(sz+(*iter)->getNumberOfComponents());
1443       _type->getValueOn(*iter,_mesh,spaceLoc,&res2[sz]);
1444     }
1445   _time_discr->getValueForTime(time,res2,res);
1446 }
1447
1448 /*!
1449  * Apply a liner function to a given component of \a this field, so that
1450  * a component value <em>(x)</em> becomes \f$ a * x + b \f$.
1451  *  \param [in] a - the first coefficient of the function.
1452  *  \param [in] b - the second coefficient of the function.
1453  *  \param [in] compoId - the index of component to modify.
1454  *  \throw If the data array(s) is(are) not set.
1455  */
1456 void MEDCouplingFieldDouble::applyLin(double a, double b, int compoId)
1457 {
1458   _time_discr->applyLin(a,b,compoId);
1459 }
1460
1461 /*!
1462  * This method sets \a this to a uniform scalar field with one component.
1463  * All tuples will have the same value 'value'.
1464  * An exception is thrown if no underlying mesh is defined.
1465  */
1466 MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator=(double value) throw(INTERP_KERNEL::Exception)
1467 {
1468   if(!_mesh)
1469     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::operator= : no mesh defined !");
1470   if(!((const MEDCouplingFieldDiscretization *)_type))
1471     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform operator = !");
1472   int nbOfTuple=_type->getNumberOfTuples(_mesh);
1473   _time_discr->setOrCreateUniformValueOnAllComponents(nbOfTuple,value);
1474   return *this;
1475 }
1476
1477 /*!
1478  * Creates data array(s) of \a this field by using a C function for value generation.
1479  *  \param [in] nbOfComp - the number of components for \a this field to have.
1480  *  \param [in] func - the function used to compute values of \a this field.
1481  *         This function is to compute a field value basing on coordinates of value
1482  *         location point.
1483  *  \throw If the mesh is not set.
1484  *  \throw If \a func returns \c false.
1485  *  \throw If the spatial discretization of \a this field is NULL.
1486  *
1487  *  \if ENABLE_EXAMPLES
1488  *  \ref cpp_mcfielddouble_fillFromAnalytic_c_func "Here is a C++ example".
1489  *  \endif
1490  */
1491 void MEDCouplingFieldDouble::fillFromAnalytic(int nbOfComp, FunctionToEvaluate func)
1492 {
1493   if(!_mesh)
1494     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic : no mesh defined !");
1495   if(!((const MEDCouplingFieldDiscretization *)_type))
1496     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform fillFromAnalytic !");
1497   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
1498   _time_discr->fillFromAnalytic(loc,nbOfComp,func);
1499 }
1500
1501 /*!
1502  * Creates data array(s) of \a this field by using a function for value generation.<br>
1503  * The function is applied to coordinates of value location points. For example, if
1504  * \a this field is on cells, the function is applied to cell barycenters.
1505  * For more info on supported expressions that can be used in the function, see \ref
1506  * MEDCouplingArrayApplyFuncExpr. <br>
1507  * The function can include arbitrary named variables
1508  * (e.g. "x","y" or "va44") to refer to components of point coordinates. Names of
1509  * variables are sorted in \b alphabetical \b order to associate a variable name with a
1510  * component. For example, in the expression "2*x+z", "x" stands for the component #0
1511  * and "z" stands for the component #1 (\b not #2)!<br>
1512  * In a general case, a value resulting from the function evaluation is assigned to all
1513  * components of a field value. But there is a possibility to have its own expression for
1514  * each component within one function. For this purpose, there are predefined variable
1515  * names (IVec, JVec, KVec, LVec etc) each dedicated to a certain component (IVec, to
1516  * the component #0 etc). A factor of such a variable is added to the
1517  * corresponding component only.<br>
1518  * For example, \a nbOfComp == 4, coordinates of a 3D point are (1.,3.,7.), then
1519  *   - "2*x + z"               produces (5.,5.,5.,5.)
1520  *   - "2*x + 0*y + z"         produces (9.,9.,9.,9.)
1521  *   - "2*x*IVec + (x+z)*LVec" produces (2.,0.,0.,4.)
1522  *   - "2*y*IVec + z*KVec + x" produces (7.,1.,1.,4.)
1523  *
1524  *  \param [in] nbOfComp - the number of components for \a this field to have.
1525  *  \param [in] func - the function used to compute values of \a this field.
1526  *         This function is used to compute a field value basing on coordinates of value
1527  *         location point. For example, if \a this field is on cells, the function
1528  *         is applied to cell barycenters.
1529  *  \throw If the mesh is not set.
1530  *  \throw If the spatial discretization of \a this field is NULL.
1531  *  \throw If computing \a func fails.
1532  *
1533  *  \if ENABLE_EXAMPLES
1534  *  \ref cpp_mcfielddouble_fillFromAnalytic "Here is a C++ example".<br>
1535  *  \ref  py_mcfielddouble_fillFromAnalytic "Here is a Python example".
1536  *  \endif
1537  */
1538 void MEDCouplingFieldDouble::fillFromAnalytic(int nbOfComp, const std::string& func)
1539 {
1540   if(!_mesh)
1541     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic : no mesh defined !");
1542   if(!((const MEDCouplingFieldDiscretization *)_type))
1543     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform fillFromAnalytic !");
1544   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
1545   _time_discr->fillFromAnalytic(loc,nbOfComp,func);
1546 }
1547
1548 /*!
1549  * Creates data array(s) of \a this field by using a function for value generation.<br>
1550  * The function is applied to coordinates of value location points. For example, if
1551  * \a this field is on cells, the function is applied to cell barycenters.<br>
1552  * This method differs from
1553  * \ref ParaMEDMEM::MEDCouplingFieldDouble::fillFromAnalytic(int nbOfComp, const std::string& func) "fillFromAnalytic()"
1554  * by the way how variable
1555  * names, used in the function, are associated with components of coordinates of field
1556  * location points; here, a variable name corresponding to a component is retrieved from
1557  * a corresponding node coordinates array (where it is set via
1558  * DataArrayDouble::setInfoOnComponent()).<br>
1559  * For more info on supported expressions that can be used in the function, see \ref
1560  * MEDCouplingArrayApplyFuncExpr. <br> 
1561  * In a general case, a value resulting from the function evaluation is assigned to all
1562  * components of a field value. But there is a possibility to have its own expression for
1563  * each component within one function. For this purpose, there are predefined variable
1564  * names (IVec, JVec, KVec, LVec etc) each dedicated to a certain component (IVec, to
1565  * the component #0 etc). A factor of such a variable is added to the
1566  * corresponding component only.<br>
1567  * For example, \a nbOfComp == 4, names of spatial components are "x", "y" and "z",
1568  * coordinates of a 3D point are (1.,3.,7.), then
1569  *   - "2*x + z"               produces (9.,9.,9.,9.)
1570  *   - "2*x*IVec + (x+z)*LVec" produces (2.,0.,0.,8.)
1571  *   - "2*y*IVec + z*KVec + x" produces (7.,1.,1.,8.)
1572  *
1573  *  \param [in] nbOfComp - the number of components for \a this field to have.
1574  *  \param [in] func - the function used to compute values of \a this field.
1575  *         This function is used to compute a field value basing on coordinates of value
1576  *         location point. For example, if \a this field is on cells, the function
1577  *         is applied to cell barycenters.
1578  *  \throw If the mesh is not set.
1579  *  \throw If the spatial discretization of \a this field is NULL.
1580  *  \throw If computing \a func fails.
1581  *
1582  *  \if ENABLE_EXAMPLES
1583  *  \ref cpp_mcfielddouble_fillFromAnalytic2 "Here is a C++ example".<br>
1584  *  \ref  py_mcfielddouble_fillFromAnalytic2 "Here is a Python example".
1585  *  \endif
1586  */
1587 void MEDCouplingFieldDouble::fillFromAnalytic2(int nbOfComp, const std::string& func)
1588 {
1589   if(!_mesh)
1590     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic2 : no mesh defined !");
1591   if(!((const MEDCouplingFieldDiscretization *)_type))
1592     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform fillFromAnalytic2 !");
1593   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
1594   _time_discr->fillFromAnalytic2(loc,nbOfComp,func);
1595 }
1596
1597 /*!
1598  * Creates data array(s) of \a this field by using a function for value generation.<br>
1599  * The function is applied to coordinates of value location points. For example, if
1600  * \a this field is on cells, the function is applied to cell barycenters.<br>
1601  * This method differs from
1602  * \ref ParaMEDMEM::MEDCouplingFieldDouble::fillFromAnalytic(int nbOfComp, const std::string& func) "fillFromAnalytic()"
1603  * by the way how variable
1604  * names, used in the function, are associated with components of coordinates of field
1605  * location points; here, a component index of a variable is defined by a
1606  * rank of the variable within the input array \a varsOrder.<br>
1607  * For more info on supported expressions that can be used in the function, see \ref
1608  * MEDCouplingArrayApplyFuncExpr.
1609  * In a general case, a value resulting from the function evaluation is assigned to all
1610  * components of a field value. But there is a possibility to have its own expression for
1611  * each component within one function. For this purpose, there are predefined variable
1612  * names (IVec, JVec, KVec, LVec etc) each dedicated to a certain component (IVec, to
1613  * the component #0 etc). A factor of such a variable is added to the
1614  * corresponding component only.<br>
1615  * For example, \a nbOfComp == 4, names of
1616  * spatial components are given in \a varsOrder: ["x", "y","z"], coordinates of a
1617  * 3D point are (1.,3.,7.), then
1618  *   - "2*x + z"               produces (9.,9.,9.,9.)
1619  *   - "2*x*IVec + (x+z)*LVec" produces (2.,0.,0.,8.)
1620  *   - "2*y*IVec + z*KVec + x" produces (7.,1.,1.,8.)
1621  *
1622  *  \param [in] nbOfComp - the number of components for \a this field to have.
1623  *  \param [in] func - the function used to compute values of \a this field.
1624  *         This function is used to compute a field value basing on coordinates of value
1625  *         location point. For example, if \a this field is on cells, the function
1626  *         is applied to cell barycenters.
1627  *  \throw If the mesh is not set.
1628  *  \throw If the spatial discretization of \a this field is NULL.
1629  *  \throw If computing \a func fails.
1630  *
1631  *  \if ENABLE_EXAMPLES
1632  *  \ref cpp_mcfielddouble_fillFromAnalytic3 "Here is a C++ example".<br>
1633  *  \ref  py_mcfielddouble_fillFromAnalytic3 "Here is a Python example".
1634  *  \endif
1635  */
1636 void MEDCouplingFieldDouble::fillFromAnalytic3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func)
1637 {
1638   if(!_mesh)
1639     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::fillFromAnalytic2 : no mesh defined !");
1640   if(!((const MEDCouplingFieldDiscretization *)_type))
1641     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform fillFromAnalytic3 !");
1642   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> loc=_type->getLocalizationOfDiscValues(_mesh);
1643   _time_discr->fillFromAnalytic3(loc,nbOfComp,varsOrder,func);
1644 }
1645
1646 /*!
1647  * Modifies values of \a this field by applying a C function to each tuple of all
1648  * data arrays.
1649  *  \param [in] nbOfComp - the number of components for \a this field to have.
1650  *  \param [in] func - the function used to compute values of \a this field.
1651  *         This function is to compute a field value basing on a current field value.
1652  *  \throw If \a func returns \c false.
1653  *
1654  *  \if ENABLE_EXAMPLES
1655  *  \ref cpp_mcfielddouble_applyFunc_c_func "Here is a C++ example".
1656  *  \endif
1657  */
1658 void MEDCouplingFieldDouble::applyFunc(int nbOfComp, FunctionToEvaluate func)
1659 {
1660   _time_discr->applyFunc(nbOfComp,func);
1661 }
1662
1663 /*!
1664  * Fill \a this field with a given value.<br>
1665  * This method is a specialization of other overloaded methods. When \a nbOfComp == 1
1666  * this method is equivalent to ParaMEDMEM::MEDCouplingFieldDouble::operator=().
1667  *  \param [in] nbOfComp - the number of components for \a this field to have.
1668  *  \param [in] val - the value to assign to every atomic value of \a this field.
1669  *  \throw If the spatial discretization of \a this field is NULL.
1670  *  \throw If the mesh is not set.
1671  *
1672  *  \if ENABLE_EXAMPLES
1673  *  \ref cpp_mcfielddouble_applyFunc_val "Here is a C++ example".<br>
1674  *  \ref  py_mcfielddouble_applyFunc_val "Here is a Python example".
1675  *  \endif
1676  */
1677 void MEDCouplingFieldDouble::applyFunc(int nbOfComp, double val)
1678 {
1679   if(!_mesh)
1680     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::applyFunc : no mesh defined !");
1681   if(!((const MEDCouplingFieldDiscretization *)_type))
1682     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform applyFunc !");
1683   int nbOfTuple=_type->getNumberOfTuples(_mesh);
1684   _time_discr->setUniformValue(nbOfTuple,nbOfComp,val);
1685 }
1686
1687 /*!
1688  * Modifies values of \a this field by applying a function to each tuple of all
1689  * data arrays.
1690  * For more info on supported expressions that can be used in the function, see \ref
1691  * MEDCouplingArrayApplyFuncExpr. <br>
1692  * The function can include arbitrary named variables
1693  * (e.g. "x","y" or "va44") to refer to components of a field value. Names of
1694  * variables are sorted in \b alphabetical \b order to associate a variable name with a
1695  * component. For example, in the expression "2*x+z", "x" stands for the component #0
1696  * and "z" stands for the component #1 (\b not #2)!<br>
1697  * In a general case, a value resulting from the function evaluation is assigned to all
1698  * components of a field value. But there is a possibility to have its own expression for
1699  * each component within one function. For this purpose, there are predefined variable
1700  * names (IVec, JVec, KVec, LVec etc) each dedicated to a certain component (IVec, to
1701  * the component #0 etc). A factor of such a variable is added to the
1702  * corresponding component only.<br>
1703  * For example, \a nbOfComp == 4, components of a field value are (1.,3.,7.), then
1704  *   - "2*x + z"               produces (5.,5.,5.,5.)
1705  *   - "2*x + 0*y + z"         produces (9.,9.,9.,9.)
1706  *   - "2*x*IVec + (x+z)*LVec" produces (2.,0.,0.,4.)
1707  *   - "2*y*IVec + z*KVec + x" produces (7.,1.,1.,4.)
1708  *
1709  *  \param [in] nbOfComp - the number of components for \a this field to have.
1710  *  \param [in] func - the function used to compute values of \a this field.
1711  *         This function is to compute a field value basing on a current field value.
1712  *  \throw If computing \a func fails.
1713  *
1714  *  \if ENABLE_EXAMPLES
1715  *  \ref cpp_mcfielddouble_applyFunc "Here is a C++ example".<br>
1716  *  \ref  py_mcfielddouble_applyFunc "Here is a Python example".
1717  *  \endif
1718  */
1719 void MEDCouplingFieldDouble::applyFunc(int nbOfComp, const std::string& func)
1720 {
1721   _time_discr->applyFunc(nbOfComp,func);
1722 }
1723
1724
1725 /*!
1726  * Modifies values of \a this field by applying a function to each tuple of all
1727  * data arrays.
1728  * For more info on supported expressions that can be used in the function, see \ref
1729  * MEDCouplingArrayApplyFuncExpr. <br>
1730  * This method differs from
1731  * \ref ParaMEDMEM::MEDCouplingFieldDouble::applyFunc(int nbOfComp, const std::string& func) "applyFunc()"
1732  * by the way how variable
1733  * names, used in the function, are associated with components of field values;
1734  * here, a variable name corresponding to a component is retrieved from
1735  * component information of an array (where it is set via
1736  * DataArrayDouble::setInfoOnComponent()).<br>
1737  * In a general case, a value resulting from the function evaluation is assigned to all
1738  * components of a field value. But there is a possibility to have its own expression for
1739  * each component within one function. For this purpose, there are predefined variable
1740  * names (IVec, JVec, KVec, LVec etc) each dedicated to a certain component (IVec, to
1741  * the component #0 etc). A factor of such a variable is added to the
1742  * corresponding component only.<br>
1743  * For example, \a nbOfComp == 4, components of a field value are (1.,3.,7.), then
1744  *   - "2*x + z"               produces (5.,5.,5.,5.)
1745  *   - "2*x + 0*y + z"         produces (9.,9.,9.,9.)
1746  *   - "2*x*IVec + (x+z)*LVec" produces (2.,0.,0.,4.)
1747  *   - "2*y*IVec + z*KVec + x" produces (7.,1.,1.,4.)
1748  *
1749  *  \param [in] nbOfComp - the number of components for \a this field to have.
1750  *  \param [in] func - the function used to compute values of \a this field.
1751  *         This function is to compute a new field value basing on a current field value.
1752  *  \throw If computing \a func fails.
1753  *
1754  *  \if ENABLE_EXAMPLES
1755  *  \ref cpp_mcfielddouble_applyFunc2 "Here is a C++ example".<br>
1756  *  \ref  py_mcfielddouble_applyFunc2 "Here is a Python example".
1757  *  \endif
1758  */
1759 void MEDCouplingFieldDouble::applyFunc2(int nbOfComp, const std::string& func)
1760 {
1761   _time_discr->applyFunc2(nbOfComp,func);
1762 }
1763
1764 /*!
1765  * Modifies values of \a this field by applying a function to each tuple of all
1766  * data arrays.
1767  * This method differs from
1768  * \ref ParaMEDMEM::MEDCouplingFieldDouble::applyFunc(int nbOfComp, const std::string& func) "applyFunc()"
1769  * by the way how variable
1770  * names, used in the function, are associated with components of field values;
1771  * here, a component index of a variable is defined by a
1772  * rank of the variable within the input array \a varsOrder.<br>
1773  * For more info on supported expressions that can be used in the function, see \ref
1774  * MEDCouplingArrayApplyFuncExpr.
1775  * In a general case, a value resulting from the function evaluation is assigned to all
1776  * components of a field value. But there is a possibility to have its own expression for
1777  * each component within one function. For this purpose, there are predefined variable
1778  * names (IVec, JVec, KVec, LVec etc) each dedicated to a certain component (IVec, to
1779  * the component #0 etc). A factor of such a variable is added to the
1780  * corresponding component only.<br>
1781  * For example, \a nbOfComp == 4, names of
1782  * components are given in \a varsOrder: ["x", "y","z"], components of a
1783  * 3D vector are (1.,3.,7.), then
1784  *   - "2*x + z"               produces (9.,9.,9.,9.)
1785  *   - "2*x*IVec + (x+z)*LVec" produces (2.,0.,0.,8.)
1786  *   - "2*y*IVec + z*KVec + x" produces (7.,1.,1.,8.)
1787  *
1788  *  \param [in] nbOfComp - the number of components for \a this field to have.
1789  *  \param [in] func - the function used to compute values of \a this field.
1790  *         This function is to compute a new field value basing on a current field value.
1791  *  \throw If computing \a func fails.
1792  *
1793  *  \if ENABLE_EXAMPLES
1794  *  \ref cpp_mcfielddouble_applyFunc3 "Here is a C++ example".<br>
1795  *  \ref  py_mcfielddouble_applyFunc3 "Here is a Python example".
1796  *  \endif
1797  */
1798 void MEDCouplingFieldDouble::applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func)
1799 {
1800   _time_discr->applyFunc3(nbOfComp,varsOrder,func);
1801 }
1802
1803 /*!
1804  * Modifies values of \a this field by applying a function to each atomic value of all
1805  * data arrays. The function computes a new single value basing on an old single value.
1806  * For more info on supported expressions that can be used in the function, see \ref
1807  * MEDCouplingArrayApplyFuncExpr. <br>
1808  * The function can include **only one** arbitrary named variable
1809  * (e.g. "x","y" or "va44") to refer to a field atomic value. <br>
1810  * In a general case, a value resulting from the function evaluation is assigned to 
1811  * a single field value. But there is a possibility to have its own expression for
1812  * each component within one function. For this purpose, there are predefined variable
1813  * names (IVec, JVec, KVec, LVec etc) each dedicated to a certain component (IVec, to
1814  * the component #0 etc). A factor of such a variable is added to the
1815  * corresponding component only.<br>
1816  * For example, components of a field value are (1.,3.,7.), then
1817  *   - "2*x - 1"               produces (1.,5.,13.)
1818  *   - "2*x*IVec + (x+3)*KVec" produces (2.,0.,10.)
1819  *   - "2*x*IVec + (x+3)*KVec + 1" produces (3.,1.,11.)
1820  *
1821  *  \param [in] func - the function used to compute values of \a this field.
1822  *         This function is to compute a field value basing on a current field value.
1823  *  \throw If computing \a func fails.
1824  *
1825  *  \if ENABLE_EXAMPLES
1826  *  \ref cpp_mcfielddouble_applyFunc_same_nb_comp "Here is a C++ example".<br>
1827  *  \ref  py_mcfielddouble_applyFunc_same_nb_comp "Here is a Python example".
1828  *  \endif
1829  */
1830 void MEDCouplingFieldDouble::applyFunc(const std::string& func)
1831 {
1832   _time_discr->applyFunc(func);
1833 }
1834
1835 /*!
1836  * Applyies the function specified by the string repr 'func' on each tuples on all arrays contained in _time_discr.
1837  * The field will contain exactly the same number of components after the call.
1838  * Use is not warranted for the moment !
1839  */
1840 void MEDCouplingFieldDouble::applyFuncFast32(const std::string& func)
1841 {
1842   _time_discr->applyFuncFast32(func);
1843 }
1844
1845 /*!
1846  * Applyies the function specified by the string repr 'func' on each tuples on all arrays contained in _time_discr.
1847  * The field will contain exactly the same number of components after the call.
1848  * Use is not warranted for the moment !
1849  */
1850 void MEDCouplingFieldDouble::applyFuncFast64(const std::string& func)
1851 {
1852   _time_discr->applyFuncFast64(func);
1853 }
1854
1855 /*!
1856  * Returns number of components in the data array. For more info on the data arrays,
1857  * see \ref MEDCouplingArrayPage.
1858  *  \return int - the number of components in the data array.
1859  *  \throw If the data array is not set.
1860  */
1861 int MEDCouplingFieldDouble::getNumberOfComponents() const
1862 {
1863   if(getArray()==0)
1864     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getNumberOfComponents : No array specified !");
1865   return getArray()->getNumberOfComponents();
1866 }
1867
1868 /*!
1869  * Returns number of tuples in \a this field, that depends on 
1870  * - the number of entities in the underlying mesh
1871  * - \ref MEDCouplingSpatialDisc "spatial discretization" of \a this field (e.g. number
1872  * of Gauss points if \a this->getTypeOfField() == 
1873  * \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT").
1874  *
1875  * The returned value does **not depend** on the number of tuples in the data array
1876  * (which has to be equal to the returned value), \b contrary to
1877  * getNumberOfComponents() and getNumberOfValues() that retrieve information from the
1878  * data array.
1879  * \warning No checkCoherency() is done here.
1880  * For more info on the data arrays, see \ref MEDCouplingArrayPage.
1881  *  \return int - the number of tuples.
1882  *  \throw If the mesh is not set.
1883  *  \throw If the spatial discretization of \a this field is NULL.
1884  *  \throw If the spatial discretization is not fully defined.
1885  */
1886 int MEDCouplingFieldDouble::getNumberOfTuples() const
1887 {
1888   if(!_mesh)
1889     throw INTERP_KERNEL::Exception("Impossible to retrieve number of tuples because no mesh specified !");
1890   if(!((const MEDCouplingFieldDiscretization *)_type))
1891     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getNumberOfTuples !");
1892   return _type->getNumberOfTuples(_mesh);
1893 }
1894
1895 /*!
1896  * Returns number of atomic double values in the data array of \a this field.
1897  * For more info on the data arrays, see \ref MEDCouplingArrayPage.
1898  *  \return int - (number of tuples) * (number of components) of the
1899  *  data array.
1900  *  \throw If the data array is not set.
1901  */
1902 int MEDCouplingFieldDouble::getNumberOfValues() const
1903 {
1904   if(getArray()==0)
1905     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getNumberOfValues : No array specified !");
1906   return getArray()->getNbOfElems();
1907 }
1908
1909 /*!
1910  * Sets own modification time by the most recently modified element of data (the mesh,
1911  * the data array etc). For more info, see \ref MEDCouplingTimeLabelPage.
1912  */
1913 void MEDCouplingFieldDouble::updateTime() const
1914 {
1915   MEDCouplingField::updateTime();
1916   updateTimeWith(*_time_discr);
1917 }
1918
1919 std::size_t MEDCouplingFieldDouble::getHeapMemorySizeWithoutChildren() const
1920 {
1921   return MEDCouplingField::getHeapMemorySizeWithoutChildren();
1922 }
1923
1924 std::vector<const BigMemoryObject *> MEDCouplingFieldDouble::getDirectChildren() const
1925 {
1926   std::vector<const BigMemoryObject *> ret(MEDCouplingField::getDirectChildren());
1927   if(_time_discr)
1928     {
1929       std::vector<const BigMemoryObject *> ret2(_time_discr->getDirectChildren());
1930       ret.insert(ret.end(),ret2.begin(),ret2.end());
1931     }
1932   return ret;
1933 }
1934
1935 /*!
1936  * Sets \ref NatureOfField.
1937  *  \param [in] nat - an item of enum ParaMEDMEM::NatureOfField.
1938  */
1939 void MEDCouplingFieldDouble::setNature(NatureOfField nat)
1940 {
1941   MEDCouplingField::setNature(nat);
1942   if(_type)
1943     _type->checkCompatibilityWithNature(nat);
1944 }
1945
1946 /*!
1947  * This method synchronizes time information (time, iteration, order, time unit) regarding the information in \c this->_mesh.
1948  * \throw If no mesh is set in this. Or if \a this is not compatible with time setting (typically NO_TIME)
1949  */
1950 void MEDCouplingFieldDouble::synchronizeTimeWithMesh()
1951 {
1952   if(!_mesh)
1953     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::synchronizeTimeWithMesh : no mesh set in this !");
1954   int it=-1,ordr=-1;
1955   double val=_mesh->getTime(it,ordr);
1956   std::string timeUnit(_mesh->getTimeUnit());
1957   setTime(val,it,ordr);
1958   setTimeUnit(timeUnit);
1959 }
1960
1961 /*!
1962  * Returns a value of \a this field of type either
1963  * \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT" or
1964  * \ref ParaMEDMEM::ON_GAUSS_NE "ON_GAUSS_NE".
1965  *  \param [in] cellId - an id of cell of interest.
1966  *  \param [in] nodeIdInCell - a node index within the cell.
1967  *  \param [in] compoId - an index of component.
1968  *  \return double - the field value corresponding to the specified parameters.
1969  *  \throw If the data array is not set.
1970  *  \throw If the mesh is not set.
1971  *  \throw If the spatial discretization of \a this field is NULL.
1972  *  \throw If \a this field if of type other than 
1973  *         \ref ParaMEDMEM::ON_GAUSS_PT "ON_GAUSS_PT" or
1974  *         \ref ParaMEDMEM::ON_GAUSS_NE "ON_GAUSS_NE".
1975  */
1976 double MEDCouplingFieldDouble::getIJK(int cellId, int nodeIdInCell, int compoId) const
1977 {
1978   if(!((const MEDCouplingFieldDiscretization *)_type))
1979     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getIJK !");
1980   return _type->getIJK(_mesh,getArray(),cellId,nodeIdInCell,compoId);
1981 }
1982
1983 /*!
1984  * Sets the data array. 
1985  *  \param [in] array - the data array holding values of \a this field. It's size
1986  *         should correspond to the mesh and
1987  *         \ref MEDCouplingSpatialDisc "spatial discretization" of \a this field
1988  *         (see getNumberOfTuples()), but this size is not checked here.
1989  */
1990 void MEDCouplingFieldDouble::setArray(DataArrayDouble *array)
1991 {
1992   _time_discr->setArray(array,this);
1993 }
1994
1995 /*!
1996  * Sets the data array holding values corresponding to an end of a time interval
1997  * for which \a this field is defined.
1998  *  \param [in] array - the data array holding values of \a this field. It's size
1999  *         should correspond to the mesh and
2000  *         \ref MEDCouplingSpatialDisc "spatial discretization" of \a this field
2001  *         (see getNumberOfTuples()), but this size is not checked here.
2002  */
2003 void MEDCouplingFieldDouble::setEndArray(DataArrayDouble *array)
2004 {
2005   _time_discr->setEndArray(array,this);
2006 }
2007
2008 /*!
2009  * Sets all data arrays needed to define the field values.
2010  *  \param [in] arrs - a vector of DataArrayDouble's holding values of \a this
2011  *         field. Size of each array should correspond to the mesh and
2012  *         \ref MEDCouplingSpatialDisc "spatial discretization" of \a this field
2013  *         (see getNumberOfTuples()), but this size is not checked here.
2014  *  \throw If number of arrays in \a arrs does not correspond to type of
2015  *         \ref MEDCouplingTemporalDisc "temporal discretization" of \a this field.
2016  */
2017 void MEDCouplingFieldDouble::setArrays(const std::vector<DataArrayDouble *>& arrs)
2018 {
2019   _time_discr->setArrays(arrs,this);
2020 }
2021
2022 void MEDCouplingFieldDouble::getTinySerializationStrInformation(std::vector<std::string>& tinyInfo) const
2023 {
2024   tinyInfo.clear();
2025   _time_discr->getTinySerializationStrInformation(tinyInfo);
2026   tinyInfo.push_back(_name);
2027   tinyInfo.push_back(_desc);
2028   tinyInfo.push_back(getTimeUnit());
2029 }
2030
2031 /*!
2032  * This method retrieves some critical values to resize and prepare remote instance.
2033  * The first two elements returned in tinyInfo correspond to the parameters to give in constructor.
2034  * @param tinyInfo out parameter resized correctly after the call. The length of this vector is tiny.
2035  */
2036 void MEDCouplingFieldDouble::getTinySerializationIntInformation(std::vector<int>& tinyInfo) const
2037 {
2038   if(!((const MEDCouplingFieldDiscretization *)_type))
2039     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getTinySerializationIntInformation !");
2040   tinyInfo.clear();
2041   tinyInfo.push_back((int)_type->getEnum());
2042   tinyInfo.push_back((int)_time_discr->getEnum());
2043   tinyInfo.push_back((int)_nature);
2044   _time_discr->getTinySerializationIntInformation(tinyInfo);
2045   std::vector<int> tinyInfo2;
2046   _type->getTinySerializationIntInformation(tinyInfo2);
2047   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
2048   tinyInfo.push_back((int)tinyInfo2.size());
2049 }
2050
2051 /*!
2052  * This method retrieves some critical values to resize and prepare remote instance.
2053  * @param tinyInfo out parameter resized correctly after the call. The length of this vector is tiny.
2054  */
2055 void MEDCouplingFieldDouble::getTinySerializationDbleInformation(std::vector<double>& tinyInfo) const
2056 {
2057   if(!((const MEDCouplingFieldDiscretization *)_type))
2058     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform getTinySerializationDbleInformation !");
2059   tinyInfo.clear();
2060   _time_discr->getTinySerializationDbleInformation(tinyInfo);
2061   std::vector<double> tinyInfo2;
2062   _type->getTinySerializationDbleInformation(tinyInfo2);
2063   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
2064   tinyInfo.push_back((int)tinyInfo2.size());//very bad, lack of time to improve it
2065 }
2066
2067 /*!
2068  * This method has to be called to the new instance filled by CORBA, MPI, File...
2069  * @param tinyInfoI is the value retrieves from distant result of getTinySerializationIntInformation on source instance to be copied.
2070  * @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.
2071  * @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.
2072  *               No decrRef must be applied to every instances in returned vector.
2073  */
2074 void MEDCouplingFieldDouble::resizeForUnserialization(const std::vector<int>& tinyInfoI, DataArrayInt *&dataInt, std::vector<DataArrayDouble *>& arrays)
2075 {
2076   if(!((const MEDCouplingFieldDiscretization *)_type))
2077     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform resizeForUnserialization !");
2078   dataInt=0;
2079   std::vector<int> tinyInfoITmp(tinyInfoI);
2080   int sz=tinyInfoITmp.back();
2081   tinyInfoITmp.pop_back();
2082   std::vector<int> tinyInfoITmp2(tinyInfoITmp.begin(),tinyInfoITmp.end()-sz);
2083   std::vector<int> tinyInfoI2(tinyInfoITmp2.begin()+3,tinyInfoITmp2.end());
2084   _time_discr->resizeForUnserialization(tinyInfoI2,arrays);
2085   std::vector<int> tinyInfoITmp3(tinyInfoITmp.end()-sz,tinyInfoITmp.end());
2086   _type->resizeForUnserialization(tinyInfoITmp3,dataInt);
2087 }
2088
2089 void MEDCouplingFieldDouble::finishUnserialization(const std::vector<int>& tinyInfoI, const std::vector<double>& tinyInfoD, const std::vector<std::string>& tinyInfoS)
2090 {
2091   if(!((const MEDCouplingFieldDiscretization *)_type))
2092     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform finishUnserialization !");
2093   std::vector<int> tinyInfoI2(tinyInfoI.begin()+3,tinyInfoI.end());
2094   //
2095   std::vector<double> tmp(tinyInfoD);
2096   int sz=(int)tinyInfoD.back();//very bad, lack of time to improve it
2097   tmp.pop_back();
2098   std::vector<double> tmp1(tmp.begin(),tmp.end()-sz);
2099   std::vector<double> tmp2(tmp.end()-sz,tmp.end());
2100   //
2101   _time_discr->finishUnserialization(tinyInfoI2,tmp1,tinyInfoS);
2102   _nature=(NatureOfField)tinyInfoI[2];
2103   _type->finishUnserialization(tmp2);
2104   int nbOfElemS=(int)tinyInfoS.size();
2105   _name=tinyInfoS[nbOfElemS-3];
2106   _desc=tinyInfoS[nbOfElemS-2];
2107   setTimeUnit(tinyInfoS[nbOfElemS-1]);
2108 }
2109
2110 /*!
2111  * Contrary to MEDCouplingPointSet class the returned arrays are \b not the responsabilities of the caller.
2112  * The values returned must be consulted only in readonly mode.
2113  */
2114 void MEDCouplingFieldDouble::serialize(DataArrayInt *&dataInt, std::vector<DataArrayDouble *>& arrays) const
2115 {
2116   if(!((const MEDCouplingFieldDiscretization *)_type))
2117     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform serialize !");
2118   _time_discr->getArrays(arrays);
2119   _type->getSerializationIntArray(dataInt);
2120 }
2121
2122 /*!
2123  * Tries to set an \a other mesh as the support of \a this field. An attempt fails, if
2124  * a current and the \a other meshes are different with use of specified equality
2125  * criteria, and then an exception is thrown.
2126  *  \param [in] other - the mesh to use as the field support if this mesh can be
2127  *         considered equal to the current mesh.
2128  *  \param [in] levOfCheck - defines equality criteria used for mesh comparison. For
2129  *         it's meaning explanation, see MEDCouplingMesh::checkGeoEquivalWith() which
2130  *         is used for mesh comparison.
2131  *  \param [in] precOnMesh - a precision used to compare nodes of the two meshes.
2132  *         It is used as \a prec parameter of MEDCouplingMesh::checkGeoEquivalWith().
2133  *  \param [in] eps - a precision used at node renumbering (if needed) to compare field
2134  *         values at merged nodes. If the values differ more than \a eps, an
2135  *         exception is thrown.
2136  *  \throw If the mesh is not set.
2137  *  \throw If \a other == NULL.
2138  *  \throw If any of the meshes is not well defined.
2139  *  \throw If the two meshes do not match.
2140  *  \throw If field values at merged nodes (if any) deffer more than \a eps.
2141  *
2142  *  \if ENABLE_EXAMPLES
2143  *  \ref cpp_mcfielddouble_changeUnderlyingMesh "Here is a C++ example".<br>
2144  *  \ref  py_mcfielddouble_changeUnderlyingMesh "Here is a Python example".
2145  *  \endif
2146  */
2147 void MEDCouplingFieldDouble::changeUnderlyingMesh(const MEDCouplingMesh *other, int levOfCheck, double precOnMesh, double eps)
2148 {
2149   if(_mesh==0 || other==0)
2150     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::changeUnderlyingMesh : is expected to operate on not null meshes !");
2151   DataArrayInt *cellCor=0,*nodeCor=0;
2152   other->checkGeoEquivalWith(_mesh,levOfCheck,precOnMesh,cellCor,nodeCor);
2153   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellCor2(cellCor),nodeCor2(nodeCor);
2154   if(cellCor)
2155     renumberCellsWithoutMesh(cellCor->getConstPointer(),false);
2156   if(nodeCor)
2157     renumberNodesWithoutMesh(nodeCor->getConstPointer(),nodeCor->getMaxValueInArray()+1,eps);
2158   setMesh(other);
2159 }
2160
2161 /*!
2162  * Subtracts another field from \a this one in case when the two fields have different
2163  * supporting meshes. The subtraction is performed provided that the tho meshes can be
2164  * considered equal with use of specified equality criteria, else an exception is thrown.
2165  * If the meshes match, the mesh of \a f is set to \a this field (\a this is permuted if 
2166  * necessary) and field values are subtracted. No interpolation is done here, only an
2167  * analysis of two underlying mesh is done to see if the meshes are geometrically
2168  * equivalent.<br>
2169  * The job of this method consists in calling
2170  * \a this->changeUnderlyingMesh() with \a f->getMesh() as the first parameter, and then
2171  * \a this -= \a f.<br>
2172  * This method requires that \a f and \a this are coherent (checkCoherency()) and that \a f
2173  * and \a this are coherent for a merge.<br>
2174  * "DM" in the method name stands for "different meshes".
2175  *  \param [in] f - the field to subtract from this.
2176  *  \param [in] levOfCheck - defines equality criteria used for mesh comparison. For
2177  *         it's meaning explanation, see MEDCouplingMesh::checkGeoEquivalWith() which
2178  *         is used for mesh comparison.
2179  *  \param [in] precOnMesh - a precision used to compare nodes of the two meshes.
2180  *         It is used as \a prec parameter of MEDCouplingMesh::checkGeoEquivalWith().
2181  *  \param [in] eps - a precision used at node renumbering (if needed) to compare field
2182  *         values at merged nodes. If the values differ more than \a eps, an
2183  *         exception is thrown.
2184  *  \throw If \a f == NULL.
2185  *  \throw If any of the meshes is not set or is not well defined.
2186  *  \throw If the two meshes do not match.
2187  *  \throw If the two fields are not coherent for merge.
2188  *  \throw If field values at merged nodes (if any) deffer more than \a eps.
2189  *
2190  *  \if ENABLE_EXAMPLES
2191  *  \ref cpp_mcfielddouble_substractInPlaceDM "Here is a C++ example".<br>
2192  *  \ref  py_mcfielddouble_substractInPlaceDM "Here is a Python example".
2193  *  \endif
2194  *  \sa changeUnderlyingMesh().
2195  */
2196 void MEDCouplingFieldDouble::substractInPlaceDM(const MEDCouplingFieldDouble *f, int levOfCheck, double precOnMesh, double eps)
2197 {
2198   checkCoherency();
2199   if(!f)
2200     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::substractInPlaceDM : input field is NULL !");
2201   f->checkCoherency();
2202   if(!areCompatibleForMerge(f))
2203     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::substractInPlaceDM : Fields are not compatible ; unable to apply mergeFields on them !");
2204   changeUnderlyingMesh(f->getMesh(),levOfCheck,precOnMesh,eps);
2205   operator-=(*f);
2206 }
2207
2208 /*!
2209  * Merges coincident nodes of the underlying mesh. If some nodes are coincident, the
2210  * underlying mesh is replaced by a new mesh instance where the coincident nodes are merged.
2211  *  \param [in] eps - a precision used to compare nodes of the two meshes.
2212  *  \param [in] epsOnVals - a precision used to compare field
2213  *         values at merged nodes. If the values differ more than \a epsOnVals, an
2214  *         exception is thrown.
2215  *  \return bool - \c true if some nodes have been merged and hence \a this field lies
2216  *         on another mesh.
2217  *  \throw If the mesh is of type not inheriting from MEDCouplingPointSet.
2218  *  \throw If the mesh is not well defined.
2219  *  \throw If the spatial discretization of \a this field is NULL.
2220  *  \throw If the data array is not set.
2221  *  \throw If field values at merged nodes (if any) deffer more than \a epsOnVals.
2222  */
2223 bool MEDCouplingFieldDouble::mergeNodes(double eps, double epsOnVals)
2224 {
2225   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
2226   if(!meshC)
2227     throw INTERP_KERNEL::Exception("Invalid support mesh to apply mergeNodes on it : must be a MEDCouplingPointSet one !");
2228   if(!((const MEDCouplingFieldDiscretization *)_type))
2229     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform mergeNodes !");
2230   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
2231   bool ret;
2232   int ret2;
2233   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->mergeNodes(eps,ret,ret2);
2234   if(!ret)//no nodes have been merged.
2235     return ret;
2236   std::vector<DataArrayDouble *> arrays;
2237   _time_discr->getArrays(arrays);
2238   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
2239     if(*iter)
2240       _type->renumberValuesOnNodes(epsOnVals,arr->getConstPointer(),meshC2->getNumberOfNodes(),*iter);
2241   setMesh(meshC2);
2242   return true;
2243 }
2244
2245 /*!
2246  * Merges coincident nodes of the underlying mesh. If some nodes are coincident, the
2247  * underlying mesh is replaced by a new mesh instance where the coincident nodes are
2248  * merged.<br>
2249  * In contrast to mergeNodes(), location of merged nodes is changed to be at their barycenter.
2250  *  \param [in] eps - a precision used to compare nodes of the two meshes.
2251  *  \param [in] epsOnVals - a precision used to compare field
2252  *         values at merged nodes. If the values differ more than \a epsOnVals, an
2253  *         exception is thrown.
2254  *  \return bool - \c true if some nodes have been merged and hence \a this field lies
2255  *         on another mesh.
2256  *  \throw If the mesh is of type not inheriting from MEDCouplingPointSet.
2257  *  \throw If the mesh is not well defined.
2258  *  \throw If the spatial discretization of \a this field is NULL.
2259  *  \throw If the data array is not set.
2260  *  \throw If field values at merged nodes (if any) deffer more than \a epsOnVals.
2261  */
2262 bool MEDCouplingFieldDouble::mergeNodes2(double eps, double epsOnVals)
2263 {
2264   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
2265   if(!meshC)
2266     throw INTERP_KERNEL::Exception("Invalid support mesh to apply mergeNodes on it : must be a MEDCouplingPointSet one !");
2267   if(!((const MEDCouplingFieldDiscretization *)_type))
2268     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform mergeNodes2 !");
2269   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
2270   bool ret;
2271   int ret2;
2272   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->mergeNodes2(eps,ret,ret2);
2273   if(!ret)//no nodes have been merged.
2274     return ret;
2275   std::vector<DataArrayDouble *> arrays;
2276   _time_discr->getArrays(arrays);
2277   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
2278     if(*iter)
2279       _type->renumberValuesOnNodes(epsOnVals,arr->getConstPointer(),meshC2->getNumberOfNodes(),*iter);
2280   setMesh(meshC2);
2281   return true;
2282 }
2283
2284 /*!
2285  * Removes from the underlying mesh nodes not used in any cell. If some nodes are
2286  * removed, the underlying mesh is replaced by a new mesh instance where the unused
2287  * nodes are removed.<br>
2288  *  \param [in] epsOnVals - a precision used to compare field
2289  *         values at merged nodes. If the values differ more than \a epsOnVals, an
2290  *         exception is thrown.
2291  *  \return bool - \c true if some nodes have been removed and hence \a this field lies
2292  *         on another mesh.
2293  *  \throw If the mesh is of type not inheriting from MEDCouplingPointSet.
2294  *  \throw If the mesh is not well defined.
2295  *  \throw If the spatial discretization of \a this field is NULL.
2296  *  \throw If the data array is not set.
2297  *  \throw If field values at merged nodes (if any) deffer more than \a epsOnVals.
2298  */
2299 bool MEDCouplingFieldDouble::zipCoords(double epsOnVals)
2300 {
2301   const MEDCouplingPointSet *meshC=dynamic_cast<const MEDCouplingPointSet *>(_mesh);
2302   if(!meshC)
2303     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::zipCoords : Invalid support mesh to apply zipCoords on it : must be a MEDCouplingPointSet one !");
2304   if(!((const MEDCouplingFieldDiscretization *)_type))
2305     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform zipCoords !");
2306   MEDCouplingAutoRefCountObjectPtr<MEDCouplingPointSet> meshC2((MEDCouplingPointSet *)meshC->deepCpy());
2307   int oldNbOfNodes=meshC2->getNumberOfNodes();
2308   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->zipCoordsTraducer();
2309   if(meshC2->getNumberOfNodes()!=oldNbOfNodes)
2310     {
2311       std::vector<DataArrayDouble *> arrays;
2312       _time_discr->getArrays(arrays);
2313       for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
2314         if(*iter)
2315           _type->renumberValuesOnNodes(epsOnVals,arr->getConstPointer(),meshC2->getNumberOfNodes(),*iter);
2316       setMesh(meshC2);
2317       return true;
2318     }
2319   return false;
2320 }
2321
2322 /*!
2323  * Removes duplicates of cells from the understanding mesh. If some cells are
2324  * removed, the underlying mesh is replaced by a new mesh instance where the cells
2325  * duplicates are removed.<br>
2326  *  \param [in] compType - specifies a cell comparison technique. Meaning of its
2327  *          valid values [0,1,2] is explained in the description of
2328  *          MEDCouplingPointSet::zipConnectivityTraducer() which is called by this method.
2329  *  \param [in] epsOnVals - a precision used to compare field
2330  *         values at merged cells. If the values differ more than \a epsOnVals, an
2331  *         exception is thrown.
2332  *  \return bool - \c true if some cells have been removed and hence \a this field lies
2333  *         on another mesh.
2334  *  \throw If the mesh is not an instance of MEDCouplingUMesh.
2335  *  \throw If the mesh is not well defined.
2336  *  \throw If the spatial discretization of \a this field is NULL.
2337  *  \throw If the data array is not set.
2338  *  \throw If field values at merged cells (if any) deffer more than \a epsOnVals.
2339  */
2340 bool MEDCouplingFieldDouble::zipConnectivity(int compType, double epsOnVals)
2341 {
2342   const MEDCouplingUMesh *meshC=dynamic_cast<const MEDCouplingUMesh *>(_mesh);
2343   if(!meshC)
2344     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::zipConnectivity : Invalid support mesh to apply zipCoords on it : must be a MEDCouplingPointSet one !");
2345   if(!((const MEDCouplingFieldDiscretization *)_type))
2346     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform zipConnectivity !");
2347   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> meshC2((MEDCouplingUMesh *)meshC->deepCpy());
2348   int oldNbOfCells=meshC2->getNumberOfCells();
2349   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->zipConnectivityTraducer(compType);
2350   if(meshC2->getNumberOfCells()!=oldNbOfCells)
2351     {
2352       std::vector<DataArrayDouble *> arrays;
2353       _time_discr->getArrays(arrays);
2354       for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
2355         if(*iter)
2356           _type->renumberValuesOnCells(epsOnVals,meshC,arr->getConstPointer(),meshC2->getNumberOfCells(),*iter);
2357       setMesh(meshC2);
2358       return true;
2359     }
2360   return false;
2361 }
2362
2363 /*!
2364  * This method calls MEDCouplingUMesh::buildSlice3D method. So this method makes the assumption that underlying mesh exists.
2365  * For the moment, this method is implemented for fields on cells.
2366  * 
2367  * \return a newly allocated field double containing the result that the user should deallocate.
2368  */
2369 MEDCouplingFieldDouble *MEDCouplingFieldDouble::extractSlice3D(const double *origin, const double *vec, double eps) const
2370 {
2371   const MEDCouplingMesh *mesh=getMesh();
2372   if(!mesh)
2373     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::extractSlice3D : underlying mesh is null !");
2374   if(getTypeOfField()!=ON_CELLS)
2375     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::extractSlice3D : only implemented for fields on cells !");
2376   const MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> umesh(mesh->buildUnstructured());
2377   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=clone(false);
2378   ret->setMesh(umesh);
2379   DataArrayInt *cellIds=0;
2380   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh2=umesh->buildSlice3D(origin,vec,eps,cellIds);
2381   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellIds2=cellIds;
2382   ret->setMesh(mesh2);
2383   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tupleIds=computeTupleIdsToSelectFromCellIds(cellIds->begin(),cellIds->end());
2384   std::vector<DataArrayDouble *> arrays;
2385   _time_discr->getArrays(arrays);
2386   int i=0;
2387   std::vector<DataArrayDouble *> newArr(arrays.size());
2388   std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > newArr2(arrays.size());
2389   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++,i++)
2390     {
2391       if(*iter)
2392         {
2393           newArr2[i]=(*iter)->selectByTupleIdSafe(cellIds->begin(),cellIds->end());
2394           newArr[i]=newArr2[i];
2395         }
2396     }
2397   ret->setArrays(newArr);
2398   return ret.retn();
2399 }
2400
2401 /*!
2402  * Divides every cell of the underlying mesh into simplices (triangles in 2D and
2403  * tetrahedra in 3D). If some cells are divided, the underlying mesh is replaced by a new
2404  * mesh instance containing the simplices.<br> 
2405  *  \param [in] policy - specifies a pattern used for splitting. For its description, see
2406  *          MEDCouplingUMesh::simplexize().
2407  *  \return bool - \c true if some cells have been divided and hence \a this field lies
2408  *         on another mesh.
2409  *  \throw If \a policy has an invalid value. For valid values, see the description of 
2410  *         MEDCouplingUMesh::simplexize().
2411  *  \throw If MEDCouplingMesh::simplexize() is not applicable to the underlying mesh.
2412  *  \throw If the mesh is not well defined.
2413  *  \throw If the spatial discretization of \a this field is NULL.
2414  *  \throw If the data array is not set.
2415  */
2416 bool MEDCouplingFieldDouble::simplexize(int policy)
2417 {
2418   if(!_mesh)
2419     throw INTERP_KERNEL::Exception("No underlying mesh on this field to perform simplexize !");
2420   if(!((const MEDCouplingFieldDiscretization *)_type))
2421     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform simplexize !");
2422   int oldNbOfCells=_mesh->getNumberOfCells();
2423   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> meshC2(_mesh->deepCpy());
2424   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr=meshC2->simplexize(policy);
2425   int newNbOfCells=meshC2->getNumberOfCells();
2426   if(oldNbOfCells==newNbOfCells)
2427     return false;
2428   std::vector<DataArrayDouble *> arrays;
2429   _time_discr->getArrays(arrays);
2430   for(std::vector<DataArrayDouble *>::const_iterator iter=arrays.begin();iter!=arrays.end();iter++)
2431     if(*iter)
2432       _type->renumberValuesOnCellsR(_mesh,arr->getConstPointer(),arr->getNbOfElems(),*iter);
2433   setMesh(meshC2);
2434   return true;
2435 }
2436
2437 /*!
2438  * Creates a new MEDCouplingFieldDouble filled with the doubly contracted product of
2439  * every tensor of \a this 6-componental field.
2440  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, whose
2441  *          each tuple is calculated from the tuple <em>(t)</em> of \a this field as
2442  *          follows: \f$ t[0]^2+t[1]^2+t[2]^2+2*t[3]^2+2*t[4]^2+2*t[5]^2\f$. 
2443  *          This new field lies on the same mesh as \a this one. The caller is to delete
2444  *          this field using decrRef() as it is no more needed.
2445  *  \throw If \a this->getNumberOfComponents() != 6.
2446  *  \throw If the spatial discretization of \a this field is NULL.
2447  */
2448 MEDCouplingFieldDouble *MEDCouplingFieldDouble::doublyContractedProduct() const
2449 {
2450   if(!((const MEDCouplingFieldDiscretization *)_type))
2451     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform doublyContractedProduct !");
2452   MEDCouplingTimeDiscretization *td=_time_discr->doublyContractedProduct();
2453   td->copyTinyAttrFrom(*_time_discr);
2454   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2455   ret->setName("DoublyContractedProduct");
2456   ret->setMesh(getMesh());
2457   return ret.retn();
2458 }
2459
2460 /*!
2461  * Creates a new MEDCouplingFieldDouble filled with the determinant of a square
2462  * matrix defined by every tuple of \a this field, having either 4, 6 or 9 components.
2463  * The case of 6 components corresponds to that of the upper triangular matrix. 
2464  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, whose
2465  *          each tuple is the determinant of matrix of the corresponding tuple of \a this 
2466  *          field. This new field lies on the same mesh as \a this one. The caller is to 
2467  *          delete this field using decrRef() as it is no more needed.
2468  *  \throw If \a this->getNumberOfComponents() is not in [4,6,9].
2469  *  \throw If the spatial discretization of \a this field is NULL.
2470  */
2471 MEDCouplingFieldDouble *MEDCouplingFieldDouble::determinant() const
2472 {
2473   if(!((const MEDCouplingFieldDiscretization *)_type))
2474     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform determinant !");
2475   MEDCouplingTimeDiscretization *td=_time_discr->determinant();
2476   td->copyTinyAttrFrom(*_time_discr);
2477   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2478   ret->setName("Determinant");
2479   ret->setMesh(getMesh());
2480   return ret.retn();
2481 }
2482
2483
2484 /*!
2485  * Creates a new MEDCouplingFieldDouble with 3 components filled with 3 eigenvalues of
2486  * an upper triangular matrix defined by every tuple of \a this 6-componental field.
2487  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, 
2488  *          having 3 components, whose each tuple contains the eigenvalues of the matrix of
2489  *          corresponding tuple of \a this field. This new field lies on the same mesh as
2490  *          \a this one. The caller is to delete this field using decrRef() as it is no
2491  *          more needed.  
2492  *  \throw If \a this->getNumberOfComponents() != 6.
2493  *  \throw If the spatial discretization of \a this field is NULL.
2494  */
2495 MEDCouplingFieldDouble *MEDCouplingFieldDouble::eigenValues() const
2496 {
2497   if(!((const MEDCouplingFieldDiscretization *)_type))
2498     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform eigenValues !");
2499   MEDCouplingTimeDiscretization *td=_time_discr->eigenValues();
2500   td->copyTinyAttrFrom(*_time_discr);
2501   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2502   ret->setName("EigenValues");
2503   ret->setMesh(getMesh());
2504   return ret.retn();
2505 }
2506
2507 /*!
2508  * Creates a new MEDCouplingFieldDouble with 9 components filled with 3 eigenvectors of
2509  * an upper triangular matrix defined by every tuple of \a this 6-componental field.
2510  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, 
2511  *          having 9 components, whose each tuple contains the eigenvectors of the matrix of
2512  *          corresponding tuple of \a this field. This new field lies on the same mesh as
2513  *          \a this one. The caller is to delete this field using decrRef() as it is no
2514  *          more needed.  
2515  *  \throw If \a this->getNumberOfComponents() != 6.
2516  *  \throw If the spatial discretization of \a this field is NULL.
2517  */
2518 MEDCouplingFieldDouble *MEDCouplingFieldDouble::eigenVectors() const
2519 {
2520   if(!((const MEDCouplingFieldDiscretization *)_type))
2521     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform eigenVectors !");
2522   MEDCouplingTimeDiscretization *td=_time_discr->eigenVectors();
2523   td->copyTinyAttrFrom(*_time_discr);
2524   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2525   ret->setName("EigenVectors");
2526   ret->setMesh(getMesh());
2527   return ret.retn();
2528 }
2529
2530 /*!
2531  * Creates a new MEDCouplingFieldDouble filled with the inverse matrix of
2532  * a matrix defined by every tuple of \a this field having either 4, 6 or 9
2533  * components. The case of 6 components corresponds to that of the upper triangular
2534  * matrix.
2535  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, 
2536  *          having the same number of components as \a this one, whose each tuple
2537  *          contains the inverse matrix of the matrix of corresponding tuple of \a this
2538  *          field. This new field lies on the same mesh as \a this one. The caller is to
2539  *          delete this field using decrRef() as it is no more needed.  
2540  *  \throw If \a this->getNumberOfComponents() is not in [4,6,9].
2541  *  \throw If the spatial discretization of \a this field is NULL.
2542  */
2543 MEDCouplingFieldDouble *MEDCouplingFieldDouble::inverse() const
2544 {
2545   if(!((const MEDCouplingFieldDiscretization *)_type))
2546     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform inverse !");
2547   MEDCouplingTimeDiscretization *td=_time_discr->inverse();
2548   td->copyTinyAttrFrom(*_time_discr);
2549   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2550   ret->setName("Inversion");
2551   ret->setMesh(getMesh());
2552   return ret.retn();
2553 }
2554
2555 /*!
2556  * Creates a new MEDCouplingFieldDouble filled with the trace of
2557  * a matrix defined by every tuple of \a this field having either 4, 6 or 9
2558  * components. The case of 6 components corresponds to that of the upper triangular
2559  * matrix.
2560  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, 
2561  *          having 1 component, whose each tuple is the trace of the matrix of
2562  *          corresponding tuple of \a this field.
2563  *          This new field lies on the same mesh as \a this one. The caller is to
2564  *          delete this field using decrRef() as it is no more needed.  
2565  *  \throw If \a this->getNumberOfComponents() is not in [4,6,9].
2566  *  \throw If the spatial discretization of \a this field is NULL.
2567  */
2568 MEDCouplingFieldDouble *MEDCouplingFieldDouble::trace() const
2569 {
2570   if(!((const MEDCouplingFieldDiscretization *)_type))
2571     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform trace !");
2572   MEDCouplingTimeDiscretization *td=_time_discr->trace();
2573   td->copyTinyAttrFrom(*_time_discr);
2574   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2575   ret->setName("Trace");
2576   ret->setMesh(getMesh());
2577   return ret.retn();
2578 }
2579
2580 /*!
2581  * Creates a new MEDCouplingFieldDouble filled with the stress deviator tensor of
2582  * a stress tensor defined by every tuple of \a this 6-componental field.
2583  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, 
2584  *          having same number of components and tuples as \a this field,
2585  *          whose each tuple contains the stress deviator tensor of the stress tensor of
2586  *          corresponding tuple of \a this field. This new field lies on the same mesh as
2587  *          \a this one. The caller is to delete this field using decrRef() as it is no
2588  *          more needed.  
2589  *  \throw If \a this->getNumberOfComponents() != 6.
2590  *  \throw If the spatial discretization of \a this field is NULL.
2591  */
2592 MEDCouplingFieldDouble *MEDCouplingFieldDouble::deviator() const
2593 {
2594   if(!((const MEDCouplingFieldDiscretization *)_type))
2595     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform deviator !");
2596   MEDCouplingTimeDiscretization *td=_time_discr->deviator();
2597   td->copyTinyAttrFrom(*_time_discr);
2598   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2599   ret->setName("Deviator");
2600   ret->setMesh(getMesh());
2601   return ret.retn();
2602 }
2603
2604 /*!
2605  * Creates a new MEDCouplingFieldDouble filled with the magnitude of
2606  * every vector of \a this field.
2607  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble, 
2608  *          having one component, whose each tuple is the magnitude of the vector
2609  *          of corresponding tuple of \a this field. This new field lies on the
2610  *          same mesh as \a this one. The caller is to
2611  *          delete this field using decrRef() as it is no more needed.  
2612  *  \throw If the spatial discretization of \a this field is NULL.
2613  */
2614 MEDCouplingFieldDouble *MEDCouplingFieldDouble::magnitude() const
2615 {
2616   if(!((const MEDCouplingFieldDiscretization *)_type))
2617     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform magnitude !");
2618   MEDCouplingTimeDiscretization *td=_time_discr->magnitude();
2619   td->copyTinyAttrFrom(*_time_discr);
2620   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2621   ret->setName("Magnitude");
2622   ret->setMesh(getMesh());
2623   return ret.retn();
2624 }
2625
2626 /*!
2627  * Creates a new scalar MEDCouplingFieldDouble filled with the maximal value among
2628  * values of every tuple of \a this field.
2629  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
2630  *          This new field lies on the same mesh as \a this one. The caller is to
2631  *          delete this field using decrRef() as it is no more needed.  
2632  *  \throw If the spatial discretization of \a this field is NULL.
2633  */
2634 MEDCouplingFieldDouble *MEDCouplingFieldDouble::maxPerTuple() const
2635 {
2636   if(!((const MEDCouplingFieldDiscretization *)_type))
2637     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform maxPerTuple !");
2638   MEDCouplingTimeDiscretization *td=_time_discr->maxPerTuple();
2639   td->copyTinyAttrFrom(*_time_discr);
2640   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2641   std::ostringstream oss;
2642   oss << "Max_" << getName();
2643   ret->setName(oss.str());
2644   ret->setMesh(getMesh());
2645   return ret.retn();
2646 }
2647
2648 /*!
2649  * Changes number of components in \a this field. If \a newNbOfComp is less
2650  * than \a this->getNumberOfComponents() then each tuple
2651  * is truncated to have \a newNbOfComp components, keeping first components. If \a
2652  * newNbOfComp is more than \a this->getNumberOfComponents() then 
2653  * each tuple is populated with \a dftValue to have \a newNbOfComp components.  
2654  *  \param [in] newNbOfComp - number of components for the new field to have.
2655  *  \param [in] dftValue - value assigned to new values added to \a this field.
2656  *  \throw If \a this is not allocated.
2657  */
2658 void MEDCouplingFieldDouble::changeNbOfComponents(int newNbOfComp, double dftValue)
2659 {
2660   _time_discr->changeNbOfComponents(newNbOfComp,dftValue);
2661 }
2662
2663 /*!
2664  * Creates a new MEDCouplingFieldDouble composed of selected components of \a this field.
2665  * The new MEDCouplingFieldDouble has the same number of tuples but includes components
2666  * specified by \a compoIds parameter. So that getNbOfElems() of the result field
2667  * can be either less, same or more than \a this->getNumberOfValues().
2668  *  \param [in] compoIds - sequence of zero based indices of components to include
2669  *              into the new field.
2670  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble that the caller
2671  *          is to delete using decrRef() as it is no more needed.
2672  *  \throw If a component index (\a i) is not valid: 
2673  *         \a i < 0 || \a i >= \a this->getNumberOfComponents().
2674  */
2675 MEDCouplingFieldDouble *MEDCouplingFieldDouble::keepSelectedComponents(const std::vector<int>& compoIds) const
2676 {
2677   if(!((const MEDCouplingFieldDiscretization *)_type))
2678     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform keepSelectedComponents !");
2679   MEDCouplingTimeDiscretization *td=_time_discr->keepSelectedComponents(compoIds);
2680   td->copyTinyAttrFrom(*_time_discr);
2681   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2682   ret->setName(getName());
2683   ret->setMesh(getMesh());
2684   return ret.retn();
2685 }
2686
2687
2688 /*!
2689  * Copy all components in a specified order from another field.
2690  * The number of tuples in \a this and the other field can be different.
2691  *  \param [in] f - the field to copy data from.
2692  *  \param [in] compoIds - sequence of zero based indices of components, data of which is
2693  *              to be copied.
2694  *  \throw If the two fields have different number of data arrays.
2695  *  \throw If a data array is set in one of fields and is not set in the other.
2696  *  \throw If \a compoIds.size() != \a a->getNumberOfComponents().
2697  *  \throw If \a compoIds[i] < 0 or \a compoIds[i] > \a this->getNumberOfComponents().
2698  */
2699 void MEDCouplingFieldDouble::setSelectedComponents(const MEDCouplingFieldDouble *f, const std::vector<int>& compoIds)
2700 {
2701   _time_discr->setSelectedComponents(f->_time_discr,compoIds);
2702 }
2703
2704 /*!
2705  * Sorts value within every tuple of \a this field.
2706  *  \param [in] asc - if \a true, the values are sorted in ascending order, else,
2707  *              in descending order.
2708  *  \throw If a data array is not allocated.
2709  */
2710 void MEDCouplingFieldDouble::sortPerTuple(bool asc)
2711 {
2712   _time_discr->sortPerTuple(asc);
2713 }
2714
2715 /*!
2716  * Creates a new MEDCouplingFieldDouble by concatenating two given fields.
2717  * Values of
2718  * the first field precede values of the second field within the result field.
2719  *  \param [in] f1 - the first field.
2720  *  \param [in] f2 - the second field.
2721  *  \return MEDCouplingFieldDouble * - the result field. It is a new instance of
2722  *          MEDCouplingFieldDouble. The caller is to delete this mesh using decrRef() 
2723  *          as it is no more needed.
2724  *  \throw If the fields are not compatible for the merge.
2725  *  \throw If the spatial discretization of \a f1 is NULL.
2726  *  \throw If the time discretization of \a f1 is NULL.
2727  *
2728  *  \if ENABLE_EXAMPLES
2729  *  \ref cpp_mcfielddouble_MergeFields "Here is a C++ example".<br>
2730  *  \ref  py_mcfielddouble_MergeFields "Here is a Python example".
2731  *  \endif
2732  */
2733 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MergeFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
2734 {
2735   if(!f1->areCompatibleForMerge(f2))
2736     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MergeFields on them !");
2737   const MEDCouplingMesh *m1(f1->getMesh()),*m2(f2->getMesh());
2738   if(!f1->_time_discr)
2739     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::MergeFields : no time discr of f1 !");
2740   if(!f1->_type)
2741     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::MergeFields : no spatial discr of f1 !");
2742   MEDCouplingTimeDiscretization *td=f1->_time_discr->aggregate(f2->_time_discr);
2743   td->copyTinyAttrFrom(*f1->_time_discr);
2744   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
2745   ret->setName(f1->getName());
2746   ret->setDescription(f1->getDescription());
2747   if(m1)
2748     {
2749       MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> m=m1->mergeMyselfWith(m2);
2750       ret->setMesh(m);
2751     }
2752   return ret.retn();
2753 }
2754
2755 /*!
2756  * Creates a new MEDCouplingFieldDouble by concatenating all given fields.
2757  * Values of the *i*-th field precede values of the (*i*+1)-th field within the result.
2758  * If there is only one field in \a a, a deepCopy() (except time information of mesh and
2759  * field) of the field is returned. 
2760  * Generally speaking the first field in \a a is used to assign tiny attributes of the
2761  * returned field. 
2762  *  \param [in] a - a vector of fields (MEDCouplingFieldDouble) to concatenate.
2763  *  \return MEDCouplingFieldDouble * - the result field. It is a new instance of
2764  *          MEDCouplingFieldDouble. The caller is to delete this mesh using decrRef() 
2765  *          as it is no more needed.
2766  *  \throw If \a a is empty.
2767  *  \throw If the fields are not compatible for the merge.
2768  *
2769  *  \if ENABLE_EXAMPLES
2770  *  \ref cpp_mcfielddouble_MergeFields "Here is a C++ example".<br>
2771  *  \ref  py_mcfielddouble_MergeFields "Here is a Python example".
2772  *  \endif
2773  */
2774 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MergeFields(const std::vector<const MEDCouplingFieldDouble *>& a)
2775 {
2776   if(a.size()<1)
2777     throw INTERP_KERNEL::Exception("FieldDouble::MergeFields : size of array must be >= 1 !");
2778   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> > ms(a.size());
2779   std::vector< const MEDCouplingUMesh *> ms2(a.size());
2780   std::vector< const MEDCouplingTimeDiscretization *> tds(a.size());
2781   std::vector<const MEDCouplingFieldDouble *>::const_iterator it=a.begin();
2782   const MEDCouplingFieldDouble *ref=(*it++);
2783   if(!ref)
2784     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::MergeFields : presence of NULL instance in first place of input vector !");
2785   for(;it!=a.end();it++)
2786     if(!ref->areCompatibleForMerge(*it))
2787       throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MergeFields on them !");
2788   for(int i=0;i<(int)a.size();i++)
2789     {
2790       if(a[i]->getMesh())
2791         { ms[i]=a[i]->getMesh()->buildUnstructured(); ms2[i]=ms[i]; }
2792       else
2793         { ms[i]=0; ms2[i]=0; }
2794       tds[i]=a[i]->_time_discr;
2795     }
2796   MEDCouplingTimeDiscretization *td=tds[0]->aggregate(tds);
2797   td->copyTinyAttrFrom(*(a[0]->_time_discr));
2798   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(a[0]->getNature(),td,a[0]->_type->clone());
2799   ret->setName(a[0]->getName());
2800   ret->setDescription(a[0]->getDescription());
2801   if(ms2[0])
2802     {
2803       MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=MEDCouplingUMesh::MergeUMeshes(ms2);
2804       m->copyTinyInfoFrom(ms2[0]);
2805       ret->setMesh(m);
2806     }
2807   return ret.retn();
2808 }
2809
2810 /*!
2811  * Creates a new MEDCouplingFieldDouble by concatenating components of two given fields.
2812  * The number of components in the result field is a sum of the number of components of
2813  * given fields. The number of tuples in the result field is same as that of each of given
2814  * arrays.
2815  * Number of tuples in the given fields must be the same.
2816  *  \param [in] f1 - a field to include in the result field.
2817  *  \param [in] f2 - another field to include in the result field.
2818  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
2819  *          The caller is to delete this result field using decrRef() as it is no more
2820  *          needed.
2821  *  \throw If the fields are not compatible for a meld (areCompatibleForMeld()).
2822  *  \throw If any of data arrays is not allocated.
2823  *  \throw If \a f1->getNumberOfTuples() != \a f2->getNumberOfTuples()
2824  */
2825 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MeldFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
2826 {
2827   if(!f1->areCompatibleForMeld(f2))
2828     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MeldFields on them !");
2829   MEDCouplingTimeDiscretization *td=f1->_time_discr->meld(f2->_time_discr);
2830   td->copyTinyAttrFrom(*f1->_time_discr);
2831   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
2832   ret->setMesh(f1->getMesh());
2833   return ret.retn();
2834 }
2835
2836 /*!
2837  * Returns a new MEDCouplingFieldDouble containing a dot product of two given fields, 
2838  * so that the i-th tuple of the result field is a sum of products of j-th components of
2839  * i-th tuples of given fields (\f$ f_i = \sum_{j=1}^n f1_j * f2_j \f$). 
2840  * Number of tuples and components in the given fields must be the same.
2841  *  \param [in] f1 - a given field.
2842  *  \param [in] f2 - another given field.
2843  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
2844  *          The caller is to delete this result field using decrRef() as it is no more
2845  *          needed.
2846  *  \throw If either \a f1 or \a f2 is NULL.
2847  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
2848  *         differ not only in values.
2849  */
2850 MEDCouplingFieldDouble *MEDCouplingFieldDouble::DotFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
2851 {
2852   if(!f1)
2853     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::DotFields : input field is NULL !");
2854   if(!f1->areStrictlyCompatible(f2))
2855     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply DotFields on them !");
2856   MEDCouplingTimeDiscretization *td=f1->_time_discr->dot(f2->_time_discr);
2857   td->copyTinyAttrFrom(*f1->_time_discr);
2858   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
2859   ret->setMesh(f1->getMesh());
2860   return ret;
2861 }
2862
2863 /*!
2864  * Returns a new MEDCouplingFieldDouble containing a cross product of two given fields, 
2865  * so that
2866  * the i-th tuple of the result field is a 3D vector which is a cross
2867  * product of two vectors defined by the i-th tuples of given fields.
2868  * Number of tuples in the given fields must be the same.
2869  * Number of components in the given fields must be 3.
2870  *  \param [in] f1 - a given field.
2871  *  \param [in] f2 - another given field.
2872  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
2873  *          The caller is to delete this result field using decrRef() as it is no more
2874  *          needed.
2875  *  \throw If either \a f1 or \a f2 is NULL.
2876  *  \throw If \a f1->getNumberOfComponents() != 3
2877  *  \throw If \a f2->getNumberOfComponents() != 3
2878  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
2879  *         differ not only in values.
2880  */
2881 MEDCouplingFieldDouble *MEDCouplingFieldDouble::CrossProductFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
2882 {
2883   if(!f1)
2884     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::CrossProductFields : input field is NULL !");
2885   if(!f1->areStrictlyCompatible(f2))
2886     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply CrossProductFields on them !");
2887   MEDCouplingTimeDiscretization *td=f1->_time_discr->crossProduct(f2->_time_discr);
2888   td->copyTinyAttrFrom(*f1->_time_discr);
2889   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
2890   ret->setMesh(f1->getMesh());
2891   return ret.retn();
2892 }
2893
2894 /*!
2895  * Returns a new MEDCouplingFieldDouble containing maximal values of two given fields.
2896  * Number of tuples and components in the given fields must be the same.
2897  *  \param [in] f1 - a field to compare values with another one.
2898  *  \param [in] f2 - another field to compare values with the first one.
2899  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
2900  *          The caller is to delete this result field using decrRef() as it is no more
2901  *          needed.
2902  *  \throw If either \a f1 or \a f2 is NULL.
2903  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
2904  *         differ not only in values.
2905  *
2906  *  \if ENABLE_EXAMPLES
2907  *  \ref cpp_mcfielddouble_MaxFields "Here is a C++ example".<br>
2908  *  \ref  py_mcfielddouble_MaxFields "Here is a Python example".
2909  *  \endif
2910  */
2911 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MaxFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
2912 {
2913   if(!f1)
2914     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::MaxFields : input field is NULL !");
2915   if(!f1->areStrictlyCompatible(f2))
2916     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MaxFields on them !");
2917   MEDCouplingTimeDiscretization *td=f1->_time_discr->max(f2->_time_discr);
2918   td->copyTinyAttrFrom(*f1->_time_discr);
2919   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
2920   ret->setMesh(f1->getMesh());
2921   return ret.retn();
2922 }
2923
2924 /*!
2925  * Returns a new MEDCouplingFieldDouble containing minimal values of two given fields.
2926  * Number of tuples and components in the given fields must be the same.
2927  *  \param [in] f1 - a field to compare values with another one.
2928  *  \param [in] f2 - another field to compare values with the first one.
2929  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
2930  *          The caller is to delete this result field using decrRef() as it is no more
2931  *          needed.
2932  *  \throw If either \a f1 or \a f2 is NULL.
2933  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
2934  *         differ not only in values.
2935  *
2936  *  \if ENABLE_EXAMPLES
2937  *  \ref cpp_mcfielddouble_MaxFields "Here is a C++ example".<br>
2938  *  \ref  py_mcfielddouble_MaxFields "Here is a Python example".
2939  *  \endif
2940  */
2941 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MinFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
2942 {
2943   if(!f1)
2944     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::MinFields : input field is NULL !");
2945   if(!f1->areStrictlyCompatible(f2))
2946     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MinFields on them !");
2947   MEDCouplingTimeDiscretization *td=f1->_time_discr->min(f2->_time_discr);
2948   td->copyTinyAttrFrom(*f1->_time_discr);
2949   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
2950   ret->setMesh(f1->getMesh());
2951   return ret.retn();
2952 }
2953
2954 /*!
2955  * Returns a copy of \a this field in which sign of all values is reversed.
2956  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble
2957  *         containing the same number of tuples and components as \a this field. 
2958  *         The caller is to delete this result field using decrRef() as it is no more
2959  *         needed. 
2960  *  \throw If the spatial discretization of \a this field is NULL.
2961  *  \throw If a data array is not allocated.
2962  */
2963 MEDCouplingFieldDouble *MEDCouplingFieldDouble::negate() const
2964 {
2965   if(!((const MEDCouplingFieldDiscretization *)_type))
2966     throw INTERP_KERNEL::Exception("No spatial discretization underlying this field to perform negate !");
2967   MEDCouplingTimeDiscretization *td=_time_discr->negate();
2968   td->copyTinyAttrFrom(*_time_discr);
2969   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(getNature(),td,_type->clone());
2970   ret->setMesh(getMesh());
2971   return ret.retn();
2972 }
2973
2974 /*!
2975  * Returns a new MEDCouplingFieldDouble containing sum values of corresponding values of
2976  * two given fields ( _f_ [ i, j ] = _f1_ [ i, j ] + _f2_ [ i, j ] ).
2977  * Number of tuples and components in the given fields must be the same.
2978  *  \param [in] f1 - a field to sum up.
2979  *  \param [in] f2 - another field to sum up.
2980  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
2981  *          The caller is to delete this result field using decrRef() as it is no more
2982  *          needed.
2983  *  \throw If either \a f1 or \a f2 is NULL.
2984  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
2985  *         differ not only in values.
2986  */
2987 MEDCouplingFieldDouble *MEDCouplingFieldDouble::AddFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
2988 {
2989   if(!f1)
2990     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::AddFields : input field is NULL !");
2991   if(!f1->areStrictlyCompatible(f2))
2992     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply AddFields on them !");
2993   MEDCouplingTimeDiscretization *td=f1->_time_discr->add(f2->_time_discr);
2994   td->copyTinyAttrFrom(*f1->_time_discr);
2995   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
2996   ret->setMesh(f1->getMesh());
2997   return ret.retn();
2998 }
2999
3000 /*!
3001  * Adds values of another MEDCouplingFieldDouble to values of \a this one
3002  * ( _this_ [ i, j ] += _other_ [ i, j ] ) using DataArrayDouble::addEqual().
3003  * The two fields must have same number of tuples, components and same underlying mesh.
3004  *  \param [in] other - the field to add to \a this one.
3005  *  \return const MEDCouplingFieldDouble & - a reference to \a this field.
3006  *  \throw If \a other is NULL.
3007  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
3008  *         differ not only in values.
3009  */
3010 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator+=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
3011 {
3012   if(!areStrictlyCompatible(&other))
3013     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply += on them !");
3014   _time_discr->addEqual(other._time_discr);
3015   return *this;
3016 }
3017
3018 /*!
3019  * Returns a new MEDCouplingFieldDouble containing subtraction of corresponding values of
3020  * two given fields ( _f_ [ i, j ] = _f1_ [ i, j ] - _f2_ [ i, j ] ).
3021  * Number of tuples and components in the given fields must be the same.
3022  *  \param [in] f1 - a field to subtract from.
3023  *  \param [in] f2 - a field to subtract.
3024  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
3025  *          The caller is to delete this result field using decrRef() as it is no more
3026  *          needed.
3027  *  \throw If either \a f1 or \a f2 is NULL.
3028  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
3029  *         differ not only in values.
3030  */
3031 MEDCouplingFieldDouble *MEDCouplingFieldDouble::SubstractFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
3032 {
3033   if(!f1)
3034     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::SubstractFields : input field is NULL !");
3035   if(!f1->areStrictlyCompatible(f2))
3036     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply SubstractFields on them !");
3037   MEDCouplingTimeDiscretization *td=f1->_time_discr->substract(f2->_time_discr);
3038   td->copyTinyAttrFrom(*f1->_time_discr);
3039   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
3040   ret->setMesh(f1->getMesh());
3041   return ret.retn();
3042 }
3043
3044 /*!
3045  * Subtract values of another MEDCouplingFieldDouble from values of \a this one
3046  * ( _this_ [ i, j ] -= _other_ [ i, j ] ) using DataArrayDouble::substractEqual().
3047  * The two fields must have same number of tuples, components and same underlying mesh.
3048  *  \param [in] other - the field to subtract from \a this one.
3049  *  \return const MEDCouplingFieldDouble & - a reference to \a this field.
3050  *  \throw If \a other is NULL.
3051  *  \throw If the fields are not strictly compatible (areStrictlyCompatible()), i.e. they
3052  *         differ not only in values.
3053  */
3054 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator-=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
3055 {
3056   if(!areStrictlyCompatible(&other))
3057     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply -= on them !");
3058   _time_discr->substractEqual(other._time_discr);
3059   return *this;
3060 }
3061
3062 /*!
3063  * Returns a new MEDCouplingFieldDouble containing product values of
3064  * two given fields. There are 2 valid cases.
3065  * 1.  The fields have same number of tuples and components. Then each value of
3066  *   the result field (_f_) is a product of the corresponding values of _f1_ and
3067  *   _f2_, i.e. _f_ [ i, j ] = _f1_ [ i, j ] * _f2_ [ i, j ].
3068  * 2.  The fields have same number of tuples and one field, say _f2_, has one
3069  *   component. Then
3070  *   _f_ [ i, j ] = _f1_ [ i, j ] * _f2_ [ i, 0 ].
3071  *
3072  * The two fields must have same number of tuples and same underlying mesh.
3073  *  \param [in] f1 - a factor field.
3074  *  \param [in] f2 - another factor field.
3075  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
3076  *          The caller is to delete this result field using decrRef() as it is no more
3077  *          needed.
3078  *  \throw If either \a f1 or \a f2 is NULL.
3079  *  \throw If the fields are not compatible for production (areCompatibleForMul()),
3080  *         i.e. they differ not only in values and possibly number of components.
3081  */
3082 MEDCouplingFieldDouble *MEDCouplingFieldDouble::MultiplyFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
3083 {
3084   if(!f1)
3085     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::MultiplyFields : input field is NULL !");
3086   if(!f1->areCompatibleForMul(f2))
3087     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply MultiplyFields on them !");
3088   MEDCouplingTimeDiscretization *td=f1->_time_discr->multiply(f2->_time_discr);
3089   td->copyTinyAttrFrom(*f1->_time_discr);
3090   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
3091   ret->setMesh(f1->getMesh());
3092   return ret.retn();
3093 }
3094
3095 /*!
3096  * Multiply values of another MEDCouplingFieldDouble to values of \a this one
3097  * using DataArrayDouble::multiplyEqual().
3098  * The two fields must have same number of tuples and same underlying mesh.
3099  * There are 2 valid cases.
3100  * 1.  The fields have same number of components. Then each value of
3101  *   \a other is multiplied to the corresponding value of \a this field, i.e.
3102  *   _this_ [ i, j ] *= _other_ [ i, j ].
3103  * 2. The _other_ field has one component. Then
3104  *   _this_ [ i, j ] *= _other_ [ i, 0 ].
3105  *
3106  * The two fields must have same number of tuples and same underlying mesh.
3107  *  \param [in] other - an field to multiply to \a this one.
3108  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
3109  *          The caller is to delete this result field using decrRef() as it is no more
3110  *          needed.
3111  *  \throw If \a other is NULL.
3112  *  \throw If the fields are not strictly compatible for production
3113  *         (areCompatibleForMul()),
3114  *         i.e. they differ not only in values and possibly in number of components.
3115  */
3116 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator*=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
3117 {
3118   if(!areCompatibleForMul(&other))
3119     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply *= on them !");
3120   _time_discr->multiplyEqual(other._time_discr);
3121   return *this;
3122 }
3123
3124 /*!
3125  * Returns a new MEDCouplingFieldDouble containing division of two given fields.
3126  * There are 2 valid cases.
3127  * 1.  The fields have same number of tuples and components. Then each value of
3128  *   the result field (_f_) is a division of the corresponding values of \a f1 and
3129  *   \a f2, i.e. _f_ [ i, j ] = _f1_ [ i, j ] / _f2_ [ i, j ].
3130  * 2.  The fields have same number of tuples and _f2_ has one component. Then
3131  *   _f_ [ i, j ] = _f1_ [ i, j ] / _f2_ [ i, 0 ].
3132  *
3133  *  \param [in] f1 - a numerator field.
3134  *  \param [in] f2 - a denominator field.
3135  *  \return MEDCouplingFieldDouble * - the new instance of MEDCouplingFieldDouble.
3136  *          The caller is to delete this result field using decrRef() as it is no more
3137  *          needed.
3138  *  \throw If either \a f1 or \a f2 is NULL.
3139  *  \throw If the fields are not compatible for division (areCompatibleForDiv()),
3140  *         i.e. they differ not only in values and possibly in number of components.
3141  */
3142 MEDCouplingFieldDouble *MEDCouplingFieldDouble::DivideFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
3143 {
3144   if(!f1)
3145     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::DivideFields : input field is NULL !");
3146   if(!f1->areCompatibleForDiv(f2))
3147     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply DivideFields on them !");
3148   MEDCouplingTimeDiscretization *td=f1->_time_discr->divide(f2->_time_discr);
3149   td->copyTinyAttrFrom(*f1->_time_discr);
3150   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
3151   ret->setMesh(f1->getMesh());
3152   return ret.retn();
3153 }
3154
3155 /*!
3156  * Divide values of \a this field by values of another MEDCouplingFieldDouble
3157  * using DataArrayDouble::divideEqual().
3158  * The two fields must have same number of tuples and same underlying mesh.
3159  * There are 2 valid cases.
3160  * 1.  The fields have same number of components. Then each value of
3161  *    \a this field is divided by the corresponding value of \a other one, i.e.
3162  *   _this_ [ i, j ] /= _other_ [ i, j ].
3163  * 2.  The \a other field has one component. Then
3164  *   _this_ [ i, j ] /= _other_ [ i, 0 ].
3165  *
3166  *  \warning No check of division by zero is performed!
3167  *  \param [in] other - an field to divide \a this one by.
3168  *  \throw If \a other is NULL.
3169  *  \throw If the fields are not compatible for division (areCompatibleForDiv()),
3170  *         i.e. they differ not only in values and possibly in number of components.
3171  */
3172 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator/=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
3173 {
3174   if(!areCompatibleForDiv(&other))
3175     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply /= on them !");
3176   _time_discr->divideEqual(other._time_discr);
3177   return *this;
3178 }
3179
3180 /*!
3181  * Directly called by MEDCouplingFieldDouble::operator^.
3182  * 
3183  * \sa MEDCouplingFieldDouble::operator^
3184  */
3185 MEDCouplingFieldDouble *MEDCouplingFieldDouble::PowFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
3186 {
3187   if(!f1)
3188     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::PowFields : input field is NULL !");
3189   if(!f1->areCompatibleForMul(f2))
3190     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply PowFields on them !");
3191   MEDCouplingTimeDiscretization *td=f1->_time_discr->pow(f2->_time_discr);
3192   td->copyTinyAttrFrom(*f1->_time_discr);
3193   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->_type->clone());
3194   ret->setMesh(f1->getMesh());
3195   return ret.retn();
3196 }
3197
3198 /*!
3199  * Directly call MEDCouplingFieldDouble::PowFields static method.
3200  * 
3201  * \sa MEDCouplingFieldDouble::PowFields
3202  */
3203 MEDCouplingFieldDouble *MEDCouplingFieldDouble::operator^(const MEDCouplingFieldDouble& other) const throw(INTERP_KERNEL::Exception)
3204 {
3205   return PowFields(this,&other);
3206 }
3207
3208 const MEDCouplingFieldDouble &MEDCouplingFieldDouble::operator^=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception)
3209 {
3210   if(!areCompatibleForDiv(&other))
3211     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to apply /= on them !");
3212   _time_discr->powEqual(other._time_discr);
3213   return *this;
3214 }
3215
3216 /*!
3217  * Writes the field series \a fs and the mesh the fields lie on in the VTK file \a fileName.
3218  * If \a fs is empty no file is written.
3219  * The result file is valid provided that no exception is thrown.
3220  * \warning All the fields must be named and lie on the same non NULL mesh.
3221  *  \param [in] fileName - the name of a VTK file to write in.
3222  *  \param [in] fs - the fields to write.
3223  *  \param [in] isBinary - specifies the VTK format of the written file. By default true (Binary mode)
3224  *  \throw If \a fs[ 0 ] == NULL.
3225  *  \throw If the fields lie not on the same mesh.
3226  *  \throw If the mesh is not set.
3227  *  \throw If any of the fields has no name.
3228  *
3229  *  \if ENABLE_EXAMPLES
3230  *  \ref cpp_mcfielddouble_WriteVTK "Here is a C++ example".<br>
3231  *  \ref  py_mcfielddouble_WriteVTK "Here is a Python example".
3232  *  \endif
3233  */
3234 void MEDCouplingFieldDouble::WriteVTK(const std::string& fileName, const std::vector<const MEDCouplingFieldDouble *>& fs, bool isBinary)
3235 {
3236   if(fs.empty())
3237     return;
3238   std::size_t nfs=fs.size();
3239   if(!fs[0])
3240     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::WriteVTK : 1st instance of field is NULL !");
3241   const MEDCouplingMesh *m=fs[0]->getMesh();
3242   if(!m)
3243     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::WriteVTK : 1st instance of field lies on NULL mesh !");
3244   for(std::size_t i=1;i<nfs;i++)
3245     if(fs[i]->getMesh()!=m)
3246       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.");
3247   if(!m)
3248     throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::WriteVTK : Fields are lying on a same mesh but it is empty !");
3249   MEDCouplingAutoRefCountObjectPtr<DataArrayByte> byteArr;
3250   if(isBinary)
3251     { byteArr=DataArrayByte::New(); byteArr->alloc(0,1); }
3252   std::ostringstream coss,noss;
3253   for(std::size_t i=0;i<nfs;i++)
3254     {
3255       const MEDCouplingFieldDouble *cur=fs[i];
3256       std::string name(cur->getName());
3257       if(name.empty())
3258         {
3259           std::ostringstream oss; oss << "MEDCouplingFieldDouble::WriteVTK : Field in pos #" << i << " has no name !";
3260           throw INTERP_KERNEL::Exception(oss.str().c_str());
3261         }
3262       TypeOfField typ=cur->getTypeOfField();
3263       if(typ==ON_CELLS)
3264         cur->getArray()->writeVTK(coss,8,cur->getName(),byteArr);
3265       else if(typ==ON_NODES)
3266         cur->getArray()->writeVTK(noss,8,cur->getName(),byteArr);
3267       else
3268         throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::WriteVTK : only node and cell fields supported for the moment !");
3269     }
3270   m->writeVTKAdvanced(fileName,coss.str(),noss.str(),byteArr);
3271 }
3272
3273 void MEDCouplingFieldDouble::reprQuickOverview(std::ostream& stream) const
3274 {
3275   stream << "MEDCouplingFieldDouble C++ instance at " << this << ". Name : \"" << _name << "\"." << std::endl;
3276   const char *nat=0;
3277   try
3278     {
3279       nat=MEDCouplingNatureOfField::GetRepr(_nature);
3280       stream << "Nature of field : " << nat << ".\n";
3281     }
3282   catch(INTERP_KERNEL::Exception& /*e*/)
3283     {  }
3284   const MEDCouplingFieldDiscretization *fd(_type);
3285   if(!fd)
3286     stream << "No spatial discretization set !";
3287   else
3288     fd->reprQuickOverview(stream);
3289   stream << std::endl;
3290   if(!_mesh)
3291     stream << "\nNo mesh support defined !";
3292   else
3293     {
3294       std::ostringstream oss;
3295       _mesh->reprQuickOverview(oss);
3296       std::string tmp(oss.str());
3297       stream << "\nMesh info : " << tmp.substr(0,tmp.find('\n'));
3298     }
3299   if(_time_discr)
3300     {
3301       const DataArrayDouble *arr=_time_discr->getArray();
3302       if(arr)
3303         {
3304           stream << "\n\nArray info : ";
3305           arr->reprQuickOverview(stream);
3306         }
3307       else
3308         {
3309           stream << "\n\nNo data array set !";
3310         }
3311     }
3312 }