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 public:
 588   FIELD();
 589   FIELD(const FIELD &m);
 590   FIELD & operator=(const FIELD &m);    // A FAIRE
 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
 781 */
 782 template <class T> FIELD<T> & FIELD<T>::operator=(const FIELD &m)
 783 {
 784   MESSAGE("Appel de FIELD<T>::operator=");
 785 }
 786
 787 /*!
 788      Overload addition operator.
 789      This operation is authorized only for compatible fields that have the same support.
 790      The compatibility checking includes equality tests of the folowing data members:/n
 791      - _support
 792      - _numberOfComponents
 793      - _numberOfValues
 794      - _componentsTypes
 795      - _MEDComponentsUnits.
 796
 797      The data members of the returned field are initialized, based on the first field, except for the name,
 798      which is the combination of the two field's names and the operator.
 799      Advised Utilisation in C++ : <tt> FIELD<T> c = a + b; </tt> /n
 800      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
 801      When using python, this operator calls the copy constructor in any case.
 802      The user has to be aware that when using operator + in associatives expressions like
 803      <tt> a = b + c + d +e; </tt> /n
 804      no optimisation is performed : the evaluation of last expression requires the construction of
 805      3 temporary fields.
 806 */
 807 template <class T>
 808 const FIELD<T> FIELD<T>::operator+(const FIELD & m) const
 809 {
 810     BEGIN_OF("FIELD<T>::operator+(const FIELD & m)");
 811     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
 812
 813     // Select mode : avoid if possible any calculation of other mode for fields m or *this
 814     medModeSwitch mode;
 815     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
 816         mode=m.getvalue()->getMode();
 817     else
 818         mode=this->getvalue()->getMode();
 819
 820     // Creation of the result - memory is allocated by FIELD constructor
 821     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
 822     //result._operation(*this,m,mode,"+"); // perform Atribute's initialization & addition
 823     result._operationInitialize(*this,m,"+"); // perform Atribute's initialization
 824     result._add_in_place(*this,m,mode); // perform addition
 825
 826     END_OF("FIELD<T>::operator+(const FIELD & m)");
 827     return result;
 828 }
 829
 830 /*!  Overloaded Operator +=
 831  *   Operations are directly performed in the first field's array.
 832  *   This operation is authorized only for compatible fields that have the same support.
 833  */
 834 template <class T>
 835 FIELD<T>& FIELD<T>::operator+=(const FIELD & m)
 836 {
 837     BEGIN_OF("FIELD<T>::operator+=(const FIELD & m)");
 838     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
 839
 840     // We choose to keep *this mode, even if it may cost a re-calculation for m
 841     medModeSwitch mode=this->getvalue()->getMode();
 842     const T* value1=m.getValue(mode); // get pointers to the values we are adding
 843
 844     // get a non const pointer to the inside array of values and perform operation
 845     T * value=const_cast<T *> (getValue(mode));
 846     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
 847     const T* endV=value+size; // pointer to the end of value
 848     for(;value!=endV; value1++,value++)
 849         *value += *value1;
 850     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
 851     this->getvalue()->clearOtherMode();
 852     END_OF("FIELD<T>::operator+=(const FIELD & m)");
 853     return *this;
 854 }
 855
 856
 857 /*! Addition of fields. Static member function.
 858  *  The function return a pointer to a new created field that holds the addition.
 859  *  Data members are checked for compatibility and initialized.
 860  *  The user is in charge of memory deallocation.
 861  */
 862 template <class T>
 863 FIELD<T>* FIELD<T>::add(const FIELD& m, const FIELD& n)
 864 {
 865     BEGIN_OF("FIELD<T>::add(const FIELD & m, const FIELD& n)");
 866     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
 867
 868     // Select mode : avoid if possible any calculation of other mode for fields m or *this
 869     medModeSwitch mode;
 870     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
 871         mode=m.getvalue()->getMode();
 872     else
 873         mode=n.getvalue()->getMode();
 874
 875     // Creation of a new field
 876     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
 877     result->_operationInitialize(m,n,"+"); // perform Atribute's initialization
 878     result->_add_in_place(m,n,mode); // perform addition
 879
 880     END_OF("FIELD<T>::add(const FIELD & m, const FIELD& n)");
 881     return result;
 882 }
 883
 884 /*!
 885      Overload substraction operator.
 886      This operation is authorized only for compatible fields that have the same support.
 887      The compatibility checking includes equality tests of the folowing data members:/n
 888      - _support
 889      - _numberOfComponents
 890      - _numberOfValues
 891      - _componentsTypes
 892      - _MEDComponentsUnits.
 893
 894      The data members of the returned field are initialized, based on the first field, except for the name,
 895      which is the combination of the two field's names and the operator.
 896      Advised Utilisation in C++ : <tt> FIELD<T> c = a - b; </tt> /n
 897      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
 898      When using python, this operator calls the copy constructor in any case.
 899      The user has to be aware that when using operator - in associatives expressions like
 900      <tt> a = b - c - d -e; </tt> /n
 901      no optimisation is performed : the evaluation of last expression requires the construction of
 902      3 temporary fields.
 903 */
 904 template <class T>
 905 const FIELD<T> FIELD<T>::operator-(const FIELD & m) const
 906 {
 907     BEGIN_OF("FIELD<T>::operator-(const FIELD & m)");
 908     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
 909
 910     // Select mode : avoid if possible any calculation of other mode for fields m or *this
 911     medModeSwitch mode;
 912     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
 913         mode=m.getvalue()->getMode();
 914     else
 915         mode=this->getvalue()->getMode();
 916
 917     // Creation of the result - memory is allocated by FIELD constructor
 918     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
 919     //result._operation(*this,m,mode,"-"); // perform Atribute's initialization & substraction
 920     result._operationInitialize(*this,m,"-"); // perform Atribute's initialization
 921     result._sub_in_place(*this,m,mode); // perform substracion
 922
 923     END_OF("FIELD<T>::operator-(const FIELD & m)");
 924     return result;
 925 }
 926
 927 template <class T>
 928 const FIELD<T> FIELD<T>::operator-() const
 929 {
 930     BEGIN_OF("FIELD<T>::operator-()");
 931
 932     medModeSwitch mode=this->getvalue()->getMode();
 933     // Creation of the result - memory is allocated by FIELD constructor
 934     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
 935     // Atribute's initialization
 936     result.setName("- "+getName());
 937     result.setComponentsNames(getComponentsNames());
 938     // not yet implemented    setComponentType(getComponentType());
 939     result.setComponentsDescriptions(getComponentsDescriptions());
 940     result.setMEDComponentsUnits(getMEDComponentsUnits());
 941     result.setComponentsUnits(getComponentsUnits());
 942     result.setIterationNumber(getIterationNumber());
 943     result.setTime(getTime());
 944     result.setOrderNumber(getOrderNumber());
 945     result.setValueType(getValueType());
 946
 947     const T* value1=getValue(mode);
 948     // get a non const pointer to the inside array of values and perform operation
 949     T * value=const_cast<T *> (result.getValue(mode));
 950     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
 951     const T* endV=value+size; // pointer to the end of value
 952
 953     for(;value!=endV; value1++,value++)
 954         *value = -(*value1);
 955     END_OF("FIELD<T>::operator-=(const FIELD & m)");
 956     return result;
 957 }
 958
 959 /*!  Overloaded Operator -=
 960  *   Operations are directly performed in the first field's array.
 961  *   This operation is authorized only for compatible fields that have the same support.
 962  */
 963 template <class T>
 964 FIELD<T>& FIELD<T>::operator-=(const FIELD & m)
 965 {
 966     BEGIN_OF("FIELD<T>::operator-=(const FIELD & m)");
 967     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
 968
 969     // We choose to keep *this mode, even if it may cost a re-calculation for m
 970     medModeSwitch mode=this->getvalue()->getMode();
 971     const T* value1=m.getValue(mode); // get pointers to the values we are adding
 972
 973     // get a non const pointer to the inside array of values and perform operation
 974     T * value=const_cast<T *> (getValue(mode));
 975     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
 976     const T* endV=value+size; // pointer to the end of value
 977
 978     for(;value!=endV; value1++,value++)
 979         *value -= *value1;
 980     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
 981     this->getvalue()->clearOtherMode();
 982     END_OF("FIELD<T>::operator-=(const FIELD & m)");
 983     return *this;
 984 }
 985
 986
 987 /*! Substraction of fields. Static member function.
 988  *  The function return a pointer to a new created field that holds the substraction.
 989  *  Data members are checked for compatibility and initialized.
 990  *  The user is in charge of memory deallocation.
 991  */
 992 template <class T>
 993 FIELD<T>* FIELD<T>::sub(const FIELD& m, const FIELD& n)
 994 {
 995     BEGIN_OF("FIELD<T>::sub(const FIELD & m, const FIELD& n)");
 996     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
 997
 998     // Select mode : avoid if possible any calculation of other mode for fields m or *this
 999     medModeSwitch mode;
1000     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
1001         mode=m.getvalue()->getMode();
1002     else
1003         mode=n.getvalue()->getMode();
1004
1005     // Creation of a new field
1006     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
1007     result->_operationInitialize(m,n,"-"); // perform Atribute's initialization
1008     result->_sub_in_place(m,n,mode); // perform substraction
1009
1010     END_OF("FIELD<T>::sub(const FIELD & m, const FIELD& n)");
1011     return result;
1012 }
1013
1014 /*!
1015      Overload multiplication operator.
1016      This operation is authorized only for compatible fields that have the same support.
1017      The compatibility checking includes equality tests of the folowing data members:/n
1018      - _support
1019      - _numberOfComponents
1020      - _numberOfValues
1021      - _componentsTypes
1022      - _MEDComponentsUnits.
1023
1024      The data members of the returned field are initialized, based on the first field, except for the name,
1025      which is the combination of the two field's names and the operator.
1026      Advised Utilisation in C++ : <tt> FIELD<T> c = a * b; </tt> /n
1027      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
1028      When using python, this operator calls the copy constructor in any case.
1029      The user has to be aware that when using operator * in associatives expressions like
1030      <tt> a = b * c * d *e; </tt> /n
1031      no optimisation is performed : the evaluation of last expression requires the construction of
1032      3 temporary fields.
1033 */
1034 template <class T>
1035 const FIELD<T> FIELD<T>::operator*(const FIELD & m) const
1036 {
1037     BEGIN_OF("FIELD<T>::operator*(const FIELD & m)");
1038     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1039
1040     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1041     medModeSwitch mode;
1042     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
1043         mode=m.getvalue()->getMode();
1044     else
1045         mode=this->getvalue()->getMode();
1046
1047     // Creation of the result - memory is allocated by FIELD constructor
1048     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
1049     //result._operation(*this,m,mode,"*"); // perform Atribute's initialization & multiplication
1050     result._operationInitialize(*this,m,"*"); // perform Atribute's initialization
1051     result._mul_in_place(*this,m,mode); // perform multiplication
1052
1053     END_OF("FIELD<T>::operator*(const FIELD & m)");
1054     return result;
1055 }
1056
1057 /*!  Overloaded Operator *=
1058  *   Operations are directly performed in the first field's array.
1059  *   This operation is authorized only for compatible fields that have the same support.
1060  */
1061 template <class T>
1062 FIELD<T>& FIELD<T>::operator*=(const FIELD & m)
1063 {
1064     BEGIN_OF("FIELD<T>::operator*=(const FIELD & m)");
1065     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1066
1067     // We choose to keep *this mode, even if it may cost a re-calculation for m
1068     medModeSwitch mode=this->getvalue()->getMode();
1069     const T* value1=m.getValue(mode); // get pointers to the values we are adding
1070
1071     // get a non const pointer to the inside array of values and perform operation
1072     T * value=const_cast<T *> (getValue(mode));
1073     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1074     const T* endV=value+size; // pointer to the end of value
1075
1076     for(;value!=endV; value1++,value++)
1077         *value *= *value1;
1078     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1079     this->getvalue()->clearOtherMode();
1080     END_OF("FIELD<T>::operator*=(const FIELD & m)");
1081     return *this;
1082 }
1083
1084
1085 /*! Multiplication of fields. Static member function.
1086  *  The function return a pointer to a new created field that holds the multiplication.
1087  *  Data members are checked for compatibility and initialized.
1088  *  The user is in charge of memory deallocation.
1089  */
1090 template <class T>
1091 FIELD<T>* FIELD<T>::mul(const FIELD& m, const FIELD& n)
1092 {
1093     BEGIN_OF("FIELD<T>::mul(const FIELD & m, const FIELD& n)");
1094     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
1095
1096     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1097     medModeSwitch mode;
1098     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
1099         mode=m.getvalue()->getMode();
1100     else
1101         mode=n.getvalue()->getMode();
1102
1103     // Creation of a new field
1104     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
1105     result->_operationInitialize(m,n,"*"); // perform Atribute's initialization
1106     result->_mul_in_place(m,n,mode); // perform multiplication
1107
1108     END_OF("FIELD<T>::mul(const FIELD & m, const FIELD& n)");
1109     return result;
1110 }
1111
1112
1113 /*!
1114      Overload division operator.
1115      This operation is authorized only for compatible fields that have the same support.
1116      The compatibility checking includes equality tests of the folowing data members:/n
1117      - _support
1118      - _numberOfComponents
1119      - _numberOfValues
1120      - _componentsTypes
1121      - _MEDComponentsUnits.
1122
1123      The data members of the returned field are initialized, based on the first field, except for the name,
1124      which is the combination of the two field's names and the operator.
1125      Advised Utilisation in C++ : <tt> FIELD<T> c = a / b; </tt> /n
1126      In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
1127      When using python, this operator calls the copy constructor in any case.
1128      The user has to be aware that when using operator / in associatives expressions like
1129      <tt> a = b / c / d /e; </tt> /n
1130      no optimisation is performed : the evaluation of last expression requires the construction of
1131      3 temporary fields.
1132 */
1133 template <class T>
1134 const FIELD<T> FIELD<T>::operator/(const FIELD & m) const
1135 {
1136     BEGIN_OF("FIELD<T>::operator/(const FIELD & m)");
1137     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1138
1139     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1140     medModeSwitch mode;
1141     if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
1142         mode=m.getvalue()->getMode();
1143     else
1144         mode=this->getvalue()->getMode();
1145
1146     // Creation of the result - memory is allocated by FIELD constructor
1147     FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
1148     //result._operation(*this,m,mode,"/"); // perform Atribute's initialization & division
1149     result._operationInitialize(*this,m,"/"); // perform Atribute's initialization
1150     result._div_in_place(*this,m,mode); // perform division
1151
1152     END_OF("FIELD<T>::operator/(const FIELD & m)");
1153     return result;
1154 }
1155
1156
1157 /*!  Overloaded Operator /=
1158  *   Operations are directly performed in the first field's array.
1159  *   This operation is authorized only for compatible fields that have the same support.
1160  */
1161 template <class T>
1162 FIELD<T>& FIELD<T>::operator/=(const FIELD & m)
1163 {
1164     BEGIN_OF("FIELD<T>::operator/=(const FIELD & m)");
1165     FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1166
1167     // We choose to keep *this mode, even if it may cost a re-calculation for m
1168     medModeSwitch mode=this->getvalue()->getMode();
1169     const T* value1=m.getValue(mode); // get pointers to the values we are adding
1170
1171     // get a non const pointer to the inside array of values and perform operation
1172     T * value=const_cast<T *> (getValue(mode));
1173     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1174     const T* endV=value+size; // pointer to the end of value
1175
1176     for(;value!=endV; value1++,value++)
1177         *value /= *value1;
1178     // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1179     this->getvalue()->clearOtherMode();
1180     END_OF("FIELD<T>::operator/=(const FIELD & m)");
1181     return *this;
1182 }
1183
1184
1185 /*! Division of fields. Static member function.
1186  *  The function return a pointer to a new created field that holds the division.
1187  *  Data members are checked for compatibility and initialized.
1188  *  The user is in charge of memory deallocation.
1189  */
1190 template <class T>
1191 FIELD<T>* FIELD<T>::div(const FIELD& m, const FIELD& n)
1192 {
1193     BEGIN_OF("FIELD<T>::div(const FIELD & m, const FIELD& n)");
1194     FIELD_::_checkFieldCompatibility(m, n); // may throw exception
1195
1196     // Select mode : avoid if possible any calculation of other mode for fields m or *this
1197     medModeSwitch mode;
1198     if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
1199         mode=m.getvalue()->getMode();
1200     else
1201         mode=n.getvalue()->getMode();
1202
1203     // Creation of a new field
1204     FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
1205     result->_operationInitialize(m,n,"/"); // perform Atribute's initialization
1206     result->_div_in_place(m,n,mode); // perform division
1207
1208     END_OF("FIELD<T>::div(const FIELD & m, const FIELD& n)");
1209     return result;
1210 }
1211
1212
1213 /*!
1214   \if developper
1215   This internal method initialize the members of a new field created to hold the result of the operation Op .
1216   Initialization is based on the first field, except for the name, which is the combination of the two field's names
1217   and the operator.
1218   \endif
1219 */
1220 template <class T>
1221 void FIELD<T>::_operationInitialize(const FIELD& m,const FIELD& n, char* Op)
1222 {
1223     MESSAGE("Appel methode interne _add" << Op);
1224
1225     // Atribute's initialization (copy of the first field's attributes)
1226     // Other data members (_support, _numberOfValues) are initialized in the field's constr.
1227     setName(m.getName()+" "+Op+" "+n.getName());
1228     setComponentsNames(m.getComponentsNames());
1229     // not yet implemented    setComponentType(m.getComponentType());
1230     setComponentsDescriptions(m.getComponentsDescriptions());
1231     setMEDComponentsUnits(m.getMEDComponentsUnits());
1232
1233     // The following data member may differ from field m to n.
1234     // The initialization is done based on the first field.
1235     setComponentsUnits(m.getComponentsUnits());
1236     setIterationNumber(m.getIterationNumber());
1237     setTime(m.getTime());
1238     setOrderNumber(m.getOrderNumber());
1239     setValueType(m.getValueType());
1240 }
1241
1242
1243 /*!
1244   \if developper
1245   Internal method called by FIELD<T>::operator+ and FIELD<T>::add to perform addition "in place".
1246   This method is applied to a just created field with medModeSwitch mode.
1247   For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1248   it doesn't exist!
1249   \endif
1250 */
1251 template <class T>
1252 void FIELD<T>::_add_in_place(const FIELD& m,const FIELD& n, const medModeSwitch mode)
1253 {
1254     // get pointers to the values we are adding
1255     const T* value1=m.getValue(mode);
1256     const T* value2=n.getValue(mode);
1257     // get a non const pointer to the inside array of values and perform operation
1258     T * value=const_cast<T *> (getValue(mode));
1259
1260     const int size=getNumberOfValues()*getNumberOfComponents();
1261     SCRUTE(size);
1262     const T* endV1=value1+size;
1263     for(;value1!=endV1; value1++,value2++,value++)
1264         *value=(*value1)+(*value2);
1265 }
1266
1267 /*!
1268   \if developper
1269   Internal method called by FIELD<T>::operator- and FIELD<T>::sub to perform substraction "in place".
1270   This method is applied to a just created field with medModeSwitch mode.
1271   For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1272   it doesn't exist!
1273   \endif
1274 */
1275 template <class T>
1276 void FIELD<T>::_sub_in_place(const FIELD& m,const FIELD& n, const medModeSwitch mode)
1277 {
1278     // get pointers to the values we are adding
1279     const T* value1=m.getValue(mode);
1280     const T* value2=n.getValue(mode);
1281     // get a non const pointer to the inside array of values and perform operation
1282     T * value=const_cast<T *> (getValue(mode));
1283
1284     const int size=getNumberOfValues()*getNumberOfComponents();
1285     SCRUTE(size);
1286     const T* endV1=value1+size;
1287     for(;value1!=endV1; value1++,value2++,value++)
1288         *value=(*value1)-(*value2);
1289 }
1290
1291 /*!
1292   \if developper
1293   Internal method called by FIELD<T>::operator* and FIELD<T>::mul to perform multiplication "in place".
1294   This method is applied to a just created field with medModeSwitch mode.
1295   For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1296   it doesn't exist!
1297   \endif
1298 */
1299 template <class T>
1300 void FIELD<T>::_mul_in_place(const FIELD& m,const FIELD& n, const medModeSwitch mode)
1301 {
1302     // get pointers to the values we are adding
1303     const T* value1=m.getValue(mode);
1304     const T* value2=n.getValue(mode);
1305     // get a non const pointer to the inside array of values and perform operation
1306     T * value=const_cast<T *> (getValue(mode));
1307
1308     const int size=getNumberOfValues()*getNumberOfComponents();
1309     SCRUTE(size);
1310     const T* endV1=value1+size;
1311     for(;value1!=endV1; value1++,value2++,value++)
1312         *value=(*value1)*(*value2);
1313 }
1314
1315 /*!
1316   \if developper
1317   Internal method called by FIELD<T>::operator/ and FIELD<T>::div to perform division "in place".
1318   This method is applied to a just created field with medModeSwitch mode.
1319   For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1320   it doesn't exist!
1321   \endif
1322 */
1323 template <class T>
1324 void FIELD<T>::_div_in_place(const FIELD& m,const FIELD& n, const medModeSwitch mode)
1325 {
1326     // get pointers to the values we are adding
1327     const T* value1=m.getValue(mode);
1328     const T* value2=n.getValue(mode);
1329     // get a non const pointer to the inside array of values and perform operation
1330     T * value=const_cast<T *> (getValue(mode));
1331
1332     const int size=getNumberOfValues()*getNumberOfComponents();
1333     SCRUTE(size);
1334     const T* endV1=value1+size;
1335     for(;value1!=endV1; value1++,value2++,value++)
1336         *value=(*value1)/(*value2);
1337 }
1338
1339 /*!  Return Max Norm
1340  */
1341 template <class T> double FIELD<T>::normMax() const throw (MEDEXCEPTION)
1342 {
1343     const T* value=getValue(getvalue()->getMode()); // get pointer to the values
1344     const int size=getNumberOfValues()*getNumberOfComponents();
1345     if (size <= 0) // Size of array has to be strictly positive
1346     {
1347         string diagnosis;
1348         diagnosis="FIELD<T>::normMax() : cannot compute the norm of "+getName()+
1349             " : it size is non positive!";
1350         throw MEDEXCEPTION(diagnosis.c_str());
1351     }
1352     const T* lastvalue=value+size; // get pointer just after last value
1353     const T* pMax=value; // pointer to the max value
1354     const T* pMin=value; // pointer to the min value
1355
1356     // get pointers to the max & min value of array
1357     while ( ++value != lastvalue )
1358     {
1359         if ( *pMin > *value )
1360             pMin=value;
1361         if ( *pMax < *value )
1362             pMax=value;
1363     }
1364
1365     T Max= *pMax>(T) 0 ? *pMax : -*pMax; // Max=abs(*pMax)
1366     T Min= *pMin>(T) 0 ? *pMin : -*pMin; // Min=abs(*pMin)
1367
1368     return Max>Min ? static_cast<double>(Max) : static_cast<double>(Min);
1369 }
1370
1371 /*!  Return Euclidien norm
1372  */
1373 template <class T> double FIELD<T>::norm2() const throw (MEDEXCEPTION)
1374 {
1375     const T* value=this->getValue(this->getvalue()->getMode()); // get const pointer to the values
1376     const int size=getNumberOfValues()*getNumberOfComponents(); // get size of array
1377     if (size <= 0) // Size of array has to be strictly positive
1378     {
1379         string diagnosis;
1380         diagnosis="FIELD<T>::norm2() : cannot compute the norm of "+getName()+
1381             " : it size is non positive!";
1382         throw MEDEXCEPTION(diagnosis.c_str());
1383     }
1384     const T* lastvalue=value+size; // point just after last value
1385
1386     T result((T)0); // init
1387     for( ; value!=lastvalue ; ++value)
1388         result += (*value) * (*value);
1389
1390     return std::sqrt(static_cast<double> (result));
1391 }
1392
1393
1394 /*!  Apply to each (scalar) field component the template parameter T_function,
1395  *   which is a pointer to function.
1396  *   Since the pointer is known at compile time, the function is inlined into the inner loop!
1397  *   calculation is done "in place".
1398  *   Use examples :
1399  *
1400  *   \code  myField.applyFunc<std::sqrt>();  // apply sqare root function \endcode
1401  *     \code myField.applyFunc<myFunction>(); // apply your own created function \endcode
1402  */
1403 template <class T> template <T T_function(T)>
1404 void FIELD<T>::applyFunc()
1405 {
1406     medModeSwitch mode=getvalue()->getMode();
1407
1408     // get a non const pointer to the inside array of values and perform operation
1409     T * value=const_cast<T *> (getValue(mode));
1410     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1411
1412     if (size>0) // for a negative size, there is nothing to do
1413     {
1414         const T* lastvalue=value+size; // pointer to the end of value
1415         for(;value!=lastvalue; ++value) // apply linear transformation
1416             *value = T_function(*value);
1417         // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1418         getvalue()->clearOtherMode();
1419     }
1420 }
1421
1422
1423 /*!  Apply to each (scalar) field component the linear function x -> ax+b.
1424  *   calculation is done "in place".
1425  */
1426 template <class T> void FIELD<T>::applyLin(T a, T b)
1427 {
1428     medModeSwitch mode=getvalue()->getMode();
1429
1430     // get a non const pointer to the inside array of values and perform operation in place
1431     T * value=const_cast<T *> (getValue(mode));
1432     const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1433
1434     if (size>0) // for a negative size, there is nothing to do
1435     {
1436         const T* lastvalue=value+size; // pointer to the end of value
1437         for(;value!=lastvalue; ++value) // apply linear transformation
1438             *value = a*(*value)+b;
1439         // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1440         getvalue()->clearOtherMode();
1441     }
1442 }
1443
1444
1445 /*!
1446  *   Return a pointer to a new field that holds the scalar product. Static member function.
1447  *   This operation is authorized only for compatible fields that have the same support.
1448  *   The compatibility checking includes equality tests of the folowing data members:/n
1449  *   - _support
1450  *   - _numberOfComponents
1451  *   - _numberOfValues
1452  *   - _componentsTypes
1453  *   - _MEDComponentsUnits.
1454  *   Data members are initialized.
1455  *   The new field point to the same support and has one component.
1456  *   Each value of it is the scalar product of the two argument's fields.
1457  *   The user is in charge of memory deallocation.
1458  */
1459 template <class T> FIELD<T>* FIELD<T>::scalarProduct(const FIELD & m, const FIELD & n)
1460 {
1461     FIELD_::_checkFieldCompatibility( m, n); // may throw exception
1462     // we need a MED_FULL_INTERLACE representation of m & n to compute the scalar product
1463     const medModeSwitch mode=MED_FULL_INTERLACE;
1464
1465     const int numberOfElements=m.getNumberOfValues(); // strictly positive
1466     const int NumberOfComponents=m.getNumberOfComponents(); // strictly positive
1467
1468     // Creation & init of a the result field on the same suppot, with one component
1469     FIELD<T>* result = new FIELD<T>(m.getSupport(),1,mode);
1470     result->setName( "scalarProduct ( " + m.getName() + " , " + n.getName() + " )" );
1471     result->setIterationNumber(m.getIterationNumber());
1472     result->setTime(m.getTime());
1473     result->setOrderNumber(m.getOrderNumber());
1474     result->setValueType(m.getValueType());
1475
1476     const T* value1=m.getValue(mode); // get const pointer to the values
1477     const T* value2=n.getValue(mode); // get const pointer to the values
1478     // get a non const pointer to the inside array of values and perform operation
1479     T * value=const_cast<T *> (result->getValue(mode));
1480
1481     const T* lastvalue=value+numberOfElements; // pointing just after last value of result
1482     for ( ; value!=lastvalue ; ++value ) // loop on all elements
1483     {
1484         *value=(T)0; // initialize value
1485         const T* endofRow=value1+NumberOfComponents; // pointing just after end of row
1486         for ( ; value1 != endofRow; ++value1, ++value2) // computation of dot product
1487             *value += (*value1) * (*value2);
1488     }
1489     return result;
1490 }
1491
1492 /*!  Return L2 Norm  of the field's component.
1493  *   Cannot be applied to a field with a support on nodes.
1494  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1495  */
1496 template <class T> double FIELD<T>::normL2(int component, const FIELD<double> * p_field_volume) const
1497 {
1498     _checkNormCompatibility(p_field_volume); // may throw exception
1499     if ( component<1 || component>getNumberOfComponents() )
1500         throw MEDEXCEPTION(STRING("FIELD<T>::normL2() : The component argument should be between 1 and the number of components"));
1501
1502     const FIELD<double> * p_field_size=p_field_volume;
1503     if(!p_field_volume) // if the user don't supply the volume
1504         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1505
1506     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1507     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1508     const T* value=getValueI( MED_NO_INTERLACE, component); // get pointer to the component's values
1509     const T* lastvalue=value+getNumberOfValues(); // pointing just after the end of column
1510
1511     double integrale=0.0;
1512     double totVol=0.0;
1513     for (; value!=lastvalue ; ++value ,++vol)
1514     {
1515         integrale += static_cast<double>((*value) * (*value)) * (*vol);
1516         totVol+=*vol;
1517     }
1518
1519     if(!p_field_volume) // if the user didn't supply the volume
1520         delete p_field_size; // delete temporary volume field
1521     if( totVol <= 0)
1522         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1523
1524     return integrale/totVol;
1525 }
1526
1527 /*!  Return L2 Norm  of the field.
1528  *   Cannot be applied to a field with a support on nodes.
1529  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1530  */
1531 template <class T> double FIELD<T>::normL2(const FIELD<double> * p_field_volume) const
1532 {
1533     _checkNormCompatibility(p_field_volume); // may throw exception
1534     const FIELD<double> * p_field_size=p_field_volume;
1535     if(!p_field_volume) // if the user don't supply the volume
1536         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1537
1538     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1539     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1540     const double* lastvol=vol+getNumberOfValues(); // pointing just after the end of vol
1541     const T* value=getValue( MED_NO_INTERLACE); // get pointer to the field's values
1542
1543     double totVol=0.0;
1544     const double* p_vol=vol;
1545     for (p_vol=vol; p_vol!=lastvol ; ++p_vol) // calculate total volume
1546         totVol+=*p_vol;
1547
1548     double integrale=0.0;
1549     for (int i=1; i<=getNumberOfComponents(); ++i) // compute integral on all components
1550         for (p_vol=vol; p_vol!=lastvol ; ++value ,++p_vol)
1551             integrale += static_cast<double>((*value) * (*value)) * (*p_vol);
1552
1553     if(!p_field_volume) // if the user didn't supply the volume
1554         delete p_field_size; // delete temporary volume field
1555     if( totVol <= 0)
1556         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1557
1558     return integrale/totVol;
1559 }
1560
1561 /*!  Return L1 Norm  of the field's component.
1562  *   Cannot be applied to a field with a support on nodes.
1563  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1564  */
1565 template <class T> double FIELD<T>::normL1(int component, const FIELD<double> * p_field_volume) const
1566 {
1567     _checkNormCompatibility(p_field_volume); // may throw exception
1568     if ( component<1 || component>getNumberOfComponents() )
1569         throw MEDEXCEPTION(STRING("FIELD<T>::normL2() : The component argument should be between 1 and the number of components"));
1570
1571     const FIELD<double> * p_field_size=p_field_volume;
1572     if(!p_field_volume) // if the user don't supply the volume
1573         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1574
1575     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1576     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1577     const T* value=getValueI( MED_NO_INTERLACE, component); // get pointer to the component's values
1578     const T* lastvalue=value+getNumberOfValues(); // pointing just after the end of column
1579
1580     double integrale=0.0;
1581     double totVol=0.0;
1582     for (; value!=lastvalue ; ++value ,++vol)
1583     {
1584         integrale += std::abs( static_cast<double>(*value) ) * (*vol);
1585         totVol+=*vol;
1586     }
1587
1588     if(!p_field_volume) // if the user didn't supply the volume
1589         delete p_field_size; // delete temporary volume field
1590     if( totVol <= 0)
1591         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1592
1593     return integrale/totVol;
1594 }
1595
1596 /*!  Return L1 Norm  of the field.
1597  *   Cannot be applied to a field with a support on nodes.
1598  *   If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1599  */
1600 template <class T> double FIELD<T>::normL1(const FIELD<double> * p_field_volume) const
1601 {
1602     _checkNormCompatibility(p_field_volume); // may throw exception
1603     const FIELD<double> * p_field_size=p_field_volume;
1604     if(!p_field_volume) // if the user don't supply the volume
1605         p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1606
1607     // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1608     const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1609     const double* lastvol=vol+getNumberOfValues(); // pointing just after the end of vol
1610     const T* value=getValue( MED_NO_INTERLACE); // get pointer to the field's values
1611
1612     double totVol=0.0;
1613     const double* p_vol=vol;
1614     for (p_vol=vol; p_vol!=lastvol ; ++p_vol) // calculate total volume
1615         totVol+=*p_vol;
1616
1617     double integrale=0.0;
1618     for (int i=1; i<=getNumberOfComponents(); ++i) // compute integral on all components
1619         for (p_vol=vol; p_vol!=lastvol ; ++value ,++p_vol)
1620             integrale += std::abs( static_cast<double>(*value) ) * (*p_vol);
1621
1622     if(!p_field_volume) // if the user didn't supply the volume
1623         delete p_field_size; // delete temporary volume field
1624     if( totVol <= 0)
1625         throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1626
1627     return integrale/totVol;
1628 }
1629
1630
1631
1632
1633 /*!
1634   Constructor with parameters; the object is set via a file and its associated
1635   driver. For the moment only the MED_DRIVER is considered and if the last two
1636   argument (iterationNumber and orderNumber) are not set; their default value
1637   is -1. If the field fieldDriverName with the iteration number
1638   iterationNumber and the order number orderNumber does not exist in the file
1639   fieldDriverName; the constructor raises an exception.
1640 */
1641 template <class T> FIELD<T>::FIELD(const SUPPORT * Support,
1642                                    driverTypes driverType,
1643                                    const string & fileName/*=""*/,
1644                                    const string & fieldDriverName/*=""*/,
1645                                    const int iterationNumber,
1646                                    const int orderNumber)
1647   throw (MEDEXCEPTION)
1648 {
1649   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) : ";
1650
1651   int current;
1652
1653   BEGIN_OF(LOC);
1654
1655   init();
1656
1657   _support = Support;
1658   _value = (MEDARRAY<T>*)NULL;
1659
1660   _iterationNumber = iterationNumber;
1661   _time = 0.0;
1662   _orderNumber = orderNumber;
1663
1664   switch(driverType)
1665     {
1666     case MED_DRIVER :
1667       {
1668         MED_FIELD_RDONLY_DRIVER<T> myDriver(fileName,this);
1669         myDriver.setFieldName(fieldDriverName);
1670         current = addDriver(myDriver);
1671         break;
1672       }
1673 //   current = addDriver(driverType,fileName,fieldDriverName);
1674 //   switch(_drivers[current]->getAccessMode() ) {
1675 //   case MED_WRONLY : {
1676 //     MESSAGE("FIELD<T>::FIELD(driverTypes driverType, .....) : driverType must have a MED_RDONLY or MED_RDWR accessMode");
1677 //     rmDriver(current);
1678 //     break;}
1679 //   default : {
1680 //   }
1681 //   }
1682     default :
1683       {
1684         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Driver type unknown : can't create driver!"));
1685         break;
1686       }
1687     }
1688
1689   _drivers[current]->open();
1690   _drivers[current]->read();
1691   _drivers[current]->close();
1692
1693   END_OF(LOC);
1694 }
1695
1696 /*!
1697   Destructor.
1698 */
1699 template <class T> FIELD<T>::~FIELD()
1700 {
1701   BEGIN_OF(" Destructeur FIELD<T>::~FIELD()");
1702   SCRUTE(this);
1703   if (_value) delete _value;
1704   END_OF(" Destructeur FIELD<T>::~FIELD()");
1705 }
1706
1707 /*!
1708
1709 */
1710 template <class T> void FIELD<T>::allocValue(const int NumberOfComponents)
1711 {
1712   const char* LOC = "FIELD<T>::allocValue(const int NumberOfComponents)" ;
1713   BEGIN_OF(LOC);
1714
1715   _numberOfComponents = NumberOfComponents ;
1716   if (_componentsTypes == NULL)
1717     _componentsTypes = new int[NumberOfComponents] ;
1718   if (_componentsNames == NULL)
1719     _componentsNames = new string[NumberOfComponents];
1720   if (_componentsDescriptions == NULL)
1721     _componentsDescriptions = new string[NumberOfComponents];
1722   if (_componentsUnits == NULL)
1723     _componentsUnits = new UNIT[NumberOfComponents];
1724   if (_MEDComponentsUnits == NULL)
1725     _MEDComponentsUnits = new string[NumberOfComponents];
1726   for (int i=0;i<NumberOfComponents;i++) {
1727     _componentsTypes[i] = 0 ;
1728   }
1729
1730   try {
1731     _numberOfValues = _support->getNumberOfElements(MED_ALL_ELEMENTS);
1732     MESSAGE(LOC <<" : "<<_numberOfValues <<" et "<< NumberOfComponents);
1733
1734     _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues);
1735
1736     _isRead = true ;
1737   }
1738   catch (MEDEXCEPTION &ex) {
1739     MESSAGE("No value defined, problem with NumberOfComponents (and may be _support) size of MEDARRAY<T>::_value !");
1740     _value = (MEDARRAY<T>*)NULL ;
1741   }
1742
1743   SCRUTE(_value);
1744   END_OF("void FIELD<T>::allocValue(const int NumberOfComponents)");
1745 }
1746
1747 /*!
1748
1749 */
1750 template <class T> void FIELD<T>::allocValue(const int NumberOfComponents, const int LengthValue)
1751 {
1752   BEGIN_OF("void FIELD<T>::allocValue(const int NumberOfComponents,const int LengthValue)");
1753
1754   _numberOfComponents = NumberOfComponents ;
1755   if (_componentsTypes == NULL)
1756     _componentsTypes = new int[NumberOfComponents] ;
1757   if (_componentsNames == NULL)
1758     _componentsNames = new string[NumberOfComponents];
1759   if (_componentsDescriptions == NULL)
1760     _componentsDescriptions = new string[NumberOfComponents];
1761   if (_componentsUnits == NULL)
1762     _componentsUnits = new UNIT[NumberOfComponents];
1763   if (_MEDComponentsUnits == NULL)
1764     _MEDComponentsUnits = new string[NumberOfComponents];
1765   for (int i=0;i<NumberOfComponents;i++) {
1766     _componentsTypes[i] = 0 ;
1767   }
1768
1769   MESSAGE("FIELD : constructeur : "<<LengthValue <<" et "<< NumberOfComponents);
1770   _numberOfValues = LengthValue ;
1771   _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues);
1772   _isRead = true ;
1773
1774   SCRUTE(_value);
1775   END_OF("void FIELD<T>::allocValue(const int NumberOfComponents,const int LengthValue)");
1776 }
1777
1778 /*!
1779
1780 */
1781 template <class T> void FIELD<T>::deallocValue()
1782 {
1783   BEGIN_OF("void FIELD<T>::deallocValue()");
1784   _numberOfValues = 0 ;
1785   _numberOfComponents = 0 ;
1786   if (_value != NULL)
1787     delete _value;
1788
1789   END_OF("void FIELD<T>::deallocValue()");
1790 }
1791
1792 // -----------------
1793 // Methodes Inline
1794 // -----------------
1795
1796
1797 template <class T>       FIELD<T>::INSTANCE_DE<MED_FIELD_RDWR_DRIVER<T> >     FIELD<T>::inst_med  ;
1798
1799 template <class T>       FIELD<T>::INSTANCE_DE<VTK_FIELD_DRIVER<T> >     FIELD<T>::inst_vtk  ;
1800
1801 template <class T> const typename FIELD<T>::INSTANCE * const FIELD<T>::instances[] = { &FIELD<T>::inst_med,  &FIELD<T>::inst_vtk} ;
1802
1803
1804 /*!
1805   Create the specified driver and return its index reference to path to
1806   read or write methods.
1807 */
1808
1809 template <class T> int FIELD<T>::addDriver(driverTypes driverType,
1810                                            const string & fileName/*="Default File Name.med"*/,
1811                                            const string & driverName/*="Default Field Name"*/)
1812 {
1813   const char * LOC = "FIELD<T>::addDriver(driverTypes driverType, const string & fileName=\"Default File Name.med\",const string & driverName=\"Default Field Name\") : ";
1814
1815   GENDRIVER * driver;
1816   int current;
1817   int itDriver = (int) NO_DRIVER;
1818
1819   BEGIN_OF(LOC);
1820
1821   SCRUTE(driverType);
1822
1823   SCRUTE(instances[driverType]);
1824
1825   switch(driverType)
1826     {
1827     case MED_DRIVER : {
1828       itDriver = (int) driverType ;
1829       break ;
1830     }
1831
1832     case VTK_DRIVER : {
1833       itDriver = 1 ;
1834       break ;
1835     }
1836
1837     case GIBI_DRIVER : {
1838       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"));
1839       break;
1840     }
1841
1842     case NO_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
1848   if (itDriver == ((int) NO_DRIVER))
1849     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"));
1850
1851   driver = instances[itDriver]->run(fileName, this) ;
1852
1853   _drivers.push_back(driver);
1854
1855   current = _drivers.size()-1;
1856
1857   _drivers[current]->setFieldName(driverName);
1858
1859   END_OF(LOC);
1860
1861   return current;
1862 }
1863
1864
1865 /*!
1866   Duplicate the given driver and return its index reference to path to
1867   read or write methods.
1868 */
1869 template <class T> inline int FIELD<T>::addDriver (GENDRIVER & driver )
1870 {
1871   const char * LOC = "FIELD<T>::addDriver(GENDRIVER &) : ";
1872   BEGIN_OF(LOC);
1873
1874   // duplicate driver to delete it with destructor !
1875   GENDRIVER * newDriver = driver.copy() ;
1876
1877   _drivers.push_back(newDriver);
1878   return _drivers.size() -1 ;
1879
1880   END_OF(LOC);
1881 };
1882
1883 /*!
1884   Remove the driver referenced by its index.
1885 */
1886 template <class T> void FIELD<T>::rmDriver (int index/*=0*/)
1887 {
1888   const char * LOC = "FIELD<T>::rmDriver (int index=0): ";
1889   BEGIN_OF(LOC);
1890
1891   if ( _drivers[index] ) {
1892     //_drivers.erase(&_drivers[index]);
1893     // why not ????
1894     MESSAGE ("detruire");
1895   }
1896   else
1897     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1898                                      << "The <index given is invalid, index must be between  0 and  |"
1899                                      << _drivers.size()
1900                                      )
1901                           );
1902
1903   END_OF(LOC);
1904 }
1905
1906 /*!
1907   Read FIELD in the file specified in the driver given by its index.
1908 */
1909 template <class T> inline  void FIELD<T>::read(int index/*=0*/)
1910 {
1911   const char * LOC = "FIELD<T>::read(int index=0) : ";
1912   BEGIN_OF(LOC);
1913
1914   if ( _drivers[index] ) {
1915     _drivers[index]->open();
1916     _drivers[index]->read();
1917     _drivers[index]->close();
1918   }
1919   else
1920     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1921                                      << "The index given is invalid, index must be between  0 and |"
1922                                      << _drivers.size()
1923                                      )
1924                           );
1925   END_OF(LOC);
1926 }
1927
1928 /*!
1929   Write FIELD in the file specified in the driver given by its index.
1930 */
1931 template <class T> inline void FIELD<T>::write(int index/*=0*/, const string & driverName /*= ""*/)
1932 {
1933   const char * LOC = "FIELD<T>::write(int index=0, const string & driverName = \"\") : ";
1934   BEGIN_OF(LOC);
1935
1936   if( _drivers[index] ) {
1937     _drivers[index]->open();
1938     if (driverName != "") _drivers[index]->setFieldName(driverName);
1939     _drivers[index]->write();
1940     _drivers[index]->close();
1941   }
1942   else
1943     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1944                                      << "The index given is invalid, index must be between  0 and |"
1945                                      << _drivers.size()
1946                                      )
1947                           );
1948   END_OF(LOC);
1949 }
1950
1951 /*!
1952   Write FIELD in the file specified in the driver given by its index. Use this
1953   method for ASCII drivers (e.g. VTK_DRIVER)
1954 */
1955 template <class T> inline void FIELD<T>::writeAppend(int index/*=0*/, const string & driverName /*= ""*/)
1956 {
1957   const char * LOC = "FIELD<T>::write(int index=0, const string & driverName = \"\") : ";
1958   BEGIN_OF(LOC);
1959
1960   if( _drivers[index] ) {
1961     _drivers[index]->openAppend();
1962     if (driverName != "") _drivers[index]->setFieldName(driverName);
1963     _drivers[index]->writeAppend();
1964     _drivers[index]->close();
1965   }
1966   else
1967     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1968                                      << "The index given is invalid, index must be between  0 and |"
1969                                      << _drivers.size()
1970                                      )
1971                           );
1972   END_OF(LOC);
1973 }
1974
1975 /*!
1976   \internal
1977   Write FIELD with the driver which is equal to the given driver.
1978
1979   Use by MED object.
1980 */
1981 template <class T> inline void FIELD<T>::write(const GENDRIVER & genDriver)
1982 {
1983   const char * LOC = " FIELD<T>::write(const GENDRIVER &) : ";
1984   BEGIN_OF(LOC);
1985
1986   for (unsigned int index=0; index < _drivers.size(); index++ )
1987     if ( *_drivers[index] == genDriver ) {
1988       _drivers[index]->open();
1989       _drivers[index]->write();
1990       _drivers[index]->close();
1991     }
1992
1993   END_OF(LOC);
1994
1995 }
1996
1997 /*!
1998   \internal
1999   Write FIELD with the driver which is equal to the given driver.
2000
2001   Use by MED object. Use this method for ASCII drivers (e.g. VTK_DRIVER).
2002 */
2003 template <class T> inline void FIELD<T>::writeAppend(const GENDRIVER & genDriver)
2004 {
2005   const char * LOC = " FIELD<T>::write(const GENDRIVER &) : ";
2006   BEGIN_OF(LOC);
2007
2008   for (unsigned int index=0; index < _drivers.size(); index++ )
2009     if ( *_drivers[index] == genDriver ) {
2010       _drivers[index]->openAppend();
2011       _drivers[index]->writeAppend();
2012       _drivers[index]->close();
2013     }
2014
2015   END_OF(LOC);
2016
2017 }
2018
2019 /*!
2020   \internal
2021   Read FIELD with the driver which is equal to the given driver.
2022
2023   Use by MED object.
2024 */
2025 template <class T> inline void FIELD<T>::read(const GENDRIVER & genDriver)
2026 {
2027   const char * LOC = " FIELD<T>::read(const GENDRIVER &) : ";
2028   BEGIN_OF(LOC);
2029
2030   for (unsigned int index=0; index < _drivers.size(); index++ )
2031     if ( *_drivers[index] == genDriver ) {
2032       _drivers[index]->open();
2033       _drivers[index]->read();
2034       _drivers[index]->close();
2035     }
2036
2037   END_OF(LOC);
2038
2039 }
2040
2041 /*!
2042   \if developper
2043   Destroy the MEDARRAY<T> in FIELD and put the new one without copy.
2044   \endif
2045 */
2046 template <class T> inline void FIELD<T>::setValue(MEDARRAY<T> *Value)
2047 {
2048   if (NULL != _value) delete _value ;
2049   _value=Value ;
2050 }
2051
2052 /*!
2053   \if developper
2054   Return a reference to  the MEDARRAY<T> in FIELD.
2055   \endif
2056 */
2057 template <class T> inline MEDARRAY<T>* FIELD<T>::getvalue() const
2058 {
2059   return _value ;
2060 }
2061
2062 /*!
2063   Return a reference to values array to read them.
2064 */
2065 template <class T> inline const T* FIELD<T>::getValue(medModeSwitch Mode) const
2066 {
2067   return _value->get(Mode) ;
2068 }
2069
2070 /*!
2071   Return a reference to i^{th} row or column - component - (depend on Mode value)
2072   of FIELD values array.
2073 */
2074 template <class T> inline const T* FIELD<T>::getValueI(medModeSwitch Mode,int i) const
2075 {
2076  if ( Mode == MED_FULL_INTERLACE )
2077  {
2078          return _value->getRow(i) ;
2079  }
2080  ASSERT (  Mode == MED_NO_INTERLACE);
2081  return _value->getColumn(i);
2082 }
2083
2084 /*!
2085   Return the value of i^{th} element and j^{th} component.
2086 */
2087 template <class T> inline T FIELD<T>::getValueIJ(int i,int j) const
2088 {
2089   return _value->getIJ(i,j) ;
2090 }
2091
2092 /*!
2093   Copy new values array in FIELD according to the given mode.
2094
2095   Array must have right size. If not results are unpredicable.
2096 */
2097 template <class T> inline void FIELD<T>::setValue(medModeSwitch mode, T* value)
2098 {
2099   _value->set(mode,value);
2100 }
2101
2102 /*!
2103   Update values array in FIELD with the given ones according to specified mode.
2104 */
2105 template <class T> inline void FIELD<T>::setValueI(medModeSwitch mode, int i, T* value)
2106 {
2107   // PROVISOIRE :
2108   if (MED_FULL_INTERLACE == mode)
2109     _value->setI(i,value);
2110   else if (MED_NO_INTERLACE == mode)
2111     _value->setJ(i,value);
2112   else
2113     throw MEDEXCEPTION(LOCALIZED("FIELD<T>::setValueI : bad medModeSwitch")) ;
2114 }
2115
2116 /*!
2117   Set the value of i^{th} element and j^{th} component with the given one.
2118 */
2119 template <class T> inline void FIELD<T>::setValueIJ(int i, int j, T value)
2120 {
2121   _value->setIJ(i,j,value);
2122 }
2123
2124 /*
2125   METHODS
2126 */
2127
2128 /*!
2129   Fill values array with volume values.
2130 */
2131 template <class T> void FIELD<T>::getVolume() const throw (MEDEXCEPTION)
2132 {
2133   const char * LOC = "FIELD<double>::getVolume() const : ";
2134   BEGIN_OF(LOC);
2135
2136   // The field has to be initilised by a non empty support and a
2137   // number of components = 1 and its value type has to be set to MED_REEL64
2138   // (ie a FIELD<double>)
2139
2140   if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2141       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"));
2142
2143   END_OF(LOC);
2144 }
2145
2146 /*!
2147   Fill values array with area values.
2148 */
2149 template <class T> void FIELD<T>::getArea() const throw (MEDEXCEPTION)
2150 {
2151   const char * LOC = "FIELD<double>::getArea() const : ";
2152   BEGIN_OF(LOC);
2153
2154   // The field has to be initilised by a non empty support and a
2155   // number of components = 1 and its value type has to be set to MED_REEL64
2156   // (ie a FIELD<double>)
2157
2158   if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2159       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"));
2160
2161   END_OF(LOC);
2162 }
2163
2164 /*!
2165   Fill values array with length values.
2166 */
2167 template <class T> void FIELD<T>::getLength() const throw (MEDEXCEPTION)
2168 {
2169   const char * LOC = "FIELD<double>::getLength() const : ";
2170   BEGIN_OF(LOC);
2171
2172   // The field has to be initilised by a non empty support and a
2173   // number of components = 1 and its value type has to be set to MED_REEL64
2174   // (ie a FIELD<double>)
2175
2176   if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2177       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"));
2178
2179   END_OF(LOC);
2180 }
2181
2182 /*!
2183   Fill values array with normal values.
2184 */
2185 template <class T> void FIELD<T>::getNormal() const throw (MEDEXCEPTION)
2186 {
2187   const char * LOC = "FIELD<double>::getNormal() const : ";
2188   BEGIN_OF(LOC);
2189
2190   // The field has to be initilised by a non empty support and a
2191   // number of components = 1 and its value type has to be set to MED_REEL64
2192   // (ie a FIELD<double>)
2193
2194   if (_support == (SUPPORT *) NULL)
2195       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"));
2196
2197   int dim_space = _support->getMesh()->getSpaceDimension();
2198
2199   if ((_numberOfComponents != dim_space) || (_valueType != MED_REEL64))
2200       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"));
2201
2202   END_OF(LOC);
2203 }
2204
2205 /*!
2206   Fill values array with barycenter values.
2207 */
2208 template <class T> void FIELD<T>::getBarycenter() const throw (MEDEXCEPTION)
2209 {
2210   const char * LOC = "FIELD<double>::getBarycenter() const : ";
2211   BEGIN_OF(LOC);
2212
2213   // The field has to be initilised by a non empty support and a number of
2214   //components = space dimension and its value type has to be set to MED_REEL64
2215   // (ie a FIELD<double>)
2216
2217   if (_support == (SUPPORT *) NULL)
2218       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"));
2219
2220   int dim_space = _support->getMesh()->getSpaceDimension();
2221
2222   if ((_numberOfComponents != dim_space) || (_valueType != MED_REEL64))
2223       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"));
2224
2225   END_OF(LOC);
2226 }
2227
2228 #endif /* FIELD_HXX */