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