Salome HOME
cde465b131bfe34f50d0686e77b6bc26e14a3842
[modules/med.git] / src / MEDMEM / MEDMEM_Field.hxx
1 /*
2  File Field.hxx
3  $Header$
4 */
5
6 #ifndef FIELD_HXX
7 #define FIELD_HXX
8
9 #include <vector>
10 #include <algorithm>
11 #include <cmath>
12
13 #include "utilities.h"
14 #include "MEDMEM_Exception.hxx"
15 #include "MEDMEM_define.hxx"
16
17 #include "MEDMEM_Support.hxx"
18 #include "MEDMEM_Unit.hxx"
19 #include "MEDMEM_Array.hxx"
20 #include "MEDMEM_GenDriver.hxx"
21 #include "MEDMEM_DriverFactory.hxx"
22
23 /*!
24
25   This class contains all the informations related with a template class FIELD :
26   - Components descriptions
27   - Time step description
28   - Location of the values (a SUPPORT class)
29
30 */
31
32 namespace MEDMEM {
33
34   template<class T> class FIELD;
35
36 class FIELD_    // GENERIC POINTER TO a template <class T> class FIELD
37 {
38 protected:
39
40   bool            _isRead ;
41
42   /*!
43     \if developper
44     Field name.
45     \endif
46   */
47   string          _name ;
48   /*!
49     \if developper
50     Field description.
51     \endif
52   */
53   string          _description ;
54   /*!
55     \if developper
56     Pointer to the support the field deals with.
57     \endif
58   */
59   const SUPPORT * _support ;
60
61   /*!
62     \if developper
63     Number of field's components.
64     \endif
65   */
66   int             _numberOfComponents ;
67   /*!
68     \if developper
69     Number of field's values.
70     \endif
71   */
72   int             _numberOfValues ;
73
74   /*!
75     \if developper
76     Array of size _numberOfComponents. \n
77     (constant, scalar, vector, tensor)\n
78     We could use an array of integer to store
79     numbers of values: \n
80     - 1 for scalar,\n
81     - space dimension for vector,\n
82     - space dimension square for tensor.\n
83     So numbers of values per entities would be
84     sum of _componentsTypes array.
85
86     Not implemented yet! All type are scalar !
87     \endif
88   */
89   int *           _componentsTypes ;
90   /*!
91     \if developper
92     Array of size _numberOfComponents
93     storing components names if any.
94     \endif
95   */
96   string *        _componentsNames;     
97   /*!
98     \if developper
99     Array of size _numberOfComponents
100     storing components descriptions if any.
101     \endif
102  */
103   string *        _componentsDescriptions;
104   /*!
105     \if developper
106     Array of size _numberOfComponents
107     storing components units if any.
108     \endif
109  */
110   UNIT *          _componentsUnits;
111   /*!
112     \if developper
113     Array of size _numberOfComponents
114     storing components units if any.
115     \endif
116   */
117   string *        _MEDComponentsUnits;
118   int             _iterationNumber ;
119   double          _time;
120   int             _orderNumber ;
121
122   // _valueType should be a static const. Here is an initialization exemple 
123   // template < classType T > struct SET_VALUE_TYPE { static const med_type_champ _valueType = 0; }
124   // template < > struct SET_VALUE_TYPE<double> { static const med_type_champ _valueType = MED_EN::MED_REEL64; }
125   // template < > struct SET_VALUE_TYPE<int> { static const med_type_champ _valueType = MED_EN::MED_INT32; }
126   // static const med_type_champ _valueType = SET_VALUE_TYPE <T>::_valueType;
127   MED_EN::med_type_champ _valueType ;
128
129   vector<GENDRIVER *> _drivers; // Storage of the drivers currently in use
130   static void _checkFieldCompatibility(const FIELD_& m, const FIELD_& n ) throw (MEDEXCEPTION);
131   void _checkNormCompatibility(const FIELD<double>* p_field_volume=NULL) const  throw (MEDEXCEPTION);
132   FIELD<double>* _getFieldSize() const;
133
134 public:
135
136   friend class MED_MED_RDONLY_DRIVER;
137   friend class MED_MED_WRONLY_DRIVER;
138   friend class MED_MED_RDWR_DRIVER;
139
140   friend class VTK_MED_DRIVER;
141
142   /*!
143     Constructor.
144   */
145   FIELD_ ();
146   /*!
147     Constructor.
148   */
149   FIELD_(const SUPPORT * Support, const int NumberOfComponents);
150   /*!
151     Copy constructor.
152   */
153   FIELD_(const FIELD_ &m);
154
155   /*!
156     Destructor.
157   */
158   virtual ~FIELD_();
159
160 //    virtual  void     setIterationNumber (int IterationNumber);
161 //    virtual  void     setOrderNumber     (int OrderNumber);
162 //    virtual  void     setFieldName       (string& fieldName);
163
164   virtual  void     rmDriver(int index);
165   virtual   int     addDriver(driverTypes driverType,
166                               const string & fileName="Default File Name.med",
167                               const string & driverFieldName="Default Field Nam",
168                               MED_EN::med_mode_acces access=MED_EN::MED_REMP) ;
169
170   virtual  int      addDriver( GENDRIVER & driver);
171   virtual  void     read (const GENDRIVER &);
172   virtual  void     read(int index=0);
173   virtual void openAppend( void );
174   virtual  void     write(const GENDRIVER &);
175   virtual  void     write(int index=0, const string & driverName="");
176
177   virtual  void     writeAppend(const GENDRIVER &);
178   virtual  void     writeAppend(int index=0, const string & driverName="");
179
180   inline void     setName(const string Name);
181   inline string   getName() const;
182   inline void     setDescription(const string Description);
183   inline string   getDescription() const;
184   inline const SUPPORT * getSupport() const;
185   inline void     setSupport(const SUPPORT * support);
186   inline void     setNumberOfComponents(const int NumberOfComponents);
187   inline int      getNumberOfComponents() const;
188   inline void     setNumberOfValues(const int NumberOfValues);
189   inline int      getNumberOfValues() const;
190   //    inline void     setComponentType(int *ComponentType);
191   //    inline int *    getComponentType() const;
192   //    inline int      getComponentTypeI(int i) const;
193   inline void     setComponentsNames(const string * ComponentsNames);
194   inline void     setComponentName(int i, const string ComponentName);
195   inline const string * getComponentsNames() const;
196   inline string   getComponentName(int i) const;
197   inline void     setComponentsDescriptions(const string * ComponentsDescriptions);
198   inline void     setComponentDescription(int i, const string ComponentDescription);
199   inline const string * getComponentsDescriptions() const;
200   inline string   getComponentDescription(int i) const;
201
202   // provisoire : en attendant de regler le probleme des unites !
203   inline void     setComponentsUnits(const UNIT * ComponentsUnits);
204   inline const UNIT *   getComponentsUnits() const;
205   inline const UNIT *   getComponentUnit(int i) const;
206   inline void     setMEDComponentsUnits(const string * MEDComponentsUnits);
207   inline void     setMEDComponentUnit(int i, const string MEDComponentUnit);
208   inline const string * getMEDComponentsUnits() const;
209   inline string   getMEDComponentUnit(int i) const;
210
211   inline void     setIterationNumber(int IterationNumber);
212   inline int      getIterationNumber() const;
213   inline void     setTime(double Time);
214   inline double   getTime() const;
215   inline void     setOrderNumber(int OrderNumber);
216   inline int      getOrderNumber() const;
217
218   inline void     setValueType (const MED_EN::med_type_champ ValueType) ;
219   inline MED_EN::med_type_champ getValueType () const;
220
221 };
222
223 // ---------------------------------
224 // Implemented Methods : constructor
225 // ---------------------------------
226
227 // -----------------
228 // Methodes Inline
229 // -----------------
230 /*!
231   Set FIELD name.
232 */
233 inline void FIELD_::setName(const string Name)
234 {
235   _name=Name;
236 }
237 /*!
238   Get FIELD name.
239 */
240 inline string FIELD_::getName() const
241 {
242   return _name;
243 }
244 /*!
245   Set FIELD description.
246 */
247 inline void FIELD_::setDescription(const string Description)
248 {
249   _description=Description;
250 }
251 /*!
252   Get FIELD description.
253 */
254 inline string FIELD_::getDescription() const
255 {
256   return _description;
257 }
258 /*!
259   Set FIELD number of components.
260 */
261 inline void FIELD_::setNumberOfComponents(const int NumberOfComponents)
262 {
263   _numberOfComponents=NumberOfComponents;
264 }
265 /*!
266   Get FIELD number of components.
267 */
268 inline int FIELD_::getNumberOfComponents() const
269 {
270   return _numberOfComponents ;
271 }
272 /*!
273   Set FIELD number of values.
274
275   It must be the same than in the associated SUPPORT object.
276 */
277 inline void FIELD_::setNumberOfValues(const int NumberOfValues)
278 {
279   _numberOfValues=NumberOfValues;
280 }
281 /*!
282   Get FIELD number of value.
283 */
284 inline int FIELD_::getNumberOfValues() const
285 {
286   return _numberOfValues ;
287 }
288
289 //  inline void FIELD_::setComponentType(int *ComponentType)
290 //  {
291 //    _componentsTypes=ComponentType ;
292 //  }
293 //  inline int * FIELD_::getComponentType() const
294 //  {
295 //    return _componentsTypes ;
296 //  }
297 //  inline int FIELD_::getComponentTypeI(int i) const
298 //  {
299 //    return _componentsTypes[i-1] ;
300 //  }
301
302 /*!
303   Set FIELD components names.
304
305   Duplicate the ComponentsNames string array to put components names in
306   FIELD. ComponentsNames size must be equal to number of components.
307 */
308 inline void FIELD_::setComponentsNames(const string * ComponentsNames)
309 {
310   if (NULL == _componentsNames)
311     _componentsNames = new string[_numberOfComponents] ;
312   for (int i=0; i<_numberOfComponents; i++)
313     _componentsNames[i]=ComponentsNames[i] ;
314 }
315 /*!
316   Set FIELD i^th component name.
317
318   i must be >=1 and <= number of components.
319 */
320 inline void FIELD_::setComponentName(int i, const string ComponentName)
321 {
322   _componentsNames[i-1]=ComponentName ;
323 }
324 /*!
325   Get a reference to the string array which contain the components names.
326
327   This Array size is equal to number of components
328 */
329 inline const string * FIELD_::getComponentsNames() const
330 {
331   return _componentsNames ;
332 }
333 /*!
334   Get the name of the i^th component.
335 */
336 inline string FIELD_::getComponentName(int i) const
337 {
338   return _componentsNames[i-1] ;
339 }
340 /*!
341   Set FIELD components descriptions.
342
343   Duplicate the ComponentsDescriptions string array to put components 
344   descriptions in FIELD.
345   ComponentsDescriptions size must be equal to number of components.
346 */
347 inline void FIELD_::setComponentsDescriptions(const string * ComponentsDescriptions)
348 {
349   if (NULL == _componentsDescriptions)
350     _componentsDescriptions = new string[_numberOfComponents] ;
351   for (int i=0; i<_numberOfComponents; i++)
352     _componentsDescriptions[i]=ComponentsDescriptions[i] ;
353 }
354 /*!
355   Set FIELD i^th component description.
356
357   i must be >=1 and <= number of components.
358 */
359 inline void FIELD_::setComponentDescription(int i,const string ComponentDescription)
360 {
361   _componentsDescriptions[i-1]=ComponentDescription ;
362 }
363 /*!
364   Get a reference to the string array which contain the components descriptions.
365
366   This Array size is equal to number of components
367 */
368 inline const string * FIELD_::getComponentsDescriptions() const
369 {
370   return _componentsDescriptions ;
371 }
372 /*!
373   Get the description of the i^th component.
374 */
375 inline string FIELD_::getComponentDescription(int i) const
376 {
377   return _componentsDescriptions[i-1];
378 }
379
380 /*!
381   \todo
382   Set FIELD components UNIT.
383
384   Duplicate the ComponentsUnits UNIT array to put components 
385   units in FIELD.
386   ComponentsUnits size must be equal to number of components.
387 */
388 inline void FIELD_::setComponentsUnits(const UNIT * ComponentsUnits)
389 {
390   if (NULL == _componentsUnits)
391     _componentsUnits = new UNIT[_numberOfComponents] ;
392   for (int i=0; i<_numberOfComponents; i++)
393     _componentsUnits[i]=ComponentsUnits[i] ;
394 }
395 /*!
396   Get a reference to the UNIT array which contain the components units.
397
398   This Array size is equal to number of components
399 */
400 inline const UNIT * FIELD_::getComponentsUnits() const
401 {
402   return _componentsUnits ;
403 }
404 /*!
405   Get the UNIT of the i^th component.
406 */
407 inline const UNIT * FIELD_::getComponentUnit(int i) const
408 {
409   return &_componentsUnits[i-1] ;
410 }
411 /*!
412   Set FIELD components unit.
413
414   Duplicate the MEDComponentsUnits string array to put components 
415   units in FIELD.
416   MEDComponentsUnits size must be equal to number of components.
417   
418 */
419 inline void FIELD_::setMEDComponentsUnits(const string * MEDComponentsUnits)
420 {
421   if (NULL == _MEDComponentsUnits)
422     _MEDComponentsUnits = new string[_numberOfComponents] ;
423   for (int i=0; i<_numberOfComponents; i++)
424     _MEDComponentsUnits[i]=MEDComponentsUnits[i] ;
425 }
426 /*!
427   Set FIELD i^th component unit.
428
429   i must be >=1 and <= number of components.
430 */
431 inline void FIELD_::setMEDComponentUnit(int i, const string MEDComponentUnit)
432 {
433   _MEDComponentsUnits[i-1]=MEDComponentUnit ;
434 }
435 /*!
436   Get a reference to the string array which contain the components units.
437
438   This Array size is equal to number of components
439 */
440 inline const string * FIELD_::getMEDComponentsUnits() const
441 {
442   return _MEDComponentsUnits ;
443 }
444 /*!
445   Get the string for unit of the i^th component.
446 */
447 inline string FIELD_::getMEDComponentUnit(int i) const
448 {
449   return _MEDComponentsUnits[i-1] ;
450 }
451 /*!
452   Set the iteration number where FIELD has been calculated.
453 */
454 inline void FIELD_::setIterationNumber(int IterationNumber)
455 {
456   _iterationNumber=IterationNumber;
457 }
458 /*!
459   Get the iteration number where FIELD has been calculated.
460 */
461 inline int FIELD_::getIterationNumber() const
462 {
463   return _iterationNumber ;
464 }
465 /*!
466   Set the time (in second) where FIELD has been calculated.
467 */
468 inline void FIELD_::setTime(double Time)
469 {
470   _time=Time ;
471 }
472 /*!
473   Get the time (in second) where FIELD has been calculated.
474 */
475 inline double FIELD_::getTime() const
476 {
477   return _time ;
478 }
479 /*!
480   Set the order number where FIELD has been calculated.
481
482   It corresponds to internal iteration during one time step.
483 */
484 inline void FIELD_::setOrderNumber(int OrderNumber)
485 {
486   _orderNumber=OrderNumber ;
487 }
488 /*!
489   Get the order number where FIELD has been calculated.
490 */
491 inline int FIELD_::getOrderNumber() const
492 {
493   return _orderNumber ;
494 }
495 /*!
496   Get a reference to the SUPPORT object associated to FIELD.
497 */
498 inline  const SUPPORT * FIELD_::getSupport() const
499 {
500   return _support ;
501 }
502 /*!
503   Set the reference to the SUPPORT object associated to FIELD.
504
505   Reference is not duplicate, so it must not be deleted.
506 */
507 inline void FIELD_::setSupport(const SUPPORT * support)
508 {
509   _support = support ;
510 }
511 /*!
512   Get the FIELD med value type (MED_INT32 or MED_REEL64).
513 */
514 inline MED_EN::med_type_champ FIELD_::getValueType () const
515 {
516   return _valueType ;
517 }
518 /*!
519   Set the FIELD med value type (MED_INT32 or MED_REEL64).
520 */
521 inline void FIELD_::setValueType (const MED_EN::med_type_champ ValueType)
522 {
523   _valueType = ValueType ;
524 }
525
526 }//End namespace MEDMEM
527
528 /////////////////////////
529 // END OF CLASS FIELD_ //
530 /////////////////////////
531
532 /*!
533
534   This template class contains informations related with a FIELD :
535   - Values of the field
536
537 */
538
539 namespace MEDMEM {
540
541   template<class T2> class MED_FIELD_RDONLY_DRIVER;
542   template<class T2> class MED_FIELD_WRONLY_DRIVER;
543   template<class T2> class VTK_FIELD_DRIVER;
544   
545 template <class T> class FIELD : public FIELD_
546 {
547 protected:
548   // ------ End of Drivers Management Part
549
550   // array of value of type T
551   MEDARRAY<T> *_value ;
552
553 private:
554   void _operation(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode, char* Op);
555   void _operationInitialize(const FIELD& m,const FIELD& n, char* Op);
556   void _add_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode);
557   void _sub_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode);
558   void _mul_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode);
559   void _div_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode) throw (MEDEXCEPTION);
560   //setValueType() ;
561
562   FIELD & operator=(const FIELD &m);    // A FAIRE
563 public:
564   FIELD();
565   FIELD(const FIELD &m);
566   FIELD(const SUPPORT * Support, const int NumberOfComponents, const MED_EN::medModeSwitch Mode=MED_EN::MED_FULL_INTERLACE)  throw (MEDEXCEPTION) ; // Ajout NB Constructeur FIELD avec allocation de memoire de tous ses attribut
567   FIELD(const SUPPORT * Support, driverTypes driverType,
568         const string & fileName="", const string & fieldName="",
569         const int iterationNumber = -1, const int orderNumber = -1)
570     throw (MEDEXCEPTION);
571   ~FIELD();
572
573   const FIELD operator+(const FIELD& m) const;
574   const FIELD operator-(const FIELD& m) const;
575   const FIELD operator*(const FIELD& m) const;
576   const FIELD operator/(const FIELD& m) const;
577   const FIELD operator-() const;
578   FIELD& operator+=(const FIELD& m);
579   FIELD& operator-=(const FIELD& m);
580   FIELD& operator*=(const FIELD& m);
581   FIELD& operator/=(const FIELD& m);
582   static FIELD* add(const FIELD& m, const FIELD& n);
583   static FIELD* sub(const FIELD& m, const FIELD& n);
584   static FIELD* mul(const FIELD& m, const FIELD& n);
585   static FIELD* div(const FIELD& m, const FIELD& n);
586   double normMax() const throw (MEDEXCEPTION);
587   double norm2() const throw (MEDEXCEPTION);
588   void   applyLin(T a, T b);
589   template <T T_function(T)> void applyFunc();
590   static FIELD* scalarProduct(const FIELD& m, const FIELD& n);
591   double normL2(int component, const FIELD<double> * p_field_volume=NULL) const;
592   double normL2(const FIELD<double> * p_field_volume=NULL) const;
593   double normL1(int component, const FIELD<double> * p_field_volume=NULL) const;
594   double normL1(const FIELD<double> * p_field_volume=NULL) const;
595
596   friend class MED_FIELD_RDONLY_DRIVER<T>;
597   friend class MED_FIELD_WRONLY_DRIVER<T>;
598   friend class VTK_FIELD_DRIVER<T>;
599   //friend class MED_FIELD_RDWR_DRIVER  <T>;
600
601   void init ();
602   void rmDriver(int index=0);
603   int  addDriver(driverTypes driverType,
604                  const string & fileName="Default File Name.med",
605                  const string & driverFieldName="Default Field Name",
606                  MED_EN::med_mode_acces access=MED_EN::MED_REMP) ;
607
608   int  addDriver(GENDRIVER & driver);
609
610   void allocValue(const int NumberOfComponents);
611   void allocValue(const int NumberOfComponents, const int LengthValue);
612
613   void deallocValue();
614
615   inline void read(int index=0);
616   inline void read(const GENDRIVER & genDriver);
617   inline void write(int index=0, const string & driverName = "");
618   inline void write(const GENDRIVER &);
619
620   inline void writeAppend(int index=0, const string & driverName = "");
621   inline void writeAppend(const GENDRIVER &);
622
623   inline void     setValue(MEDARRAY<T> *Value);
624
625   inline MEDARRAY<T>* getvalue() const;
626   inline int getValueLength(MED_EN::medModeSwitch Mode) const;
627   inline const T*       getValue(MED_EN::medModeSwitch Mode) const;
628   inline const T*       getValueI(MED_EN::medModeSwitch Mode,int i) const;
629   inline T        getValueIJ(int i,int j) const;
630
631   inline void setValue(MED_EN::medModeSwitch mode, T* value);
632   inline void setValueI(MED_EN::medModeSwitch mode, int i, T* value);
633   inline void setValueIJ(int i, int j, T value);
634
635   /*!
636     This fonction feeds the FIELD<double> private attributs _value with the
637     volume of each cells belonging to the argument Support. The field has to be
638     initialised via the constructor FIELD<double>(const SUPPORT * , const int )
639     with Support as SUPPORT argument, 1 has the number of components, and Support
640     has be a SUPPORT on 3D cells. This initialisation could be done by the empty
641     constructor followed by a setSupport and setNumberOfComponents call but it has
642     be followed by a setValueType(MED_REEL64) call.
643    */
644   void getVolume() const throw (MEDEXCEPTION) ;
645   /*!
646     This fonction feeds the FIELD<double> private attributs _value with the
647     area of each cells (or faces) belonging to the attribut _support. The field
648     has to be initialised via the constructor FIELD<double>(const SUPPORT * ,
649     const int ) with 1 has the number of components, and _support has be a
650     SUPPORT on 2D cells or 3D faces. This initialisation could be done by the
651     empty constructor followed by a setSupport and setNumberOfComponents call but
652     it has be followed by a setValueType(MED_REEL64) call.
653    */
654   void getArea() const throw (MEDEXCEPTION) ;
655   /*!
656     This fonction feeds the FIELD<double> private attributs _value with the
657     length of each segments belonging to the attribut _support. The field has
658     to be initialised via the constructor FIELD<double>(const SUPPORT * ,
659     const int ) with 1 has the number of components, and _support has be a
660     SUPPORT on 3D edges or 2D faces. This initialisation could be done by the
661     empty constructor followed by a setSupport and setNumberOfComponents call but
662     it has be followed by a setValueType(MED_REEL64) call.
663    */
664   void getLength() const throw (MEDEXCEPTION) ;
665   /*!
666     This fonction feeds the FIELD<double> private attributs _value with the
667     normal vector of each faces belonging to the attribut _support. The field
668     has to be initialised via the constructor FIELD<double>(const SUPPORT * ,
669     const int ) with the space dimension has the number of components, and
670     _support has be a SUPPORT on 3D or 2D faces. This initialisation could be done
671     by the empty constructor followed by a setSupport and setNumberOfComponents
672     call but it has be followed by a setValueType(MED_REEL64) call.
673    */
674   void getNormal() const throw (MEDEXCEPTION) ;
675   /*!
676     This fonction feeds the FIELD<double> private attributs _value with the
677     barycenter of each faces or cells or edges belonging to the attribut _support.
678     The field has to be initialised via the constructor
679     FIELD<double>(const SUPPORT * ,const int ) with the space dimension has the
680     number of components, and _support has be a SUPPORT on 3D cells or 2D faces.
681     This initialisation could be done by the empty constructor followed by a
682     setSupport and setNumberOfComponents call but it has be followed by a
683     setValueType(MED_REEL64) call.
684    */
685   void getBarycenter() const throw (MEDEXCEPTION) ;
686 };
687
688 // --------------------
689 // Implemented Methods
690 // --------------------
691
692 /*!
693   Constructor with no parameter, most of the attribut members are set to NULL.
694 */
695 template <class T>  FIELD<T>::FIELD():
696  _value((MEDARRAY<T>*)NULL)
697 {
698   MESSAGE("Constructeur FIELD sans parametre");
699 }
700
701 /*!
702   Constructor with parameters such that all attrribut are set but the _value
703   attribut is allocated but not set.
704 */
705 template <class T>  FIELD<T>::FIELD(const SUPPORT * Support,
706                                     const int NumberOfComponents, const MED_EN::medModeSwitch Mode) throw (MEDEXCEPTION) :
707   FIELD_(Support, NumberOfComponents)
708 {
709   BEGIN_OF("FIELD<T>::FIELD(const SUPPORT * Support, const int NumberOfComponents, const medModeSwitch Mode)");
710   SCRUTE(this);
711
712   try {
713     _numberOfValues = Support->getNumberOfElements(MED_ALL_ELEMENTS);
714   }
715   catch (MEDEXCEPTION &ex) {
716     MESSAGE("No value defined ! ("<<ex.what()<<")");
717     _value = (MEDARRAY<T>*)NULL ;
718   }
719   MESSAGE("FIELD : constructeur : "<< _numberOfValues <<" et "<< NumberOfComponents);
720   if (0<_numberOfValues) {
721     _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues,Mode);
722     _isRead = true ;
723   } else
724     _value = (MEDARRAY<T>*)NULL ;
725
726   END_OF("FIELD<T>::FIELD(const SUPPORT * Support, const int NumberOfComponents, const medModeSwitch Mode)");
727 }
728
729 /*!
730   \if developper
731   \endif
732 */
733 template <class T> void FIELD<T>::init ()
734 {
735 }
736
737 /*!
738   Copy constructor.
739 */
740 template <class T> FIELD<T>::FIELD(const FIELD & m):
741   FIELD_((FIELD_) m)
742 {
743   MESSAGE("Constructeur FIELD de recopie");
744   if (m._value != NULL)
745     {
746       // copie only default !
747       _value = new MEDARRAY<T>(* m._value,false);
748     }
749   else
750     _value = (MEDARRAY<T> *) NULL;
751   //_drivers = m._drivers;
752 }
753
754 /*!
755   \if developper
756   Not implemented.
757   \endif
758 */
759 template <class T> FIELD<T> & FIELD<T>::operator=(const FIELD &m)
760 {
761   MESSAGE("Appel de FIELD<T>::operator=") ;
762   // operator= on FIELD_
763   // ignore driver
764   // copy values array
765 }
766
767 /*!
768      Overload addition operator.
769      This operation is authorized only for compatible fields that have the same support.
770      The compatibility checking includes equality tests of the folowing data members:\n
771      - _support
772      - _numberOfComponents
773      - _numberOfValues
774      - _componentsTypes
775      - _MEDComponentsUnits.
776
777      The data members of the returned field are initialized, based on the first field, except for the name, 
778      which is the combination of the two field's names and the operator.
779      Advised Utilisation in C++ : <tt> FIELD<T> c = a + b; </tt> \n
780      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
781      When using python, this operator calls the copy constructor in any case.
782      The user has to be aware that when using operator + in associatives expressions like
783      <tt> a = b + c + d +e; </tt> \n
784      no optimisation is performed : the evaluation of last expression requires the construction of
785      3 temporary fields.
786 */
787 template <class T>
788 const FIELD<T> FIELD<T>::operator+(const FIELD & m) const
789 {
790     BEGIN_OF("FIELD<T>::operator+(const FIELD & m)");
791     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
792
793     // Select mode : avoid if possible any calculation of other mode for fields m or *this
794     MED_EN::medModeSwitch mode;
795     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
796         mode=m.getvalue()->getMode();
797     else
798         mode=this->getvalue()->getMode();
799     
800     // Creation of the result - memory is allocated by FIELD constructor
801     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
802     //result._operation(*this,m,mode,"+"); // perform Atribute's initialization & addition
803     result._operationInitialize(*this,m,"+"); // perform Atribute's initialization
804     result._add_in_place(*this,m,mode); // perform addition
805
806     END_OF("FIELD<T>::operator+(const FIELD & m)");
807     return result;
808 }
809
810 /*!  Overloaded Operator +=
811  *   Operations are directly performed in the first field's array.
812  *   This operation is authorized only for compatible fields that have the same support.
813  */
814 template <class T>
815 FIELD<T>& FIELD<T>::operator+=(const FIELD & m)
816 {
817     BEGIN_OF("FIELD<T>::operator+=(const FIELD & m)");
818     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
819
820     // We choose to keep *this mode, even if it may cost a re-calculation for m
821     MED_EN::medModeSwitch mode=this->getvalue()->getMode();
822     const T* value1=m.getValue(mode); // get pointers to the values we are adding
823
824     // get a non const pointer to the inside array of values and perform operation
825     T * value=const_cast<T *> (getValue(mode));
826     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
827     const T* endV=value+size; // pointer to the end of value
828     for(;value!=endV; value1++,value++)
829         *value += *value1;
830     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
831     this->getvalue()->clearOtherMode();
832     END_OF("FIELD<T>::operator+=(const FIELD & m)");
833     return *this;
834 }
835
836
837 /*! Addition of fields. Static member function.
838  *  The function return a pointer to a new created field that holds the addition.
839  *  Data members are checked for compatibility and initialized.
840  *  The user is in charge of memory deallocation.
841  */
842 template <class T>
843 FIELD<T>* FIELD<T>::add(const FIELD& m, const FIELD& n)
844 {
845     BEGIN_OF("FIELD<T>::add(const FIELD & m, const FIELD& n)");
846     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
847
848     // Select mode : avoid if possible any calculation of other mode for fields m or *this
849     MED_EN::medModeSwitch mode;
850     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
851         mode=m.getvalue()->getMode();
852     else
853         mode=n.getvalue()->getMode();
854     
855     // Creation of a new field
856     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
857     result->_operationInitialize(m,n,"+"); // perform Atribute's initialization
858     result->_add_in_place(m,n,mode); // perform addition
859
860     END_OF("FIELD<T>::add(const FIELD & m, const FIELD& n)");
861     return result;
862 }
863
864 /*!
865      Overload substraction operator.
866      This operation is authorized only for compatible fields that have the same support.
867      The compatibility checking includes equality tests of the folowing data members:\n
868      - _support
869      - _numberOfComponents
870      - _numberOfValues
871      - _componentsTypes
872      - _MEDComponentsUnits.
873
874      The data members of the returned field are initialized, based on the first field, except for the name, 
875      which is the combination of the two field's names and the operator.
876      Advised Utilisation in C++ : <tt> FIELD<T> c = a - b; </tt> \n
877      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
878      When using python, this operator calls the copy constructor in any case.
879      The user has to be aware that when using operator - in associatives expressions like
880      <tt> a = b - c - d -e; </tt> \n
881      no optimisation is performed : the evaluation of last expression requires the construction of
882      3 temporary fields.
883 */
884 template <class T>
885 const FIELD<T> FIELD<T>::operator-(const FIELD & m) const
886 {
887     BEGIN_OF("FIELD<T>::operator-(const FIELD & m)");
888     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
889
890     MED_EN::medModeSwitch mode;
891     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
892         mode=m.getvalue()->getMode();
893     else
894         mode=this->getvalue()->getMode();
895     
896     // Creation of the result - memory is allocated by FIELD constructor
897     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
898     //result._operation(*this,m,mode,"-"); // perform Atribute's initialization & substraction
899     result._operationInitialize(*this,m,"-"); // perform Atribute's initialization
900     result._sub_in_place(*this,m,mode); // perform substracion
901
902     END_OF("FIELD<T>::operator-(const FIELD & m)");
903     return result;
904 }
905
906 template <class T>
907 const FIELD<T> FIELD<T>::operator-() const
908 {
909     BEGIN_OF("FIELD<T>::operator-()");
910
911     MED_EN::medModeSwitch mode=this->getvalue()->getMode();
912     // Creation of the result - memory is allocated by FIELD constructor
913     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
914     // Atribute's initialization 
915     result.setName("- "+getName());
916     result.setComponentsNames(getComponentsNames());
917     // not yet implemented    setComponentType(getComponentType());
918     result.setComponentsDescriptions(getComponentsDescriptions());
919     result.setMEDComponentsUnits(getMEDComponentsUnits());
920     result.setComponentsUnits(getComponentsUnits());
921     result.setIterationNumber(getIterationNumber());
922     result.setTime(getTime());
923     result.setOrderNumber(getOrderNumber());
924     result.setValueType(getValueType());
925
926     const T* value1=getValue(mode); 
927     // get a non const pointer to the inside array of values and perform operation
928     T * value=const_cast<T *> (result.getValue(mode));
929     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
930     const T* endV=value+size; // pointer to the end of value
931
932     for(;value!=endV; value1++,value++)
933         *value = -(*value1);
934     END_OF("FIELD<T>::operator-=(const FIELD & m)");
935     return result;
936 }
937
938 /*!  Overloaded Operator -=
939  *   Operations are directly performed in the first field's array.
940  *   This operation is authorized only for compatible fields that have the same support.
941  */
942 template <class T>
943 FIELD<T>& FIELD<T>::operator-=(const FIELD & m)
944 {
945     BEGIN_OF("FIELD<T>::operator-=(const FIELD & m)");
946     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
947
948     // We choose to keep *this mode, even if it may cost a re-calculation for m
949     MED_EN::medModeSwitch mode=this->getvalue()->getMode();
950     const T* value1=m.getValue(mode); // get pointers to the values we are adding
951
952     // get a non const pointer to the inside array of values and perform operation
953     T * value=const_cast<T *> (getValue(mode));
954     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
955     const T* endV=value+size; // pointer to the end of value
956
957     for(;value!=endV; value1++,value++)
958         *value -= *value1;
959     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
960     this->getvalue()->clearOtherMode();
961     END_OF("FIELD<T>::operator-=(const FIELD & m)");
962     return *this;
963 }
964
965
966 /*! Substraction of fields. Static member function.
967  *  The function return a pointer to a new created field that holds the substraction.
968  *  Data members are checked for compatibility and initialized.
969  *  The user is in charge of memory deallocation.
970  */
971 template <class T>
972 FIELD<T>* FIELD<T>::sub(const FIELD& m, const FIELD& n)
973 {
974     BEGIN_OF("FIELD<T>::sub(const FIELD & m, const FIELD& n)");
975     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
976
977     // Select mode : avoid if possible any calculation of other mode for fields m or *this
978     MED_EN::medModeSwitch mode;
979     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
980         mode=m.getvalue()->getMode();
981     else
982         mode=n.getvalue()->getMode();
983     
984     // Creation of a new field
985     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
986     result->_operationInitialize(m,n,"-"); // perform Atribute's initialization
987     result->_sub_in_place(m,n,mode); // perform substraction
988
989     END_OF("FIELD<T>::sub(const FIELD & m, const FIELD& n)");
990     return result;
991 }
992
993 /*!
994      Overload multiplication operator.
995      This operation is authorized only for compatible fields that have the same support.
996      The compatibility checking includes equality tests of the folowing data members:\n
997      - _support
998      - _numberOfComponents
999      - _numberOfValues
1000      - _componentsTypes
1001      - _MEDComponentsUnits.
1002
1003      The data members of the returned field are initialized, based on the first field, except for the name, 
1004      which is the combination of the two field's names and the operator.
1005      Advised Utilisation in C++ : <tt> FIELD<T> c = a * b; </tt> \n
1006      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
1007      When using python, this operator calls the copy constructor in any case.
1008      The user has to be aware that when using operator * in associatives expressions like
1009      <tt> a = b * c * d *e; </tt> \n
1010      no optimisation is performed : the evaluation of last expression requires the construction of
1011      3 temporary fields.
1012 */
1013 template <class T>
1014 const FIELD<T> FIELD<T>::operator*(const FIELD & m) const
1015 {
1016     BEGIN_OF("FIELD<T>::operator*(const FIELD & m)");
1017     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1018
1019     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1020     MED_EN::medModeSwitch mode;
1021     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
1022         mode=m.getvalue()->getMode();
1023     else
1024         mode=this->getvalue()->getMode();
1025     
1026     // Creation of the result - memory is allocated by FIELD constructor
1027     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
1028     //result._operation(*this,m,mode,"*"); // perform Atribute's initialization & multiplication
1029     result._operationInitialize(*this,m,"*"); // perform Atribute's initialization
1030     result._mul_in_place(*this,m,mode); // perform multiplication
1031
1032     END_OF("FIELD<T>::operator*(const FIELD & m)");
1033     return result;
1034 }
1035
1036 /*!  Overloaded Operator *=
1037  *   Operations are directly performed in the first field's array.
1038  *   This operation is authorized only for compatible fields that have the same support.
1039  */
1040 template <class T>
1041 FIELD<T>& FIELD<T>::operator*=(const FIELD & m)
1042 {
1043     BEGIN_OF("FIELD<T>::operator*=(const FIELD & m)");
1044     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1045
1046     // We choose to keep *this mode, even if it may cost a re-calculation for m
1047     MED_EN::medModeSwitch mode=this->getvalue()->getMode();
1048     const T* value1=m.getValue(mode); // get pointers to the values we are adding
1049
1050     // get a non const pointer to the inside array of values and perform operation
1051     T * value=const_cast<T *> (getValue(mode));
1052     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1053     const T* endV=value+size; // pointer to the end of value
1054
1055     for(;value!=endV; value1++,value++)
1056         *value *= *value1;
1057     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1058     this->getvalue()->clearOtherMode();
1059     END_OF("FIELD<T>::operator*=(const FIELD & m)");
1060     return *this;
1061 }
1062
1063
1064 /*! Multiplication of fields. Static member function.
1065  *  The function return a pointer to a new created field that holds the multiplication.
1066  *  Data members are checked for compatibility and initialized.
1067  *  The user is in charge of memory deallocation.
1068  */
1069 template <class T>
1070 FIELD<T>* FIELD<T>::mul(const FIELD& m, const FIELD& n)
1071 {
1072     BEGIN_OF("FIELD<T>::mul(const FIELD & m, const FIELD& n)");
1073     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
1074
1075     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1076     MED_EN::medModeSwitch mode;
1077     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
1078         mode=m.getvalue()->getMode();
1079     else
1080         mode=n.getvalue()->getMode();
1081     
1082     // Creation of a new field
1083     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
1084     result->_operationInitialize(m,n,"*"); // perform Atribute's initialization
1085     result->_mul_in_place(m,n,mode); // perform multiplication
1086
1087     END_OF("FIELD<T>::mul(const FIELD & m, const FIELD& n)");
1088     return result;
1089 }
1090
1091
1092 /*!
1093      Overload division operator.
1094      This operation is authorized only for compatible fields that have the same support.
1095      The compatibility checking includes equality tests of the folowing data members:\n
1096      - _support
1097      - _numberOfComponents
1098      - _numberOfValues
1099      - _componentsTypes
1100      - _MEDComponentsUnits.
1101
1102      The data members of the returned field are initialized, based on the first field, except for the name, 
1103      which is the combination of the two field's names and the operator.
1104      Advised Utilisation in C++ : <tt> FIELD<T> c = a / b; </tt> \n
1105      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
1106      When using python, this operator calls the copy constructor in any case.
1107      The user has to be aware that when using operator / in associatives expressions like
1108      <tt> a = b / c / d /e; </tt> \n
1109      no optimisation is performed : the evaluation of last expression requires the construction of
1110      3 temporary fields.
1111 */
1112 template <class T>
1113 const FIELD<T> FIELD<T>::operator/(const FIELD & m) const
1114 {
1115     BEGIN_OF("FIELD<T>::operator/(const FIELD & m)");
1116     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1117
1118     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1119     MED_EN::medModeSwitch mode;
1120     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
1121         mode=m.getvalue()->getMode();
1122     else
1123         mode=this->getvalue()->getMode();
1124     
1125     // Creation of the result - memory is allocated by FIELD constructor
1126     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
1127     //result._operation(*this,m,mode,"/"); // perform Atribute's initialization & division
1128     result._operationInitialize(*this,m,"/"); // perform Atribute's initialization
1129     result._div_in_place(*this,m,mode); // perform division
1130
1131     END_OF("FIELD<T>::operator/(const FIELD & m)");
1132     return result;
1133 }
1134
1135
1136 /*!  Overloaded Operator /=
1137  *   Operations are directly performed in the first field's array.
1138  *   This operation is authorized only for compatible fields that have the same support.
1139  */
1140 template <class T>
1141 FIELD<T>& FIELD<T>::operator/=(const FIELD & m)
1142 {
1143     BEGIN_OF("FIELD<T>::operator/=(const FIELD & m)");
1144     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1145
1146     // We choose to keep *this mode, even if it may cost a re-calculation for m
1147     MED_EN::medModeSwitch mode=this->getvalue()->getMode();
1148     const T* value1=m.getValue(mode); // get pointers to the values we are adding
1149
1150     // get a non const pointer to the inside array of values and perform operation
1151     T * value=const_cast<T *> (getValue(mode));
1152     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1153     const T* endV=value+size; // pointer to the end of value
1154
1155     for(;value!=endV; value1++,value++)
1156         *value /= *value1;
1157     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1158     this->getvalue()->clearOtherMode();
1159     END_OF("FIELD<T>::operator/=(const FIELD & m)");
1160     return *this;
1161 }
1162
1163
1164 /*! Division of fields. Static member function.
1165  *  The function return a pointer to a new created field that holds the division.
1166  *  Data members are checked for compatibility and initialized.
1167  *  The user is in charge of memory deallocation.
1168  */
1169 template <class T>
1170 FIELD<T>* FIELD<T>::div(const FIELD& m, const FIELD& n)
1171 {
1172     BEGIN_OF("FIELD<T>::div(const FIELD & m, const FIELD& n)");
1173     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
1174
1175     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1176     MED_EN::medModeSwitch mode;
1177     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
1178         mode=m.getvalue()->getMode();
1179     else
1180         mode=n.getvalue()->getMode();
1181     
1182     // Creation of a new field
1183     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
1184     result->_operationInitialize(m,n,"/"); // perform Atribute's initialization
1185     result->_div_in_place(m,n,mode); // perform division
1186
1187     END_OF("FIELD<T>::div(const FIELD & m, const FIELD& n)");
1188     return result;
1189 }
1190
1191
1192 /*!
1193   \if developper
1194   This internal method initialize the members of a new field created to hold the result of the operation Op .
1195   Initialization is based on the first field, except for the name, which is the combination of the two field's names
1196   and the operator.
1197   \endif
1198 */
1199 template <class T>
1200 void FIELD<T>::_operationInitialize(const FIELD& m,const FIELD& n, char* Op)
1201 {
1202     MESSAGE("Appel methode interne _add" << Op);
1203
1204     // Atribute's initialization (copy of the first field's attributes)
1205     // Other data members (_support, _numberOfValues) are initialized in the field's constr.
1206     setName(m.getName()+" "+Op+" "+n.getName());
1207     setComponentsNames(m.getComponentsNames());
1208     // not yet implemented    setComponentType(m.getComponentType());
1209     setComponentsDescriptions(m.getComponentsDescriptions());
1210     setMEDComponentsUnits(m.getMEDComponentsUnits());
1211
1212     // The following data member may differ from field m to n.
1213     // The initialization is done based on the first field.
1214     setComponentsUnits(m.getComponentsUnits());
1215     setIterationNumber(m.getIterationNumber());
1216     setTime(m.getTime());
1217     setOrderNumber(m.getOrderNumber());
1218     setValueType(m.getValueType());
1219 }
1220
1221
1222 /*!
1223   \if developper
1224   Internal method called by FIELD<T>::operator+ and FIELD<T>::add to perform addition "in place".
1225   This method is applied to a just created field with medModeSwitch mode.
1226   For this reason, the alternate mode doesn't need to be set to 0 after performing operation : 
1227   it doesn't exist!
1228   \endif
1229 */
1230 template <class T>
1231 void FIELD<T>::_add_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode)
1232 {
1233     // get pointers to the values we are adding
1234     const T* value1=m.getValue(mode);
1235     const T* value2=n.getValue(mode);
1236     // get a non const pointer to the inside array of values and perform operation
1237     T * value=const_cast<T *> (getValue(mode));
1238
1239     const int size=getNumberOfValues()*getNumberOfComponents();
1240     SCRUTE(size);
1241     const T* endV1=value1+size;
1242     for(;value1!=endV1; value1++,value2++,value++)
1243         *value=(*value1)+(*value2);
1244 }
1245
1246 /*!
1247   \if developper
1248   Internal method called by FIELD<T>::operator- and FIELD<T>::sub to perform substraction "in place".
1249   This method is applied to a just created field with medModeSwitch mode.
1250   For this reason, the alternate mode doesn't need to be set to 0 after performing operation : 
1251   it doesn't exist!
1252   \endif
1253 */
1254 template <class T>
1255 void FIELD<T>::_sub_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode)
1256 {
1257     // get pointers to the values we are adding
1258     const T* value1=m.getValue(mode);
1259     const T* value2=n.getValue(mode);
1260     // get a non const pointer to the inside array of values and perform operation
1261     T * value=const_cast<T *> (getValue(mode));
1262
1263     const int size=getNumberOfValues()*getNumberOfComponents();
1264     SCRUTE(size);
1265     const T* endV1=value1+size;
1266     for(;value1!=endV1; value1++,value2++,value++)
1267         *value=(*value1)-(*value2);
1268 }
1269
1270 /*!
1271   \if developper
1272   Internal method called by FIELD<T>::operator* and FIELD<T>::mul to perform multiplication "in place".
1273   This method is applied to a just created field with medModeSwitch mode.
1274   For this reason, the alternate mode doesn't need to be set to 0 after performing operation : 
1275   it doesn't exist!
1276   \endif
1277 */
1278 template <class T>
1279 void FIELD<T>::_mul_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode)
1280 {
1281     // get pointers to the values we are adding
1282     const T* value1=m.getValue(mode);
1283     const T* value2=n.getValue(mode);
1284     // get a non const pointer to the inside array of values and perform operation
1285     T * value=const_cast<T *> (getValue(mode));
1286
1287     const int size=getNumberOfValues()*getNumberOfComponents();
1288     SCRUTE(size);
1289     const T* endV1=value1+size;
1290     for(;value1!=endV1; value1++,value2++,value++)
1291         *value=(*value1)*(*value2);
1292 }
1293
1294 /*!
1295   \if developper
1296   Internal method called by FIELD<T>::operator/ and FIELD<T>::div to perform division "in place".
1297   This method is applied to a just created field with medModeSwitch mode.
1298   For this reason, the alternate mode doesn't need to be set to 0 after performing operation : 
1299   it doesn't exist!
1300   \endif
1301 */
1302 template <class T>
1303 void FIELD<T>::_div_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode) throw (MEDEXCEPTION)
1304 {
1305     // get pointers to the values we are adding
1306     const T* value1=m.getValue(mode);
1307     const T* value2=n.getValue(mode);
1308     // get a non const pointer to the inside array of values and perform operation
1309     T * value=const_cast<T *> (getValue(mode));
1310
1311     const int size=getNumberOfValues()*getNumberOfComponents();
1312     SCRUTE(size);
1313     const T* endV1=value1+size;
1314     for(;value1!=endV1; value1++,value2++,value++){
1315       if ( *value2 == 0 ) {
1316           string diagnosis;
1317           diagnosis="FIELD<T>::_div_in_place(...) : Divide by zero !";
1318           throw MEDEXCEPTION(diagnosis.c_str());
1319         }
1320         *value=(*value1)/(*value2);
1321     }
1322 }
1323
1324 /*!  Return Max Norm 
1325  */
1326 template <class T> double FIELD<T>::normMax() const throw (MEDEXCEPTION)
1327 {
1328     const T* value=getValue(getvalue()->getMode()); // get pointer to the values
1329     const int size=getNumberOfValues()*getNumberOfComponents();
1330     if (size <= 0) // Size of array has to be strictly positive
1331     {
1332         string diagnosis;
1333         diagnosis="FIELD<T>::normMax() : cannot compute the norm of "+getName()+
1334             " : it size is non positive!";
1335         throw MEDEXCEPTION(diagnosis.c_str());
1336     }
1337     const T* lastvalue=value+size; // get pointer just after last value
1338     const T* pMax=value; // pointer to the max value
1339     const T* pMin=value; // pointer to the min value
1340
1341     // get pointers to the max & min value of array
1342     while ( ++value != lastvalue )
1343     {
1344         if ( *pMin > *value )
1345             pMin=value;
1346         if ( *pMax < *value )
1347             pMax=value;
1348     }
1349
1350     T Max= *pMax>(T) 0 ? *pMax : -*pMax; // Max=abs(*pMax)
1351     T Min= *pMin>(T) 0 ? *pMin : -*pMin; // Min=abs(*pMin)
1352
1353     return Max>Min ? static_cast<double>(Max) : static_cast<double>(Min);
1354 }
1355
1356 /*!  Return Euclidien norm 
1357  */
1358 template <class T> double FIELD<T>::norm2() const throw (MEDEXCEPTION)
1359 {
1360     const T* value=this->getValue(this->getvalue()->getMode()); // get const pointer to the values
1361     const int size=getNumberOfValues()*getNumberOfComponents(); // get size of array
1362     if (size <= 0) // Size of array has to be strictly positive
1363     {
1364         string diagnosis;
1365         diagnosis="FIELD<T>::norm2() : cannot compute the norm of "+getName()+
1366             " : it size is non positive!";
1367         throw MEDEXCEPTION(diagnosis.c_str());
1368     }
1369     const T* lastvalue=value+size; // point just after last value
1370     
1371     T result((T)0); // init
1372     for( ; value!=lastvalue ; ++value)
1373         result += (*value) * (*value);
1374     
1375     return std::sqrt(static_cast<double> (result));
1376 }
1377
1378
1379 /*!  Apply to each (scalar) field component the template parameter T_function, 
1380  *   which is a pointer to function.
1381  *   Since the pointer is known at compile time, the function is inlined into the inner loop!
1382  *   calculation is done "in place".
1383  *   Use examples : 
1384  *   
1385  *   \code  myField.applyFunc<std::sqrt>();  // apply sqare root function \endcode
1386  *     \code myField.applyFunc<myFunction>(); // apply your own created function \endcode
1387  */
1388 template <class T> template <T T_function(T)> 
1389 void FIELD<T>::applyFunc()
1390 {
1391     MED_EN::medModeSwitch mode=getvalue()->getMode();
1392
1393     // get a non const pointer to the inside array of values and perform operation
1394     T * value=const_cast<T *> (getValue(mode));
1395     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1396
1397     if (size>0) // for a negative size, there is nothing to do
1398     {
1399         const T* lastvalue=value+size; // pointer to the end of value
1400         for(;value!=lastvalue; ++value) // apply linear transformation
1401             *value = T_function(*value);
1402         // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1403         getvalue()->clearOtherMode();
1404     }
1405 }
1406     
1407   
1408 /*!  Apply to each (scalar) field component the linear function x -> ax+b.
1409  *   calculation is done "in place".
1410  */
1411 template <class T> void FIELD<T>::applyLin(T a, T b)
1412 {
1413     MED_EN::medModeSwitch mode=getvalue()->getMode();
1414
1415     // get a non const pointer to the inside array of values and perform operation in place
1416     T * value=const_cast<T *> (getValue(mode));
1417     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1418
1419     if (size>0) // for a negative size, there is nothing to do
1420     {
1421         const T* lastvalue=value+size; // pointer to the end of value
1422         for(;value!=lastvalue; ++value) // apply linear transformation
1423             *value = a*(*value)+b;
1424         // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1425         getvalue()->clearOtherMode();
1426     }
1427 }
1428
1429
1430 /*!
1431  *   Return a pointer to a new field that holds the scalar product. Static member function.
1432  *   This operation is authorized only for compatible fields that have the same support.
1433  *   The compatibility checking includes equality tests of the folowing data members:\n
1434  *   - _support
1435  *   - _numberOfComponents
1436  *   - _numberOfValues
1437  *   - _componentsTypes
1438  *   - _MEDComponentsUnits.
1439  *   Data members are initialized.
1440  *   The new field point to the same support and has one component.
1441  *   Each value of it is the scalar product of the two argument's fields.
1442  *   The user is in charge of memory deallocation.
1443  */
1444 template <class T> FIELD<T>* FIELD<T>::scalarProduct(const FIELD & m, const FIELD & n)
1445 {
1446     FIELD_::_checkFieldCompatibility( m, n); // may throw exception
1447     // we need a MED_FULL_INTERLACE representation of m & n to compute the scalar product
1448     const MED_EN::medModeSwitch mode=MED_EN::MED_FULL_INTERLACE; 
1449
1450     const int numberOfElements=m.getNumberOfValues(); // strictly positive
1451     const int NumberOfComponents=m.getNumberOfComponents(); // strictly positive
1452
1453     // Creation & init of a the result field on the same suppot, with one component
1454     FIELD<T>* result = new FIELD<T>(m.getSupport(),1,mode);
1455     result->setName( "scalarProduct ( " + m.getName() + " , " + n.getName() + " )" );
1456     result->setIterationNumber(m.getIterationNumber());
1457     result->setTime(m.getTime());
1458     result->setOrderNumber(m.getOrderNumber());
1459     result->setValueType(m.getValueType());
1460
1461     const T* value1=m.getValue(mode); // get const pointer to the values
1462     const T* value2=n.getValue(mode); // get const pointer to the values
1463     // get a non const pointer to the inside array of values and perform operation
1464     T * value=const_cast<T *> (result->getValue(mode));
1465     
1466     const T* lastvalue=value+numberOfElements; // pointing just after last value of result
1467     for ( ; value!=lastvalue ; ++value ) // loop on all elements
1468     {
1469         *value=(T)0; // initialize value
1470         const T* endofRow=value1+NumberOfComponents; // pointing just after end of row
1471         for ( ; value1 != endofRow; ++value1, ++value2) // computation of dot product
1472             *value += (*value1) * (*value2);
1473     }
1474     return result;
1475 }
1476
1477 /*!  Return L2 Norm  of the field's component.
1478  *   Cannot be applied to a field with a support on nodes.
1479  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1480  */
1481 template <class T> double FIELD<T>::normL2(int component, const FIELD<double> * p_field_volume) const
1482 {
1483     _checkNormCompatibility(p_field_volume); // may throw exception
1484     if ( component<1 || component>getNumberOfComponents() )
1485         throw MEDEXCEPTION(STRING("FIELD<T>::normL2() : The component argument should be between 1 and the number of components"));
1486
1487     const FIELD<double> * p_field_size=p_field_volume;
1488     if(!p_field_volume) // if the user don't supply the volume
1489         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1490
1491     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1492     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE); 
1493     const T* value=getValueI( MED_NO_INTERLACE, component); // get pointer to the component's values
1494     const T* lastvalue=value+getNumberOfValues(); // pointing just after the end of column
1495
1496     double integrale=0.0;
1497     double totVol=0.0;
1498     for (; value!=lastvalue ; ++value ,++vol) 
1499     {
1500         integrale += static_cast<double>((*value) * (*value)) * (*vol);
1501         totVol+=*vol;
1502     }
1503
1504     if(!p_field_volume) // if the user didn't supply the volume
1505         delete p_field_size; // delete temporary volume field
1506     if( totVol <= 0)
1507         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1508
1509     return integrale/totVol;
1510 }
1511
1512 /*!  Return L2 Norm  of the field.
1513  *   Cannot be applied to a field with a support on nodes.
1514  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1515  */
1516 template <class T> double FIELD<T>::normL2(const FIELD<double> * p_field_volume) const
1517 {
1518     _checkNormCompatibility(p_field_volume); // may throw exception
1519     const FIELD<double> * p_field_size=p_field_volume;
1520     if(!p_field_volume) // if the user don't supply the volume
1521         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1522
1523     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1524     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE); 
1525     const double* lastvol=vol+getNumberOfValues(); // pointing just after the end of vol
1526     const T* value=getValue( MED_NO_INTERLACE); // get pointer to the field's values
1527
1528     double totVol=0.0;
1529     const double* p_vol=vol;
1530     for (p_vol=vol; p_vol!=lastvol ; ++p_vol) // calculate total volume
1531         totVol+=*p_vol;
1532     
1533     double integrale=0.0;
1534     for (int i=1; i<=getNumberOfComponents(); ++i) // compute integral on all components
1535         for (p_vol=vol; p_vol!=lastvol ; ++value ,++p_vol) 
1536             integrale += static_cast<double>((*value) * (*value)) * (*p_vol);
1537
1538     if(!p_field_volume) // if the user didn't supply the volume
1539         delete p_field_size; // delete temporary volume field
1540     if( totVol <= 0)
1541         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1542
1543     return integrale/totVol;
1544 }
1545
1546 /*!  Return L1 Norm  of the field's component.
1547  *   Cannot be applied to a field with a support on nodes.
1548  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1549  */
1550 template <class T> double FIELD<T>::normL1(int component, const FIELD<double> * p_field_volume) const
1551 {
1552     _checkNormCompatibility(p_field_volume); // may throw exception
1553     if ( component<1 || component>getNumberOfComponents() )
1554         throw MEDEXCEPTION(STRING("FIELD<T>::normL2() : The component argument should be between 1 and the number of components"));
1555
1556     const FIELD<double> * p_field_size=p_field_volume;
1557     if(!p_field_volume) // if the user don't supply the volume
1558         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1559
1560     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1561     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE); 
1562     const T* value=getValueI( MED_NO_INTERLACE, component); // get pointer to the component's values
1563     const T* lastvalue=value+getNumberOfValues(); // pointing just after the end of column
1564
1565     double integrale=0.0;
1566     double totVol=0.0;
1567     for (; value!=lastvalue ; ++value ,++vol) 
1568     {
1569         integrale += std::abs( static_cast<double>(*value) ) * (*vol);
1570         totVol+=*vol;
1571     }
1572
1573     if(!p_field_volume) // if the user didn't supply the volume
1574         delete p_field_size; // delete temporary volume field
1575     if( totVol <= 0)
1576         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1577
1578     return integrale/totVol;
1579 }
1580
1581 /*!  Return L1 Norm  of the field.
1582  *   Cannot be applied to a field with a support on nodes.
1583  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1584  */
1585 template <class T> double FIELD<T>::normL1(const FIELD<double> * p_field_volume) const
1586 {
1587     _checkNormCompatibility(p_field_volume); // may throw exception
1588     const FIELD<double> * p_field_size=p_field_volume;
1589     if(!p_field_volume) // if the user don't supply the volume
1590         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1591
1592     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1593     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE); 
1594     const double* lastvol=vol+getNumberOfValues(); // pointing just after the end of vol
1595     const T* value=getValue( MED_NO_INTERLACE); // get pointer to the field's values
1596
1597     double totVol=0.0;
1598     const double* p_vol=vol;
1599     for (p_vol=vol; p_vol!=lastvol ; ++p_vol) // calculate total volume
1600         totVol+=*p_vol;
1601     
1602     double integrale=0.0;
1603     for (int i=1; i<=getNumberOfComponents(); ++i) // compute integral on all components
1604         for (p_vol=vol; p_vol!=lastvol ; ++value ,++p_vol) 
1605             integrale += std::abs( static_cast<double>(*value) ) * (*p_vol);
1606
1607     if(!p_field_volume) // if the user didn't supply the volume
1608         delete p_field_size; // delete temporary volume field
1609     if( totVol <= 0)
1610         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1611
1612     return integrale/totVol;
1613 }
1614
1615
1616
1617
1618 /*!
1619   Constructor with parameters; the object is set via a file and its associated
1620   driver. For the moment only the MED_DRIVER is considered and if the last two
1621   argument (iterationNumber and orderNumber) are not set; their default value
1622   is -1. If the field fieldDriverName with the iteration number
1623   iterationNumber and the order number orderNumber does not exist in the file
1624   fieldDriverName; the constructor raises an exception.
1625 */
1626 template <class T> FIELD<T>::FIELD(const SUPPORT * Support,
1627                                    driverTypes driverType,
1628                                    const string & fileName/*=""*/,
1629                                    const string & fieldDriverName/*=""*/,
1630                                    const int iterationNumber,
1631                                    const int orderNumber)
1632   throw (MEDEXCEPTION)
1633 {
1634   const char * LOC = "template <class T> FIELD<T>::FIELD(const SUPPORT * Support, driverTypes driverType, const string & fileName=\"\", const string & fieldName=\"\", const int iterationNumber=-1, const int orderNumber=-1) : ";
1635
1636   int current;
1637
1638   BEGIN_OF(LOC);
1639
1640   init();
1641
1642   _support = Support;
1643   _value = (MEDARRAY<T>*)NULL;
1644
1645   _iterationNumber = iterationNumber;
1646   _time = 0.0;
1647   _orderNumber = orderNumber;
1648
1649   current = addDriver(driverType,fileName,fieldDriverName,MED_LECT);
1650
1651 //   switch(driverType)
1652 //     {
1653 //     case MED_DRIVER :
1654 //       {
1655 //      MED_FIELD_RDONLY_DRIVER<T> myDriver(fileName,this);
1656 //      myDriver.setFieldName(fieldDriverName);
1657 //      current = addDriver(myDriver);
1658 //      break;
1659 //       }
1660 //   current = addDriver(driverType,fileName,fieldDriverName);
1661 //   switch(_drivers[current]->getAccessMode() ) {
1662 //   case MED_WRONLY : {
1663 //     MESSAGE("FIELD<T>::FIELD(driverTypes driverType, .....) : driverType must have a MED_RDONLY or MED_RDWR accessMode");
1664 //     rmDriver(current);
1665 //     break;}
1666 //   default : {
1667 //   }
1668 //   }
1669 //     default :
1670 //       {
1671 //      throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Driver type unknown : can't create driver!"));
1672 //      break;
1673 //       }
1674 //     }
1675
1676   _drivers[current]->open();
1677   _drivers[current]->read();
1678   _drivers[current]->close();
1679
1680   END_OF(LOC);
1681 }
1682
1683 /*!
1684   Destructor.
1685 */
1686 template <class T> FIELD<T>::~FIELD()
1687 {
1688   BEGIN_OF(" Destructeur FIELD<T>::~FIELD()");
1689   SCRUTE(this);
1690   if (_value) delete _value;
1691   END_OF(" Destructeur FIELD<T>::~FIELD()");
1692 }
1693
1694 /*!
1695   
1696 */
1697 template <class T> void FIELD<T>::allocValue(const int NumberOfComponents)
1698 {
1699   const char* LOC = "FIELD<T>::allocValue(const int NumberOfComponents)" ;
1700   BEGIN_OF(LOC);
1701
1702   _numberOfComponents = NumberOfComponents ;
1703   if (_componentsTypes == NULL)
1704     _componentsTypes = new int[NumberOfComponents] ;
1705   if (_componentsNames == NULL)
1706     _componentsNames = new string[NumberOfComponents];
1707   if (_componentsDescriptions == NULL)
1708     _componentsDescriptions = new string[NumberOfComponents];
1709   if (_componentsUnits == NULL)
1710     _componentsUnits = new UNIT[NumberOfComponents];
1711   if (_MEDComponentsUnits == NULL)
1712     _MEDComponentsUnits = new string[NumberOfComponents];
1713   for (int i=0;i<NumberOfComponents;i++) {
1714     _componentsTypes[i] = 0 ;
1715   }
1716
1717   try {
1718     _numberOfValues = _support->getNumberOfElements(MED_ALL_ELEMENTS);
1719     MESSAGE(LOC <<" : "<<_numberOfValues <<" et "<< NumberOfComponents);
1720
1721     _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues);
1722
1723     _isRead = true ;
1724   }
1725   catch (MEDEXCEPTION &ex) {
1726     MESSAGE("No value defined, problem with NumberOfComponents (and may be _support) size of MEDARRAY<T>::_value !");
1727     _value = (MEDARRAY<T>*)NULL ;
1728   }
1729
1730   SCRUTE(_value);
1731   END_OF("void FIELD<T>::allocValue(const int NumberOfComponents)");
1732 }
1733
1734 /*!
1735   
1736 */
1737 template <class T> void FIELD<T>::allocValue(const int NumberOfComponents, const int LengthValue)
1738 {
1739   BEGIN_OF("void FIELD<T>::allocValue(const int NumberOfComponents,const int LengthValue)");
1740
1741   _numberOfComponents = NumberOfComponents ;
1742   if (_componentsTypes == NULL)
1743     _componentsTypes = new int[NumberOfComponents] ;
1744   if (_componentsNames == NULL)
1745     _componentsNames = new string[NumberOfComponents];
1746   if (_componentsDescriptions == NULL)
1747     _componentsDescriptions = new string[NumberOfComponents];
1748   if (_componentsUnits == NULL)
1749     _componentsUnits = new UNIT[NumberOfComponents];
1750   if (_MEDComponentsUnits == NULL)
1751     _MEDComponentsUnits = new string[NumberOfComponents];
1752   for (int i=0;i<NumberOfComponents;i++) {
1753     _componentsTypes[i] = 0 ;
1754   }
1755
1756   MESSAGE("FIELD : constructeur : "<<LengthValue <<" et "<< NumberOfComponents);
1757   _numberOfValues = LengthValue ;
1758   _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues);
1759   _isRead = true ;
1760
1761   SCRUTE(_value);
1762   END_OF("void FIELD<T>::allocValue(const int NumberOfComponents,const int LengthValue)");
1763 }
1764
1765 /*!
1766   
1767 */
1768 template <class T> void FIELD<T>::deallocValue()
1769 {
1770   BEGIN_OF("void FIELD<T>::deallocValue()");
1771   _numberOfValues = 0 ;
1772   _numberOfComponents = 0 ;
1773   if (_value != NULL)
1774     delete _value;
1775
1776   END_OF("void FIELD<T>::deallocValue()");
1777 }
1778
1779 // -----------------
1780 // Methodes Inline
1781 // -----------------
1782
1783 /*!
1784   Create the specified driver and return its index reference to path to 
1785   read or write methods.
1786 */
1787
1788 template <class T> int FIELD<T>::addDriver(driverTypes driverType,
1789                                            const string & fileName/*="Default File Name.med"*/,
1790                                            const string & driverName/*="Default Field Name"*/,
1791                                            MED_EN::med_mode_acces access)
1792 {
1793   const char * LOC = "FIELD<T>::addDriver(driverTypes driverType, const string & fileName=\"Default File Name.med\",const string & driverName=\"Default Field Name\,MED_EN::med_mode_acces access) : ";
1794
1795   GENDRIVER * driver;
1796
1797   BEGIN_OF(LOC);
1798
1799   SCRUTE(driverType);
1800
1801   driver = DRIVERFACTORY::buildDriverForField(driverType,fileName,this,access);
1802
1803   _drivers.push_back(driver);
1804
1805   int current = _drivers.size()-1;
1806
1807   _drivers[current]->setFieldName(driverName);
1808
1809   END_OF(LOC);
1810
1811   return current;
1812 }
1813
1814
1815 /*!
1816   Duplicate the given driver and return its index reference to path to 
1817   read or write methods.
1818 */
1819 template <class T> inline int FIELD<T>::addDriver (GENDRIVER & driver )
1820 {
1821   const char * LOC = "FIELD<T>::addDriver(GENDRIVER &) : ";
1822   int current;
1823
1824   BEGIN_OF(LOC);
1825
1826   // duplicate driver to delete it with destructor !
1827   GENDRIVER * newDriver = driver.copy() ;
1828
1829   _drivers.push_back(newDriver);
1830
1831   current = _drivers.size()-1;
1832   SCRUTE(current);
1833   driver.setId(current); 
1834
1835   MESSAGE(LOC << " je suis la 1");
1836   END_OF(LOC);
1837   MESSAGE(LOC << " je suis la 2");
1838
1839   return current ;
1840 };
1841
1842 /*!
1843   Remove the driver referenced by its index.
1844 */
1845 template <class T> void FIELD<T>::rmDriver (int index/*=0*/)
1846 {
1847   const char * LOC = "FIELD<T>::rmDriver (int index=0): ";
1848   BEGIN_OF(LOC);
1849
1850   if ( _drivers[index] ) {
1851     //_drivers.erase(&_drivers[index]);
1852     // why not ????
1853     MESSAGE ("detruire");
1854   }
1855   else
1856     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1857                                      << "The <index given is invalid, index must be between  0 and  |"
1858                                      << _drivers.size()
1859                                      )
1860                           );
1861
1862   END_OF(LOC);
1863 }
1864
1865 /*!
1866   Read FIELD in the file specified in the driver given by its index.
1867 */
1868 template <class T> inline  void FIELD<T>::read(int index/*=0*/)
1869 {
1870   const char * LOC = "FIELD<T>::read(int index=0) : ";
1871   BEGIN_OF(LOC);
1872
1873   if ( _drivers[index] ) {
1874     _drivers[index]->open();
1875     _drivers[index]->read();
1876     _drivers[index]->close();
1877   }
1878   else
1879     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1880                                      << "The index given is invalid, index must be between  0 and |"
1881                                      << _drivers.size()
1882                                      )
1883                           );
1884   END_OF(LOC);
1885 }
1886
1887 /*!
1888   Write FIELD in the file specified in the driver given by its index.
1889 */
1890 template <class T> inline void FIELD<T>::write(int index/*=0*/, const string & driverName /*= ""*/)
1891 {
1892   const char * LOC = "FIELD<T>::write(int index=0, const string & driverName = \"\") : ";
1893   BEGIN_OF(LOC);
1894
1895   if( _drivers[index] ) {
1896     _drivers[index]->open();
1897     if (driverName != "") _drivers[index]->setFieldName(driverName);
1898     _drivers[index]->write();
1899     _drivers[index]->close();
1900   }
1901   else
1902     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1903                                      << "The index given is invalid, index must be between  0 and |"
1904                                      << _drivers.size()
1905                                      )
1906                           );
1907   END_OF(LOC);
1908 }
1909
1910 /*!
1911   Write FIELD in the file specified in the driver given by its index. Use this
1912   method for ASCII drivers (e.g. VTK_DRIVER)
1913 */
1914 template <class T> inline void FIELD<T>::writeAppend(int index/*=0*/, const string & driverName /*= ""*/)
1915 {
1916   const char * LOC = "FIELD<T>::write(int index=0, const string & driverName = \"\") : ";
1917   BEGIN_OF(LOC);
1918
1919   if( _drivers[index] ) {
1920     _drivers[index]->openAppend();
1921     if (driverName != "") _drivers[index]->setFieldName(driverName);
1922     _drivers[index]->writeAppend();
1923     _drivers[index]->close();
1924   }
1925   else
1926     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1927                                      << "The index given is invalid, index must be between  0 and |"
1928                                      << _drivers.size()
1929                                      )
1930                           );
1931   END_OF(LOC);
1932 }
1933
1934 /*!
1935   \internal
1936   Write FIELD with the driver which is equal to the given driver.
1937
1938   Use by MED object.
1939 */
1940 template <class T> inline void FIELD<T>::write(const GENDRIVER & genDriver)
1941 {
1942   const char * LOC = " FIELD<T>::write(const GENDRIVER &) : ";
1943   BEGIN_OF(LOC);
1944
1945   for (unsigned int index=0; index < _drivers.size(); index++ )
1946     if ( *_drivers[index] == genDriver ) {
1947       _drivers[index]->open();
1948       _drivers[index]->write();
1949       _drivers[index]->close();
1950     }
1951
1952   END_OF(LOC);
1953
1954 }
1955
1956 /*!
1957   \internal
1958   Write FIELD with the driver which is equal to the given driver.
1959
1960   Use by MED object. Use this method for ASCII drivers (e.g. VTK_DRIVER).
1961 */
1962 template <class T> inline void FIELD<T>::writeAppend(const GENDRIVER & genDriver)
1963 {
1964   const char * LOC = " FIELD<T>::write(const GENDRIVER &) : ";
1965   BEGIN_OF(LOC);
1966
1967   for (unsigned int index=0; index < _drivers.size(); index++ )
1968     if ( *_drivers[index] == genDriver ) {
1969       _drivers[index]->openAppend();
1970       _drivers[index]->writeAppend();
1971       _drivers[index]->close();
1972     }
1973
1974   END_OF(LOC);
1975
1976 }
1977
1978 /*!
1979   \internal
1980   Read FIELD with the driver which is equal to the given driver.
1981
1982   Use by MED object.
1983 */
1984 template <class T> inline void FIELD<T>::read(const GENDRIVER & genDriver)
1985 {
1986   const char * LOC = " FIELD<T>::read(const GENDRIVER &) : ";
1987   BEGIN_OF(LOC);
1988
1989   for (unsigned int index=0; index < _drivers.size(); index++ )
1990     if ( *_drivers[index] == genDriver ) {
1991       _drivers[index]->open();
1992       _drivers[index]->read();
1993       _drivers[index]->close();
1994     }
1995
1996   END_OF(LOC);
1997
1998 }
1999
2000 /*!
2001   \if developper
2002   Destroy the MEDARRAY<T> in FIELD and put the new one without copy.
2003   \endif
2004 */
2005 template <class T> inline void FIELD<T>::setValue(MEDARRAY<T> *Value)
2006 {
2007   if (NULL != _value) delete _value ;
2008   _value=Value ;
2009 }
2010
2011 /*!
2012   \if developper
2013   Return a reference to  the MEDARRAY<T> in FIELD.
2014   \endif  
2015 */
2016 template <class T> inline MEDARRAY<T>* FIELD<T>::getvalue() const
2017 {
2018   return _value ;
2019 }
2020
2021 /*!
2022   Return the actual length of the reference to values array returned by getValue.
2023 */
2024 template <class T> inline int FIELD<T>::getValueLength(MED_EN::medModeSwitch Mode) const{
2025   return _numberOfComponents*_numberOfValues;
2026 }
2027
2028 /*!
2029   Return a reference to values array to read them.
2030 */
2031 template <class T> inline const T* FIELD<T>::getValue(MED_EN::medModeSwitch Mode) const
2032 {
2033   return _value->get(Mode) ;
2034 }
2035
2036 /*!
2037   Return a reference to i^{th} row or column - component - (depend on Mode value)
2038   of FIELD values array.
2039 */
2040 template <class T> inline const T* FIELD<T>::getValueI(MED_EN::medModeSwitch Mode,int i) const
2041 {
2042  if ( Mode == MED_EN::MED_FULL_INTERLACE )
2043  {
2044          return _value->getRow(i) ;
2045  }
2046  ASSERT (  Mode == MED_EN::MED_NO_INTERLACE);
2047  return _value->getColumn(i);
2048 }
2049
2050 /*!
2051   Return the value of i^{th} element and j^{th} component.
2052 */
2053 template <class T> inline T FIELD<T>::getValueIJ(int i,int j) const
2054 {
2055   return _value->getIJ(i,j) ;
2056 }
2057
2058 /*!
2059   Copy new values array in FIELD according to the given mode.
2060
2061   Array must have right size. If not results are unpredicable.
2062 */
2063 template <class T> inline void FIELD<T>::setValue(MED_EN::medModeSwitch mode, T* value)
2064 {
2065   _value->set(mode,value);
2066 }
2067
2068 /*!
2069   Update values array in FIELD with the given ones according to specified mode.
2070 */
2071 template <class T> inline void FIELD<T>::setValueI(MED_EN::medModeSwitch mode, int i, T* value)
2072 {
2073   // PROVISOIRE :
2074   if (MED_EN::MED_FULL_INTERLACE == mode)
2075     _value->setI(i,value);
2076   else if (MED_EN::MED_NO_INTERLACE == mode)
2077     _value->setJ(i,value);
2078   else
2079     throw MEDEXCEPTION(LOCALIZED("FIELD<T>::setValueI : bad medModeSwitch")) ;
2080 }
2081
2082 /*!
2083   Set the value of i^{th} element and j^{th} component with the given one.
2084 */
2085 template <class T> inline void FIELD<T>::setValueIJ(int i, int j, T value)
2086 {
2087   _value->setIJ(i,j,value);
2088 }
2089
2090 /*
2091   METHODS
2092 */
2093
2094 /*!
2095   Fill values array with volume values.
2096 */
2097 template <class T> void FIELD<T>::getVolume() const throw (MEDEXCEPTION)
2098 {
2099   const char * LOC = "FIELD<double>::getVolume() const : ";
2100   BEGIN_OF(LOC);
2101
2102   // The field has to be initilised by a non empty support and a
2103   // number of components = 1 and its value type has to be set to MED_REEL64
2104   // (ie a FIELD<double>)
2105
2106   if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2107       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to 1 and a value type set to MED_REEL64"));
2108
2109   END_OF(LOC);
2110 }
2111
2112 /*!
2113   Fill values array with area values.
2114 */
2115 template <class T> void FIELD<T>::getArea() const throw (MEDEXCEPTION)
2116 {
2117   const char * LOC = "FIELD<double>::getArea() const : ";
2118   BEGIN_OF(LOC);
2119
2120   // The field has to be initilised by a non empty support and a
2121   // number of components = 1 and its value type has to be set to MED_REEL64
2122   // (ie a FIELD<double>)
2123
2124   if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2125       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to 1 and a value type set to MED_REEL64"));
2126
2127   END_OF(LOC);
2128 }
2129
2130 /*!
2131   Fill values array with length values.
2132 */
2133 template <class T> void FIELD<T>::getLength() const throw (MEDEXCEPTION)
2134 {
2135   const char * LOC = "FIELD<double>::getLength() const : ";
2136   BEGIN_OF(LOC);
2137
2138   // The field has to be initilised by a non empty support and a
2139   // number of components = 1 and its value type has to be set to MED_REEL64
2140   // (ie a FIELD<double>)
2141
2142   if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2143       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to 1 and a value type set to MED_REEL64"));
2144
2145   END_OF(LOC);
2146 }
2147
2148 /*!
2149   Fill values array with normal values.
2150 */
2151 template <class T> void FIELD<T>::getNormal() const throw (MEDEXCEPTION)
2152 {
2153   const char * LOC = "FIELD<double>::getNormal() const : ";
2154   BEGIN_OF(LOC);
2155
2156   // The field has to be initilised by a non empty support and a
2157   // number of components = 1 and its value type has to be set to MED_REEL64
2158   // (ie a FIELD<double>)
2159
2160   if (_support == (SUPPORT *) NULL)
2161       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2162
2163   int dim_space = _support->getMesh()->getSpaceDimension();
2164
2165   if ((_numberOfComponents != dim_space) || (_valueType != MED_REEL64))
2166       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2167
2168   END_OF(LOC);
2169 }
2170
2171 /*!
2172   Fill values array with barycenter values.
2173 */
2174 template <class T> void FIELD<T>::getBarycenter() const throw (MEDEXCEPTION)
2175 {
2176   const char * LOC = "FIELD<double>::getBarycenter() const : ";
2177   BEGIN_OF(LOC);
2178
2179   // The field has to be initilised by a non empty support and a number of
2180   //components = space dimension and its value type has to be set to MED_REEL64
2181   // (ie a FIELD<double>)
2182
2183   if (_support == (SUPPORT *) NULL)
2184       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2185
2186   int dim_space = _support->getMesh()->getSpaceDimension();
2187
2188   if ((_numberOfComponents != dim_space) || (_valueType != MED_REEL64))
2189       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2190
2191   END_OF(LOC);
2192 }
2193
2194 }//End namespace MEDMEM
2195
2196 #endif /* FIELD_HXX */