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