13 #include "utilities.h"
14 #include "MEDMEM_Exception.hxx"
15 #include "MEDMEM_define.hxx"
17 #include "MEDMEM_Support.hxx"
18 #include "MEDMEM_Unit.hxx"
19 #include "MEDMEM_Array.hxx"
20 #include "MEDMEM_GenDriver.hxx"
21 #include "MEDMEM_DriverFactory.hxx"
25 This class contains all the informations related with a template class FIELD :
26 - Components descriptions
27 - Time step description
28 - Location of the values (a SUPPORT class)
34 template<class T> class FIELD;
36 class FIELD_ // GENERIC POINTER TO a template <class T> class FIELD
56 Pointer to the support the field deals with.
59 const SUPPORT * _support ;
63 Number of field's components.
66 int _numberOfComponents ;
69 Number of field's values.
76 Array of size _numberOfComponents. \n
77 (constant, scalar, vector, tensor)\n
78 We could use an array of integer to store
81 - space dimension for vector,\n
82 - space dimension square for tensor.\n
83 So numbers of values per entities would be
84 sum of _componentsTypes array.
86 Not implemented yet! All type are scalar !
89 int * _componentsTypes ;
92 Array of size _numberOfComponents
93 storing components names if any.
96 string * _componentsNames;
99 Array of size _numberOfComponents
100 storing components descriptions if any.
103 string * _componentsDescriptions;
106 Array of size _numberOfComponents
107 storing components units if any.
110 UNIT * _componentsUnits;
113 Array of size _numberOfComponents
114 storing components units if any.
117 string * _MEDComponentsUnits;
118 int _iterationNumber ;
122 // _valueType should be a static const. Here is an initialization exemple
123 // template < classType T > struct SET_VALUE_TYPE { static const med_type_champ _valueType = 0; }
124 // template < > struct SET_VALUE_TYPE<double> { static const med_type_champ _valueType = MED_EN::MED_REEL64; }
125 // template < > struct SET_VALUE_TYPE<int> { static const med_type_champ _valueType = MED_EN::MED_INT32; }
126 // static const med_type_champ _valueType = SET_VALUE_TYPE <T>::_valueType;
127 MED_EN::med_type_champ _valueType ;
129 vector<GENDRIVER *> _drivers; // Storage of the drivers currently in use
130 static void _checkFieldCompatibility(const FIELD_& m, const FIELD_& n ) throw (MEDEXCEPTION);
131 void _checkNormCompatibility(const FIELD<double>* p_field_volume=NULL) const throw (MEDEXCEPTION);
132 FIELD<double>* _getFieldSize() const;
136 friend class MED_MED_RDONLY_DRIVER;
137 friend class MED_MED_WRONLY_DRIVER;
138 friend class MED_MED_RDWR_DRIVER;
140 friend class VTK_MED_DRIVER;
149 FIELD_(const SUPPORT * Support, const int NumberOfComponents);
153 FIELD_(const FIELD_ &m);
160 // virtual void setIterationNumber (int IterationNumber);
161 // virtual void setOrderNumber (int OrderNumber);
162 // virtual void setFieldName (string& fieldName);
164 virtual void rmDriver(int index);
165 virtual int addDriver(driverTypes driverType,
166 const string & fileName="Default File Name.med",
167 const string & driverFieldName="Default Field Nam",
168 MED_EN::med_mode_acces access=MED_EN::MED_REMP) ;
170 virtual int addDriver( GENDRIVER & driver);
171 virtual void read (const GENDRIVER &);
172 virtual void read(int index=0);
173 virtual void openAppend( void );
174 virtual void write(const GENDRIVER &);
175 virtual void write(int index=0, const string & driverName="");
177 virtual void writeAppend(const GENDRIVER &);
178 virtual void writeAppend(int index=0, const string & driverName="");
180 inline void setName(const string Name);
181 inline string getName() const;
182 inline void setDescription(const string Description);
183 inline string getDescription() const;
184 inline const SUPPORT * getSupport() const;
185 inline void setSupport(const SUPPORT * support);
186 inline void setNumberOfComponents(const int NumberOfComponents);
187 inline int getNumberOfComponents() const;
188 inline void setNumberOfValues(const int NumberOfValues);
189 inline int getNumberOfValues() const;
190 // inline void setComponentType(int *ComponentType);
191 // inline int * getComponentType() const;
192 // inline int getComponentTypeI(int i) const;
193 inline void setComponentsNames(const string * ComponentsNames);
194 inline void setComponentName(int i, const string ComponentName);
195 inline const string * getComponentsNames() const;
196 inline string getComponentName(int i) const;
197 inline void setComponentsDescriptions(const string * ComponentsDescriptions);
198 inline void setComponentDescription(int i, const string ComponentDescription);
199 inline const string * getComponentsDescriptions() const;
200 inline string getComponentDescription(int i) const;
202 // provisoire : en attendant de regler le probleme des unites !
203 inline void setComponentsUnits(const UNIT * ComponentsUnits);
204 inline const UNIT * getComponentsUnits() const;
205 inline const UNIT * getComponentUnit(int i) const;
206 inline void setMEDComponentsUnits(const string * MEDComponentsUnits);
207 inline void setMEDComponentUnit(int i, const string MEDComponentUnit);
208 inline const string * getMEDComponentsUnits() const;
209 inline string getMEDComponentUnit(int i) const;
211 inline void setIterationNumber(int IterationNumber);
212 inline int getIterationNumber() const;
213 inline void setTime(double Time);
214 inline double getTime() const;
215 inline void setOrderNumber(int OrderNumber);
216 inline int getOrderNumber() const;
218 inline void setValueType (const MED_EN::med_type_champ ValueType) ;
219 inline MED_EN::med_type_champ getValueType () const;
223 // ---------------------------------
224 // Implemented Methods : constructor
225 // ---------------------------------
233 inline void FIELD_::setName(const string Name)
240 inline string FIELD_::getName() const
245 Set FIELD description.
247 inline void FIELD_::setDescription(const string Description)
249 _description=Description;
252 Get FIELD description.
254 inline string FIELD_::getDescription() const
259 Set FIELD number of components.
261 inline void FIELD_::setNumberOfComponents(const int NumberOfComponents)
263 _numberOfComponents=NumberOfComponents;
266 Get FIELD number of components.
268 inline int FIELD_::getNumberOfComponents() const
270 return _numberOfComponents ;
273 Set FIELD number of values.
275 It must be the same than in the associated SUPPORT object.
277 inline void FIELD_::setNumberOfValues(const int NumberOfValues)
279 _numberOfValues=NumberOfValues;
282 Get FIELD number of value.
284 inline int FIELD_::getNumberOfValues() const
286 return _numberOfValues ;
289 // inline void FIELD_::setComponentType(int *ComponentType)
291 // _componentsTypes=ComponentType ;
293 // inline int * FIELD_::getComponentType() const
295 // return _componentsTypes ;
297 // inline int FIELD_::getComponentTypeI(int i) const
299 // return _componentsTypes[i-1] ;
303 Set FIELD components names.
305 Duplicate the ComponentsNames string array to put components names in
306 FIELD. ComponentsNames size must be equal to number of components.
308 inline void FIELD_::setComponentsNames(const string * ComponentsNames)
310 if (NULL == _componentsNames)
311 _componentsNames = new string[_numberOfComponents] ;
312 for (int i=0; i<_numberOfComponents; i++)
313 _componentsNames[i]=ComponentsNames[i] ;
316 Set FIELD i^th component name.
318 i must be >=1 and <= number of components.
320 inline void FIELD_::setComponentName(int i, const string ComponentName)
322 _componentsNames[i-1]=ComponentName ;
325 Get a reference to the string array which contain the components names.
327 This Array size is equal to number of components
329 inline const string * FIELD_::getComponentsNames() const
331 return _componentsNames ;
334 Get the name of the i^th component.
336 inline string FIELD_::getComponentName(int i) const
338 return _componentsNames[i-1] ;
341 Set FIELD components descriptions.
343 Duplicate the ComponentsDescriptions string array to put components
344 descriptions in FIELD.
345 ComponentsDescriptions size must be equal to number of components.
347 inline void FIELD_::setComponentsDescriptions(const string * ComponentsDescriptions)
349 if (NULL == _componentsDescriptions)
350 _componentsDescriptions = new string[_numberOfComponents] ;
351 for (int i=0; i<_numberOfComponents; i++)
352 _componentsDescriptions[i]=ComponentsDescriptions[i] ;
355 Set FIELD i^th component description.
357 i must be >=1 and <= number of components.
359 inline void FIELD_::setComponentDescription(int i,const string ComponentDescription)
361 _componentsDescriptions[i-1]=ComponentDescription ;
364 Get a reference to the string array which contain the components descriptions.
366 This Array size is equal to number of components
368 inline const string * FIELD_::getComponentsDescriptions() const
370 return _componentsDescriptions ;
373 Get the description of the i^th component.
375 inline string FIELD_::getComponentDescription(int i) const
377 return _componentsDescriptions[i-1];
382 Set FIELD components UNIT.
384 Duplicate the ComponentsUnits UNIT array to put components
386 ComponentsUnits size must be equal to number of components.
388 inline void FIELD_::setComponentsUnits(const UNIT * ComponentsUnits)
390 if (NULL == _componentsUnits)
391 _componentsUnits = new UNIT[_numberOfComponents] ;
392 for (int i=0; i<_numberOfComponents; i++)
393 _componentsUnits[i]=ComponentsUnits[i] ;
396 Get a reference to the UNIT array which contain the components units.
398 This Array size is equal to number of components
400 inline const UNIT * FIELD_::getComponentsUnits() const
402 return _componentsUnits ;
405 Get the UNIT of the i^th component.
407 inline const UNIT * FIELD_::getComponentUnit(int i) const
409 return &_componentsUnits[i-1] ;
412 Set FIELD components unit.
414 Duplicate the MEDComponentsUnits string array to put components
416 MEDComponentsUnits size must be equal to number of components.
419 inline void FIELD_::setMEDComponentsUnits(const string * MEDComponentsUnits)
421 if (NULL == _MEDComponentsUnits)
422 _MEDComponentsUnits = new string[_numberOfComponents] ;
423 for (int i=0; i<_numberOfComponents; i++)
424 _MEDComponentsUnits[i]=MEDComponentsUnits[i] ;
427 Set FIELD i^th component unit.
429 i must be >=1 and <= number of components.
431 inline void FIELD_::setMEDComponentUnit(int i, const string MEDComponentUnit)
433 _MEDComponentsUnits[i-1]=MEDComponentUnit ;
436 Get a reference to the string array which contain the components units.
438 This Array size is equal to number of components
440 inline const string * FIELD_::getMEDComponentsUnits() const
442 return _MEDComponentsUnits ;
445 Get the string for unit of the i^th component.
447 inline string FIELD_::getMEDComponentUnit(int i) const
449 return _MEDComponentsUnits[i-1] ;
452 Set the iteration number where FIELD has been calculated.
454 inline void FIELD_::setIterationNumber(int IterationNumber)
456 _iterationNumber=IterationNumber;
459 Get the iteration number where FIELD has been calculated.
461 inline int FIELD_::getIterationNumber() const
463 return _iterationNumber ;
466 Set the time (in second) where FIELD has been calculated.
468 inline void FIELD_::setTime(double Time)
473 Get the time (in second) where FIELD has been calculated.
475 inline double FIELD_::getTime() const
480 Set the order number where FIELD has been calculated.
482 It corresponds to internal iteration during one time step.
484 inline void FIELD_::setOrderNumber(int OrderNumber)
486 _orderNumber=OrderNumber ;
489 Get the order number where FIELD has been calculated.
491 inline int FIELD_::getOrderNumber() const
493 return _orderNumber ;
496 Get a reference to the SUPPORT object associated to FIELD.
498 inline const SUPPORT * FIELD_::getSupport() const
503 Set the reference to the SUPPORT object associated to FIELD.
505 Reference is not duplicate, so it must not be deleted.
507 inline void FIELD_::setSupport(const SUPPORT * support)
512 Get the FIELD med value type (MED_INT32 or MED_REEL64).
514 inline MED_EN::med_type_champ FIELD_::getValueType () const
519 Set the FIELD med value type (MED_INT32 or MED_REEL64).
521 inline void FIELD_::setValueType (const MED_EN::med_type_champ ValueType)
523 _valueType = ValueType ;
526 }//End namespace MEDMEM
528 /////////////////////////
529 // END OF CLASS FIELD_ //
530 /////////////////////////
534 This template class contains informations related with a FIELD :
535 - Values of the field
541 template<class T2> class MED_FIELD_RDONLY_DRIVER;
542 template<class T2> class MED_FIELD_WRONLY_DRIVER;
543 template<class T2> class VTK_FIELD_DRIVER;
545 template <class T> class FIELD : public FIELD_
548 // ------ End of Drivers Management Part
550 // array of value of type T
551 MEDARRAY<T> *_value ;
554 void _operation(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode, char* Op);
555 void _operationInitialize(const FIELD& m,const FIELD& n, char* Op);
556 void _add_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode);
557 void _sub_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode);
558 void _mul_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode);
559 void _div_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode) throw (MEDEXCEPTION);
562 FIELD & operator=(const FIELD &m); // A FAIRE
565 FIELD(const FIELD &m);
566 FIELD(const SUPPORT * Support, const int NumberOfComponents, const MED_EN::medModeSwitch Mode=MED_EN::MED_FULL_INTERLACE) throw (MEDEXCEPTION) ; // Ajout NB Constructeur FIELD avec allocation de memoire de tous ses attribut
567 FIELD(const SUPPORT * Support, driverTypes driverType,
568 const string & fileName="", const string & fieldName="",
569 const int iterationNumber = -1, const int orderNumber = -1)
570 throw (MEDEXCEPTION);
573 const FIELD operator+(const FIELD& m) const;
574 const FIELD operator-(const FIELD& m) const;
575 const FIELD operator*(const FIELD& m) const;
576 const FIELD operator/(const FIELD& m) const;
577 const FIELD operator-() const;
578 FIELD& operator+=(const FIELD& m);
579 FIELD& operator-=(const FIELD& m);
580 FIELD& operator*=(const FIELD& m);
581 FIELD& operator/=(const FIELD& m);
582 static FIELD* add(const FIELD& m, const FIELD& n);
583 static FIELD* sub(const FIELD& m, const FIELD& n);
584 static FIELD* mul(const FIELD& m, const FIELD& n);
585 static FIELD* div(const FIELD& m, const FIELD& n);
586 double normMax() const throw (MEDEXCEPTION);
587 double norm2() const throw (MEDEXCEPTION);
588 void applyLin(T a, T b);
589 template <T T_function(T)> void applyFunc();
590 static FIELD* scalarProduct(const FIELD& m, const FIELD& n);
591 double normL2(int component, const FIELD<double> * p_field_volume=NULL) const;
592 double normL2(const FIELD<double> * p_field_volume=NULL) const;
593 double normL1(int component, const FIELD<double> * p_field_volume=NULL) const;
594 double normL1(const FIELD<double> * p_field_volume=NULL) const;
596 friend class MED_FIELD_RDONLY_DRIVER<T>;
597 friend class MED_FIELD_WRONLY_DRIVER<T>;
598 friend class VTK_FIELD_DRIVER<T>;
599 //friend class MED_FIELD_RDWR_DRIVER <T>;
602 void rmDriver(int index=0);
603 int addDriver(driverTypes driverType,
604 const string & fileName="Default File Name.med",
605 const string & driverFieldName="Default Field Name",
606 MED_EN::med_mode_acces access=MED_EN::MED_REMP) ;
608 int addDriver(GENDRIVER & driver);
610 void allocValue(const int NumberOfComponents);
611 void allocValue(const int NumberOfComponents, const int LengthValue);
615 inline void read(int index=0);
616 inline void read(const GENDRIVER & genDriver);
617 inline void write(int index=0, const string & driverName = "");
618 inline void write(const GENDRIVER &);
620 inline void writeAppend(int index=0, const string & driverName = "");
621 inline void writeAppend(const GENDRIVER &);
623 inline void setValue(MEDARRAY<T> *Value);
625 inline MEDARRAY<T>* getvalue() const;
626 inline int getValueLength(MED_EN::medModeSwitch Mode) const;
627 inline const T* getValue(MED_EN::medModeSwitch Mode) const;
628 inline const T* getValueI(MED_EN::medModeSwitch Mode,int i) const;
629 inline T getValueIJ(int i,int j) const;
631 inline void setValue(MED_EN::medModeSwitch mode, T* value);
632 inline void setValueI(MED_EN::medModeSwitch mode, int i, T* value);
633 inline void setValueIJ(int i, int j, T value);
636 This fonction feeds the FIELD<double> private attributs _value with the
637 volume of each cells belonging to the argument Support. The field has to be
638 initialised via the constructor FIELD<double>(const SUPPORT * , const int )
639 with Support as SUPPORT argument, 1 has the number of components, and Support
640 has be a SUPPORT on 3D cells. This initialisation could be done by the empty
641 constructor followed by a setSupport and setNumberOfComponents call but it has
642 be followed by a setValueType(MED_REEL64) call.
644 void getVolume() const throw (MEDEXCEPTION) ;
646 This fonction feeds the FIELD<double> private attributs _value with the
647 area of each cells (or faces) belonging to the attribut _support. The field
648 has to be initialised via the constructor FIELD<double>(const SUPPORT * ,
649 const int ) with 1 has the number of components, and _support has be a
650 SUPPORT on 2D cells or 3D faces. This initialisation could be done by the
651 empty constructor followed by a setSupport and setNumberOfComponents call but
652 it has be followed by a setValueType(MED_REEL64) call.
654 void getArea() const throw (MEDEXCEPTION) ;
656 This fonction feeds the FIELD<double> private attributs _value with the
657 length of each segments belonging to the attribut _support. The field has
658 to be initialised via the constructor FIELD<double>(const SUPPORT * ,
659 const int ) with 1 has the number of components, and _support has be a
660 SUPPORT on 3D edges or 2D faces. This initialisation could be done by the
661 empty constructor followed by a setSupport and setNumberOfComponents call but
662 it has be followed by a setValueType(MED_REEL64) call.
664 void getLength() const throw (MEDEXCEPTION) ;
666 This fonction feeds the FIELD<double> private attributs _value with the
667 normal vector of each faces belonging to the attribut _support. The field
668 has to be initialised via the constructor FIELD<double>(const SUPPORT * ,
669 const int ) with the space dimension has the number of components, and
670 _support has be a SUPPORT on 3D or 2D faces. This initialisation could be done
671 by the empty constructor followed by a setSupport and setNumberOfComponents
672 call but it has be followed by a setValueType(MED_REEL64) call.
674 void getNormal() const throw (MEDEXCEPTION) ;
676 This fonction feeds the FIELD<double> private attributs _value with the
677 barycenter of each faces or cells or edges belonging to the attribut _support.
678 The field has to be initialised via the constructor
679 FIELD<double>(const SUPPORT * ,const int ) with the space dimension has the
680 number of components, and _support has be a SUPPORT on 3D cells or 2D faces.
681 This initialisation could be done by the empty constructor followed by a
682 setSupport and setNumberOfComponents call but it has be followed by a
683 setValueType(MED_REEL64) call.
685 void getBarycenter() const throw (MEDEXCEPTION) ;
688 // --------------------
689 // Implemented Methods
690 // --------------------
693 Constructor with no parameter, most of the attribut members are set to NULL.
695 template <class T> FIELD<T>::FIELD():
696 _value((MEDARRAY<T>*)NULL)
698 MESSAGE("Constructeur FIELD sans parametre");
702 Constructor with parameters such that all attrribut are set but the _value
703 attribut is allocated but not set.
705 template <class T> FIELD<T>::FIELD(const SUPPORT * Support,
706 const int NumberOfComponents, const MED_EN::medModeSwitch Mode) throw (MEDEXCEPTION) :
707 FIELD_(Support, NumberOfComponents)
709 BEGIN_OF("FIELD<T>::FIELD(const SUPPORT * Support, const int NumberOfComponents, const medModeSwitch Mode)");
713 _numberOfValues = Support->getNumberOfElements(MED_ALL_ELEMENTS);
715 catch (MEDEXCEPTION &ex) {
716 MESSAGE("No value defined ! ("<<ex.what()<<")");
717 _value = (MEDARRAY<T>*)NULL ;
719 MESSAGE("FIELD : constructeur : "<< _numberOfValues <<" et "<< NumberOfComponents);
720 if (0<_numberOfValues) {
721 _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues,Mode);
724 _value = (MEDARRAY<T>*)NULL ;
726 END_OF("FIELD<T>::FIELD(const SUPPORT * Support, const int NumberOfComponents, const medModeSwitch Mode)");
733 template <class T> void FIELD<T>::init ()
740 template <class T> FIELD<T>::FIELD(const FIELD & m):
743 MESSAGE("Constructeur FIELD de recopie");
744 if (m._value != NULL)
746 // copie only default !
747 _value = new MEDARRAY<T>(* m._value,false);
750 _value = (MEDARRAY<T> *) NULL;
751 //_drivers = m._drivers;
759 template <class T> FIELD<T> & FIELD<T>::operator=(const FIELD &m)
761 MESSAGE("Appel de FIELD<T>::operator=") ;
762 // operator= on FIELD_
768 Overload addition operator.
769 This operation is authorized only for compatible fields that have the same support.
770 The compatibility checking includes equality tests of the folowing data members:\n
772 - _numberOfComponents
775 - _MEDComponentsUnits.
777 The data members of the returned field are initialized, based on the first field, except for the name,
778 which is the combination of the two field's names and the operator.
779 Advised Utilisation in C++ : <tt> FIELD<T> c = a + b; </tt> \n
780 In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
781 When using python, this operator calls the copy constructor in any case.
782 The user has to be aware that when using operator + in associatives expressions like
783 <tt> a = b + c + d +e; </tt> \n
784 no optimisation is performed : the evaluation of last expression requires the construction of
788 const FIELD<T> FIELD<T>::operator+(const FIELD & m) const
790 BEGIN_OF("FIELD<T>::operator+(const FIELD & m)");
791 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
793 // Select mode : avoid if possible any calculation of other mode for fields m or *this
794 MED_EN::medModeSwitch mode;
795 if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
796 mode=m.getvalue()->getMode();
798 mode=this->getvalue()->getMode();
800 // Creation of the result - memory is allocated by FIELD constructor
801 FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
802 //result._operation(*this,m,mode,"+"); // perform Atribute's initialization & addition
803 result._operationInitialize(*this,m,"+"); // perform Atribute's initialization
804 result._add_in_place(*this,m,mode); // perform addition
806 END_OF("FIELD<T>::operator+(const FIELD & m)");
810 /*! Overloaded Operator +=
811 * Operations are directly performed in the first field's array.
812 * This operation is authorized only for compatible fields that have the same support.
815 FIELD<T>& FIELD<T>::operator+=(const FIELD & m)
817 BEGIN_OF("FIELD<T>::operator+=(const FIELD & m)");
818 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
820 // We choose to keep *this mode, even if it may cost a re-calculation for m
821 MED_EN::medModeSwitch mode=this->getvalue()->getMode();
822 const T* value1=m.getValue(mode); // get pointers to the values we are adding
824 // get a non const pointer to the inside array of values and perform operation
825 T * value=const_cast<T *> (getValue(mode));
826 const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
827 const T* endV=value+size; // pointer to the end of value
828 for(;value!=endV; value1++,value++)
830 // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
831 this->getvalue()->clearOtherMode();
832 END_OF("FIELD<T>::operator+=(const FIELD & m)");
837 /*! Addition of fields. Static member function.
838 * The function return a pointer to a new created field that holds the addition.
839 * Data members are checked for compatibility and initialized.
840 * The user is in charge of memory deallocation.
843 FIELD<T>* FIELD<T>::add(const FIELD& m, const FIELD& n)
845 BEGIN_OF("FIELD<T>::add(const FIELD & m, const FIELD& n)");
846 FIELD_::_checkFieldCompatibility(m, n); // may throw exception
848 // Select mode : avoid if possible any calculation of other mode for fields m or *this
849 MED_EN::medModeSwitch mode;
850 if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
851 mode=m.getvalue()->getMode();
853 mode=n.getvalue()->getMode();
855 // Creation of a new field
856 FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
857 result->_operationInitialize(m,n,"+"); // perform Atribute's initialization
858 result->_add_in_place(m,n,mode); // perform addition
860 END_OF("FIELD<T>::add(const FIELD & m, const FIELD& n)");
865 Overload substraction operator.
866 This operation is authorized only for compatible fields that have the same support.
867 The compatibility checking includes equality tests of the folowing data members:\n
869 - _numberOfComponents
872 - _MEDComponentsUnits.
874 The data members of the returned field are initialized, based on the first field, except for the name,
875 which is the combination of the two field's names and the operator.
876 Advised Utilisation in C++ : <tt> FIELD<T> c = a - b; </tt> \n
877 In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
878 When using python, this operator calls the copy constructor in any case.
879 The user has to be aware that when using operator - in associatives expressions like
880 <tt> a = b - c - d -e; </tt> \n
881 no optimisation is performed : the evaluation of last expression requires the construction of
885 const FIELD<T> FIELD<T>::operator-(const FIELD & m) const
887 BEGIN_OF("FIELD<T>::operator-(const FIELD & m)");
888 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
890 MED_EN::medModeSwitch mode;
891 if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
892 mode=m.getvalue()->getMode();
894 mode=this->getvalue()->getMode();
896 // Creation of the result - memory is allocated by FIELD constructor
897 FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
898 //result._operation(*this,m,mode,"-"); // perform Atribute's initialization & substraction
899 result._operationInitialize(*this,m,"-"); // perform Atribute's initialization
900 result._sub_in_place(*this,m,mode); // perform substracion
902 END_OF("FIELD<T>::operator-(const FIELD & m)");
907 const FIELD<T> FIELD<T>::operator-() const
909 BEGIN_OF("FIELD<T>::operator-()");
911 MED_EN::medModeSwitch mode=this->getvalue()->getMode();
912 // Creation of the result - memory is allocated by FIELD constructor
913 FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
914 // Atribute's initialization
915 result.setName("- "+getName());
916 result.setComponentsNames(getComponentsNames());
917 // not yet implemented setComponentType(getComponentType());
918 result.setComponentsDescriptions(getComponentsDescriptions());
919 result.setMEDComponentsUnits(getMEDComponentsUnits());
920 result.setComponentsUnits(getComponentsUnits());
921 result.setIterationNumber(getIterationNumber());
922 result.setTime(getTime());
923 result.setOrderNumber(getOrderNumber());
924 result.setValueType(getValueType());
926 const T* value1=getValue(mode);
927 // get a non const pointer to the inside array of values and perform operation
928 T * value=const_cast<T *> (result.getValue(mode));
929 const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
930 const T* endV=value+size; // pointer to the end of value
932 for(;value!=endV; value1++,value++)
934 END_OF("FIELD<T>::operator-=(const FIELD & m)");
938 /*! Overloaded Operator -=
939 * Operations are directly performed in the first field's array.
940 * This operation is authorized only for compatible fields that have the same support.
943 FIELD<T>& FIELD<T>::operator-=(const FIELD & m)
945 BEGIN_OF("FIELD<T>::operator-=(const FIELD & m)");
946 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
948 // We choose to keep *this mode, even if it may cost a re-calculation for m
949 MED_EN::medModeSwitch mode=this->getvalue()->getMode();
950 const T* value1=m.getValue(mode); // get pointers to the values we are adding
952 // get a non const pointer to the inside array of values and perform operation
953 T * value=const_cast<T *> (getValue(mode));
954 const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
955 const T* endV=value+size; // pointer to the end of value
957 for(;value!=endV; value1++,value++)
959 // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
960 this->getvalue()->clearOtherMode();
961 END_OF("FIELD<T>::operator-=(const FIELD & m)");
966 /*! Substraction of fields. Static member function.
967 * The function return a pointer to a new created field that holds the substraction.
968 * Data members are checked for compatibility and initialized.
969 * The user is in charge of memory deallocation.
972 FIELD<T>* FIELD<T>::sub(const FIELD& m, const FIELD& n)
974 BEGIN_OF("FIELD<T>::sub(const FIELD & m, const FIELD& n)");
975 FIELD_::_checkFieldCompatibility(m, n); // may throw exception
977 // Select mode : avoid if possible any calculation of other mode for fields m or *this
978 MED_EN::medModeSwitch mode;
979 if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
980 mode=m.getvalue()->getMode();
982 mode=n.getvalue()->getMode();
984 // Creation of a new field
985 FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
986 result->_operationInitialize(m,n,"-"); // perform Atribute's initialization
987 result->_sub_in_place(m,n,mode); // perform substraction
989 END_OF("FIELD<T>::sub(const FIELD & m, const FIELD& n)");
994 Overload multiplication operator.
995 This operation is authorized only for compatible fields that have the same support.
996 The compatibility checking includes equality tests of the folowing data members:\n
998 - _numberOfComponents
1001 - _MEDComponentsUnits.
1003 The data members of the returned field are initialized, based on the first field, except for the name,
1004 which is the combination of the two field's names and the operator.
1005 Advised Utilisation in C++ : <tt> FIELD<T> c = a * b; </tt> \n
1006 In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
1007 When using python, this operator calls the copy constructor in any case.
1008 The user has to be aware that when using operator * in associatives expressions like
1009 <tt> a = b * c * d *e; </tt> \n
1010 no optimisation is performed : the evaluation of last expression requires the construction of
1014 const FIELD<T> FIELD<T>::operator*(const FIELD & m) const
1016 BEGIN_OF("FIELD<T>::operator*(const FIELD & m)");
1017 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1019 // Select mode : avoid if possible any calculation of other mode for fields m or *this
1020 MED_EN::medModeSwitch mode;
1021 if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
1022 mode=m.getvalue()->getMode();
1024 mode=this->getvalue()->getMode();
1026 // Creation of the result - memory is allocated by FIELD constructor
1027 FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
1028 //result._operation(*this,m,mode,"*"); // perform Atribute's initialization & multiplication
1029 result._operationInitialize(*this,m,"*"); // perform Atribute's initialization
1030 result._mul_in_place(*this,m,mode); // perform multiplication
1032 END_OF("FIELD<T>::operator*(const FIELD & m)");
1036 /*! Overloaded Operator *=
1037 * Operations are directly performed in the first field's array.
1038 * This operation is authorized only for compatible fields that have the same support.
1041 FIELD<T>& FIELD<T>::operator*=(const FIELD & m)
1043 BEGIN_OF("FIELD<T>::operator*=(const FIELD & m)");
1044 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1046 // We choose to keep *this mode, even if it may cost a re-calculation for m
1047 MED_EN::medModeSwitch mode=this->getvalue()->getMode();
1048 const T* value1=m.getValue(mode); // get pointers to the values we are adding
1050 // get a non const pointer to the inside array of values and perform operation
1051 T * value=const_cast<T *> (getValue(mode));
1052 const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1053 const T* endV=value+size; // pointer to the end of value
1055 for(;value!=endV; value1++,value++)
1057 // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1058 this->getvalue()->clearOtherMode();
1059 END_OF("FIELD<T>::operator*=(const FIELD & m)");
1064 /*! Multiplication of fields. Static member function.
1065 * The function return a pointer to a new created field that holds the multiplication.
1066 * Data members are checked for compatibility and initialized.
1067 * The user is in charge of memory deallocation.
1070 FIELD<T>* FIELD<T>::mul(const FIELD& m, const FIELD& n)
1072 BEGIN_OF("FIELD<T>::mul(const FIELD & m, const FIELD& n)");
1073 FIELD_::_checkFieldCompatibility(m, n); // may throw exception
1075 // Select mode : avoid if possible any calculation of other mode for fields m or *this
1076 MED_EN::medModeSwitch mode;
1077 if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
1078 mode=m.getvalue()->getMode();
1080 mode=n.getvalue()->getMode();
1082 // Creation of a new field
1083 FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
1084 result->_operationInitialize(m,n,"*"); // perform Atribute's initialization
1085 result->_mul_in_place(m,n,mode); // perform multiplication
1087 END_OF("FIELD<T>::mul(const FIELD & m, const FIELD& n)");
1093 Overload division operator.
1094 This operation is authorized only for compatible fields that have the same support.
1095 The compatibility checking includes equality tests of the folowing data members:\n
1097 - _numberOfComponents
1100 - _MEDComponentsUnits.
1102 The data members of the returned field are initialized, based on the first field, except for the name,
1103 which is the combination of the two field's names and the operator.
1104 Advised Utilisation in C++ : <tt> FIELD<T> c = a / b; </tt> \n
1105 In this case, the (recent) compilators perform optimisation and don't call the copy constructor.
1106 When using python, this operator calls the copy constructor in any case.
1107 The user has to be aware that when using operator / in associatives expressions like
1108 <tt> a = b / c / d /e; </tt> \n
1109 no optimisation is performed : the evaluation of last expression requires the construction of
1113 const FIELD<T> FIELD<T>::operator/(const FIELD & m) const
1115 BEGIN_OF("FIELD<T>::operator/(const FIELD & m)");
1116 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1118 // Select mode : avoid if possible any calculation of other mode for fields m or *this
1119 MED_EN::medModeSwitch mode;
1120 if(this->getvalue()->getMode()==m.getvalue()->getMode() || this->getvalue()->isOtherCalculated())
1121 mode=m.getvalue()->getMode();
1123 mode=this->getvalue()->getMode();
1125 // Creation of the result - memory is allocated by FIELD constructor
1126 FIELD<T> result(this->getSupport(),this->getNumberOfComponents(),mode);
1127 //result._operation(*this,m,mode,"/"); // perform Atribute's initialization & division
1128 result._operationInitialize(*this,m,"/"); // perform Atribute's initialization
1129 result._div_in_place(*this,m,mode); // perform division
1131 END_OF("FIELD<T>::operator/(const FIELD & m)");
1136 /*! Overloaded Operator /=
1137 * Operations are directly performed in the first field's array.
1138 * This operation is authorized only for compatible fields that have the same support.
1141 FIELD<T>& FIELD<T>::operator/=(const FIELD & m)
1143 BEGIN_OF("FIELD<T>::operator/=(const FIELD & m)");
1144 FIELD_::_checkFieldCompatibility(*this, m); // may throw exception
1146 // We choose to keep *this mode, even if it may cost a re-calculation for m
1147 MED_EN::medModeSwitch mode=this->getvalue()->getMode();
1148 const T* value1=m.getValue(mode); // get pointers to the values we are adding
1150 // get a non const pointer to the inside array of values and perform operation
1151 T * value=const_cast<T *> (getValue(mode));
1152 const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1153 const T* endV=value+size; // pointer to the end of value
1155 for(;value!=endV; value1++,value++)
1157 // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1158 this->getvalue()->clearOtherMode();
1159 END_OF("FIELD<T>::operator/=(const FIELD & m)");
1164 /*! Division of fields. Static member function.
1165 * The function return a pointer to a new created field that holds the division.
1166 * Data members are checked for compatibility and initialized.
1167 * The user is in charge of memory deallocation.
1170 FIELD<T>* FIELD<T>::div(const FIELD& m, const FIELD& n)
1172 BEGIN_OF("FIELD<T>::div(const FIELD & m, const FIELD& n)");
1173 FIELD_::_checkFieldCompatibility(m, n); // may throw exception
1175 // Select mode : avoid if possible any calculation of other mode for fields m or *this
1176 MED_EN::medModeSwitch mode;
1177 if(m.getvalue()->getMode()==n.getvalue()->getMode() || n.getvalue()->isOtherCalculated())
1178 mode=m.getvalue()->getMode();
1180 mode=n.getvalue()->getMode();
1182 // Creation of a new field
1183 FIELD<T>* result = new FIELD<T>(m.getSupport(),m.getNumberOfComponents(),mode);
1184 result->_operationInitialize(m,n,"/"); // perform Atribute's initialization
1185 result->_div_in_place(m,n,mode); // perform division
1187 END_OF("FIELD<T>::div(const FIELD & m, const FIELD& n)");
1194 This internal method initialize the members of a new field created to hold the result of the operation Op .
1195 Initialization is based on the first field, except for the name, which is the combination of the two field's names
1200 void FIELD<T>::_operationInitialize(const FIELD& m,const FIELD& n, char* Op)
1202 MESSAGE("Appel methode interne _add" << Op);
1204 // Atribute's initialization (copy of the first field's attributes)
1205 // Other data members (_support, _numberOfValues) are initialized in the field's constr.
1206 setName(m.getName()+" "+Op+" "+n.getName());
1207 setComponentsNames(m.getComponentsNames());
1208 // not yet implemented setComponentType(m.getComponentType());
1209 setComponentsDescriptions(m.getComponentsDescriptions());
1210 setMEDComponentsUnits(m.getMEDComponentsUnits());
1212 // The following data member may differ from field m to n.
1213 // The initialization is done based on the first field.
1214 setComponentsUnits(m.getComponentsUnits());
1215 setIterationNumber(m.getIterationNumber());
1216 setTime(m.getTime());
1217 setOrderNumber(m.getOrderNumber());
1218 setValueType(m.getValueType());
1224 Internal method called by FIELD<T>::operator+ and FIELD<T>::add to perform addition "in place".
1225 This method is applied to a just created field with medModeSwitch mode.
1226 For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1231 void FIELD<T>::_add_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode)
1233 // get pointers to the values we are adding
1234 const T* value1=m.getValue(mode);
1235 const T* value2=n.getValue(mode);
1236 // get a non const pointer to the inside array of values and perform operation
1237 T * value=const_cast<T *> (getValue(mode));
1239 const int size=getNumberOfValues()*getNumberOfComponents();
1241 const T* endV1=value1+size;
1242 for(;value1!=endV1; value1++,value2++,value++)
1243 *value=(*value1)+(*value2);
1248 Internal method called by FIELD<T>::operator- and FIELD<T>::sub to perform substraction "in place".
1249 This method is applied to a just created field with medModeSwitch mode.
1250 For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1255 void FIELD<T>::_sub_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode)
1257 // get pointers to the values we are adding
1258 const T* value1=m.getValue(mode);
1259 const T* value2=n.getValue(mode);
1260 // get a non const pointer to the inside array of values and perform operation
1261 T * value=const_cast<T *> (getValue(mode));
1263 const int size=getNumberOfValues()*getNumberOfComponents();
1265 const T* endV1=value1+size;
1266 for(;value1!=endV1; value1++,value2++,value++)
1267 *value=(*value1)-(*value2);
1272 Internal method called by FIELD<T>::operator* and FIELD<T>::mul to perform multiplication "in place".
1273 This method is applied to a just created field with medModeSwitch mode.
1274 For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1279 void FIELD<T>::_mul_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode)
1281 // get pointers to the values we are adding
1282 const T* value1=m.getValue(mode);
1283 const T* value2=n.getValue(mode);
1284 // get a non const pointer to the inside array of values and perform operation
1285 T * value=const_cast<T *> (getValue(mode));
1287 const int size=getNumberOfValues()*getNumberOfComponents();
1289 const T* endV1=value1+size;
1290 for(;value1!=endV1; value1++,value2++,value++)
1291 *value=(*value1)*(*value2);
1296 Internal method called by FIELD<T>::operator/ and FIELD<T>::div to perform division "in place".
1297 This method is applied to a just created field with medModeSwitch mode.
1298 For this reason, the alternate mode doesn't need to be set to 0 after performing operation :
1303 void FIELD<T>::_div_in_place(const FIELD& m,const FIELD& n, const MED_EN::medModeSwitch mode) throw (MEDEXCEPTION)
1305 // get pointers to the values we are adding
1306 const T* value1=m.getValue(mode);
1307 const T* value2=n.getValue(mode);
1308 // get a non const pointer to the inside array of values and perform operation
1309 T * value=const_cast<T *> (getValue(mode));
1311 const int size=getNumberOfValues()*getNumberOfComponents();
1313 const T* endV1=value1+size;
1314 for(;value1!=endV1; value1++,value2++,value++){
1315 if ( *value2 == 0 ) {
1317 diagnosis="FIELD<T>::_div_in_place(...) : Divide by zero !";
1318 throw MEDEXCEPTION(diagnosis.c_str());
1320 *value=(*value1)/(*value2);
1326 template <class T> double FIELD<T>::normMax() const throw (MEDEXCEPTION)
1328 const T* value=getValue(getvalue()->getMode()); // get pointer to the values
1329 const int size=getNumberOfValues()*getNumberOfComponents();
1330 if (size <= 0) // Size of array has to be strictly positive
1333 diagnosis="FIELD<T>::normMax() : cannot compute the norm of "+getName()+
1334 " : it size is non positive!";
1335 throw MEDEXCEPTION(diagnosis.c_str());
1337 const T* lastvalue=value+size; // get pointer just after last value
1338 const T* pMax=value; // pointer to the max value
1339 const T* pMin=value; // pointer to the min value
1341 // get pointers to the max & min value of array
1342 while ( ++value != lastvalue )
1344 if ( *pMin > *value )
1346 if ( *pMax < *value )
1350 T Max= *pMax>(T) 0 ? *pMax : -*pMax; // Max=abs(*pMax)
1351 T Min= *pMin>(T) 0 ? *pMin : -*pMin; // Min=abs(*pMin)
1353 return Max>Min ? static_cast<double>(Max) : static_cast<double>(Min);
1356 /*! Return Euclidien norm
1358 template <class T> double FIELD<T>::norm2() const throw (MEDEXCEPTION)
1360 const T* value=this->getValue(this->getvalue()->getMode()); // get const pointer to the values
1361 const int size=getNumberOfValues()*getNumberOfComponents(); // get size of array
1362 if (size <= 0) // Size of array has to be strictly positive
1365 diagnosis="FIELD<T>::norm2() : cannot compute the norm of "+getName()+
1366 " : it size is non positive!";
1367 throw MEDEXCEPTION(diagnosis.c_str());
1369 const T* lastvalue=value+size; // point just after last value
1371 T result((T)0); // init
1372 for( ; value!=lastvalue ; ++value)
1373 result += (*value) * (*value);
1375 return std::sqrt(static_cast<double> (result));
1379 /*! Apply to each (scalar) field component the template parameter T_function,
1380 * which is a pointer to function.
1381 * Since the pointer is known at compile time, the function is inlined into the inner loop!
1382 * calculation is done "in place".
1385 * \code myField.applyFunc<std::sqrt>(); // apply sqare root function \endcode
1386 * \code myField.applyFunc<myFunction>(); // apply your own created function \endcode
1388 template <class T> template <T T_function(T)>
1389 void FIELD<T>::applyFunc()
1391 MED_EN::medModeSwitch mode=getvalue()->getMode();
1393 // get a non const pointer to the inside array of values and perform operation
1394 T * value=const_cast<T *> (getValue(mode));
1395 const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1397 if (size>0) // for a negative size, there is nothing to do
1399 const T* lastvalue=value+size; // pointer to the end of value
1400 for(;value!=lastvalue; ++value) // apply linear transformation
1401 *value = T_function(*value);
1402 // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1403 getvalue()->clearOtherMode();
1408 /*! Apply to each (scalar) field component the linear function x -> ax+b.
1409 * calculation is done "in place".
1411 template <class T> void FIELD<T>::applyLin(T a, T b)
1413 MED_EN::medModeSwitch mode=getvalue()->getMode();
1415 // get a non const pointer to the inside array of values and perform operation in place
1416 T * value=const_cast<T *> (getValue(mode));
1417 const int size=getNumberOfValues()*getNumberOfComponents(); // size of array
1419 if (size>0) // for a negative size, there is nothing to do
1421 const T* lastvalue=value+size; // pointer to the end of value
1422 for(;value!=lastvalue; ++value) // apply linear transformation
1423 *value = a*(*value)+b;
1424 // if it exists, the alternate mode needs to be cleared because it is not up-to-date anymore
1425 getvalue()->clearOtherMode();
1431 * Return a pointer to a new field that holds the scalar product. Static member function.
1432 * This operation is authorized only for compatible fields that have the same support.
1433 * The compatibility checking includes equality tests of the folowing data members:\n
1435 * - _numberOfComponents
1437 * - _componentsTypes
1438 * - _MEDComponentsUnits.
1439 * Data members are initialized.
1440 * The new field point to the same support and has one component.
1441 * Each value of it is the scalar product of the two argument's fields.
1442 * The user is in charge of memory deallocation.
1444 template <class T> FIELD<T>* FIELD<T>::scalarProduct(const FIELD & m, const FIELD & n)
1446 FIELD_::_checkFieldCompatibility( m, n); // may throw exception
1447 // we need a MED_FULL_INTERLACE representation of m & n to compute the scalar product
1448 const MED_EN::medModeSwitch mode=MED_EN::MED_FULL_INTERLACE;
1450 const int numberOfElements=m.getNumberOfValues(); // strictly positive
1451 const int NumberOfComponents=m.getNumberOfComponents(); // strictly positive
1453 // Creation & init of a the result field on the same suppot, with one component
1454 FIELD<T>* result = new FIELD<T>(m.getSupport(),1,mode);
1455 result->setName( "scalarProduct ( " + m.getName() + " , " + n.getName() + " )" );
1456 result->setIterationNumber(m.getIterationNumber());
1457 result->setTime(m.getTime());
1458 result->setOrderNumber(m.getOrderNumber());
1459 result->setValueType(m.getValueType());
1461 const T* value1=m.getValue(mode); // get const pointer to the values
1462 const T* value2=n.getValue(mode); // get const pointer to the values
1463 // get a non const pointer to the inside array of values and perform operation
1464 T * value=const_cast<T *> (result->getValue(mode));
1466 const T* lastvalue=value+numberOfElements; // pointing just after last value of result
1467 for ( ; value!=lastvalue ; ++value ) // loop on all elements
1469 *value=(T)0; // initialize value
1470 const T* endofRow=value1+NumberOfComponents; // pointing just after end of row
1471 for ( ; value1 != endofRow; ++value1, ++value2) // computation of dot product
1472 *value += (*value1) * (*value2);
1477 /*! Return L2 Norm of the field's component.
1478 * Cannot be applied to a field with a support on nodes.
1479 * If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1481 template <class T> double FIELD<T>::normL2(int component, const FIELD<double> * p_field_volume) const
1483 _checkNormCompatibility(p_field_volume); // may throw exception
1484 if ( component<1 || component>getNumberOfComponents() )
1485 throw MEDEXCEPTION(STRING("FIELD<T>::normL2() : The component argument should be between 1 and the number of components"));
1487 const FIELD<double> * p_field_size=p_field_volume;
1488 if(!p_field_volume) // if the user don't supply the volume
1489 p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1491 // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1492 const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1493 const T* value=getValueI( MED_NO_INTERLACE, component); // get pointer to the component's values
1494 const T* lastvalue=value+getNumberOfValues(); // pointing just after the end of column
1496 double integrale=0.0;
1498 for (; value!=lastvalue ; ++value ,++vol)
1500 integrale += static_cast<double>((*value) * (*value)) * (*vol);
1504 if(!p_field_volume) // if the user didn't supply the volume
1505 delete p_field_size; // delete temporary volume field
1507 throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1509 return integrale/totVol;
1512 /*! Return L2 Norm of the field.
1513 * Cannot be applied to a field with a support on nodes.
1514 * If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1516 template <class T> double FIELD<T>::normL2(const FIELD<double> * p_field_volume) const
1518 _checkNormCompatibility(p_field_volume); // may throw exception
1519 const FIELD<double> * p_field_size=p_field_volume;
1520 if(!p_field_volume) // if the user don't supply the volume
1521 p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1523 // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1524 const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1525 const double* lastvol=vol+getNumberOfValues(); // pointing just after the end of vol
1526 const T* value=getValue( MED_NO_INTERLACE); // get pointer to the field's values
1529 const double* p_vol=vol;
1530 for (p_vol=vol; p_vol!=lastvol ; ++p_vol) // calculate total volume
1533 double integrale=0.0;
1534 for (int i=1; i<=getNumberOfComponents(); ++i) // compute integral on all components
1535 for (p_vol=vol; p_vol!=lastvol ; ++value ,++p_vol)
1536 integrale += static_cast<double>((*value) * (*value)) * (*p_vol);
1538 if(!p_field_volume) // if the user didn't supply the volume
1539 delete p_field_size; // delete temporary volume field
1541 throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1543 return integrale/totVol;
1546 /*! Return L1 Norm of the field's component.
1547 * Cannot be applied to a field with a support on nodes.
1548 * If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1550 template <class T> double FIELD<T>::normL1(int component, const FIELD<double> * p_field_volume) const
1552 _checkNormCompatibility(p_field_volume); // may throw exception
1553 if ( component<1 || component>getNumberOfComponents() )
1554 throw MEDEXCEPTION(STRING("FIELD<T>::normL2() : The component argument should be between 1 and the number of components"));
1556 const FIELD<double> * p_field_size=p_field_volume;
1557 if(!p_field_volume) // if the user don't supply the volume
1558 p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1560 // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1561 const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1562 const T* value=getValueI( MED_NO_INTERLACE, component); // get pointer to the component's values
1563 const T* lastvalue=value+getNumberOfValues(); // pointing just after the end of column
1565 double integrale=0.0;
1567 for (; value!=lastvalue ; ++value ,++vol)
1569 integrale += std::abs( static_cast<double>(*value) ) * (*vol);
1573 if(!p_field_volume) // if the user didn't supply the volume
1574 delete p_field_size; // delete temporary volume field
1576 throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1578 return integrale/totVol;
1581 /*! Return L1 Norm of the field.
1582 * Cannot be applied to a field with a support on nodes.
1583 * If the optional p_field_volume argument is furnished, the volume is not re-calculated.
1585 template <class T> double FIELD<T>::normL1(const FIELD<double> * p_field_volume) const
1587 _checkNormCompatibility(p_field_volume); // may throw exception
1588 const FIELD<double> * p_field_size=p_field_volume;
1589 if(!p_field_volume) // if the user don't supply the volume
1590 p_field_size=_getFieldSize(); // we calculate the volume [PROVISOIRE, en attendant l'implémentation dans mesh]
1592 // get pointer to the element's volumes. MED_FULL_INTERLACE is the default mode for p_field_size
1593 const double* vol=p_field_size->getValue(MED_FULL_INTERLACE);
1594 const double* lastvol=vol+getNumberOfValues(); // pointing just after the end of vol
1595 const T* value=getValue( MED_NO_INTERLACE); // get pointer to the field's values
1598 const double* p_vol=vol;
1599 for (p_vol=vol; p_vol!=lastvol ; ++p_vol) // calculate total volume
1602 double integrale=0.0;
1603 for (int i=1; i<=getNumberOfComponents(); ++i) // compute integral on all components
1604 for (p_vol=vol; p_vol!=lastvol ; ++value ,++p_vol)
1605 integrale += std::abs( static_cast<double>(*value) ) * (*p_vol);
1607 if(!p_field_volume) // if the user didn't supply the volume
1608 delete p_field_size; // delete temporary volume field
1610 throw MEDEXCEPTION(STRING("cannot compute sobolev norm : volume is not positive!"));
1612 return integrale/totVol;
1619 Constructor with parameters; the object is set via a file and its associated
1620 driver. For the moment only the MED_DRIVER is considered and if the last two
1621 argument (iterationNumber and orderNumber) are not set; their default value
1622 is -1. If the field fieldDriverName with the iteration number
1623 iterationNumber and the order number orderNumber does not exist in the file
1624 fieldDriverName; the constructor raises an exception.
1626 template <class T> FIELD<T>::FIELD(const SUPPORT * Support,
1627 driverTypes driverType,
1628 const string & fileName/*=""*/,
1629 const string & fieldDriverName/*=""*/,
1630 const int iterationNumber,
1631 const int orderNumber)
1632 throw (MEDEXCEPTION)
1634 const char * LOC = "template <class T> FIELD<T>::FIELD(const SUPPORT * Support, driverTypes driverType, const string & fileName=\"\", const string & fieldName=\"\", const int iterationNumber=-1, const int orderNumber=-1) : ";
1643 _value = (MEDARRAY<T>*)NULL;
1645 _iterationNumber = iterationNumber;
1647 _orderNumber = orderNumber;
1649 current = addDriver(driverType,fileName,fieldDriverName,MED_LECT);
1651 // switch(driverType)
1653 // case MED_DRIVER :
1655 // MED_FIELD_RDONLY_DRIVER<T> myDriver(fileName,this);
1656 // myDriver.setFieldName(fieldDriverName);
1657 // current = addDriver(myDriver);
1660 // current = addDriver(driverType,fileName,fieldDriverName);
1661 // switch(_drivers[current]->getAccessMode() ) {
1662 // case MED_WRONLY : {
1663 // MESSAGE("FIELD<T>::FIELD(driverTypes driverType, .....) : driverType must have a MED_RDONLY or MED_RDWR accessMode");
1664 // rmDriver(current);
1671 // throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Driver type unknown : can't create driver!"));
1676 _drivers[current]->open();
1677 _drivers[current]->read();
1678 _drivers[current]->close();
1686 template <class T> FIELD<T>::~FIELD()
1688 BEGIN_OF(" Destructeur FIELD<T>::~FIELD()");
1690 if (_value) delete _value;
1691 END_OF(" Destructeur FIELD<T>::~FIELD()");
1697 template <class T> void FIELD<T>::allocValue(const int NumberOfComponents)
1699 const char* LOC = "FIELD<T>::allocValue(const int NumberOfComponents)" ;
1702 _numberOfComponents = NumberOfComponents ;
1703 if (_componentsTypes == NULL)
1704 _componentsTypes = new int[NumberOfComponents] ;
1705 if (_componentsNames == NULL)
1706 _componentsNames = new string[NumberOfComponents];
1707 if (_componentsDescriptions == NULL)
1708 _componentsDescriptions = new string[NumberOfComponents];
1709 if (_componentsUnits == NULL)
1710 _componentsUnits = new UNIT[NumberOfComponents];
1711 if (_MEDComponentsUnits == NULL)
1712 _MEDComponentsUnits = new string[NumberOfComponents];
1713 for (int i=0;i<NumberOfComponents;i++) {
1714 _componentsTypes[i] = 0 ;
1718 _numberOfValues = _support->getNumberOfElements(MED_ALL_ELEMENTS);
1719 MESSAGE(LOC <<" : "<<_numberOfValues <<" et "<< NumberOfComponents);
1721 _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues);
1725 catch (MEDEXCEPTION &ex) {
1726 MESSAGE("No value defined, problem with NumberOfComponents (and may be _support) size of MEDARRAY<T>::_value !");
1727 _value = (MEDARRAY<T>*)NULL ;
1731 END_OF("void FIELD<T>::allocValue(const int NumberOfComponents)");
1737 template <class T> void FIELD<T>::allocValue(const int NumberOfComponents, const int LengthValue)
1739 BEGIN_OF("void FIELD<T>::allocValue(const int NumberOfComponents,const int LengthValue)");
1741 _numberOfComponents = NumberOfComponents ;
1742 if (_componentsTypes == NULL)
1743 _componentsTypes = new int[NumberOfComponents] ;
1744 if (_componentsNames == NULL)
1745 _componentsNames = new string[NumberOfComponents];
1746 if (_componentsDescriptions == NULL)
1747 _componentsDescriptions = new string[NumberOfComponents];
1748 if (_componentsUnits == NULL)
1749 _componentsUnits = new UNIT[NumberOfComponents];
1750 if (_MEDComponentsUnits == NULL)
1751 _MEDComponentsUnits = new string[NumberOfComponents];
1752 for (int i=0;i<NumberOfComponents;i++) {
1753 _componentsTypes[i] = 0 ;
1756 MESSAGE("FIELD : constructeur : "<<LengthValue <<" et "<< NumberOfComponents);
1757 _numberOfValues = LengthValue ;
1758 _value = new MEDARRAY<T>(_numberOfComponents,_numberOfValues);
1762 END_OF("void FIELD<T>::allocValue(const int NumberOfComponents,const int LengthValue)");
1768 template <class T> void FIELD<T>::deallocValue()
1770 BEGIN_OF("void FIELD<T>::deallocValue()");
1771 _numberOfValues = 0 ;
1772 _numberOfComponents = 0 ;
1776 END_OF("void FIELD<T>::deallocValue()");
1779 // -----------------
1781 // -----------------
1784 Create the specified driver and return its index reference to path to
1785 read or write methods.
1788 template <class T> int FIELD<T>::addDriver(driverTypes driverType,
1789 const string & fileName/*="Default File Name.med"*/,
1790 const string & driverName/*="Default Field Name"*/,
1791 MED_EN::med_mode_acces access)
1793 const char * LOC = "FIELD<T>::addDriver(driverTypes driverType, const string & fileName=\"Default File Name.med\",const string & driverName=\"Default Field Name\,MED_EN::med_mode_acces access) : ";
1801 driver = DRIVERFACTORY::buildDriverForField(driverType,fileName,this,access);
1803 _drivers.push_back(driver);
1805 int current = _drivers.size()-1;
1807 _drivers[current]->setFieldName(driverName);
1816 Duplicate the given driver and return its index reference to path to
1817 read or write methods.
1819 template <class T> inline int FIELD<T>::addDriver (GENDRIVER & driver )
1821 const char * LOC = "FIELD<T>::addDriver(GENDRIVER &) : ";
1826 // duplicate driver to delete it with destructor !
1827 GENDRIVER * newDriver = driver.copy() ;
1829 _drivers.push_back(newDriver);
1831 current = _drivers.size()-1;
1833 driver.setId(current);
1835 MESSAGE(LOC << " je suis la 1");
1837 MESSAGE(LOC << " je suis la 2");
1843 Remove the driver referenced by its index.
1845 template <class T> void FIELD<T>::rmDriver (int index/*=0*/)
1847 const char * LOC = "FIELD<T>::rmDriver (int index=0): ";
1850 if ( _drivers[index] ) {
1851 //_drivers.erase(&_drivers[index]);
1853 MESSAGE ("detruire");
1856 throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1857 << "The <index given is invalid, index must be between 0 and |"
1866 Read FIELD in the file specified in the driver given by its index.
1868 template <class T> inline void FIELD<T>::read(int index/*=0*/)
1870 const char * LOC = "FIELD<T>::read(int index=0) : ";
1873 if ( _drivers[index] ) {
1874 _drivers[index]->open();
1875 _drivers[index]->read();
1876 _drivers[index]->close();
1879 throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1880 << "The index given is invalid, index must be between 0 and |"
1888 Write FIELD in the file specified in the driver given by its index.
1890 template <class T> inline void FIELD<T>::write(int index/*=0*/, const string & driverName /*= ""*/)
1892 const char * LOC = "FIELD<T>::write(int index=0, const string & driverName = \"\") : ";
1895 if( _drivers[index] ) {
1896 _drivers[index]->open();
1897 if (driverName != "") _drivers[index]->setFieldName(driverName);
1898 _drivers[index]->write();
1899 _drivers[index]->close();
1902 throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1903 << "The index given is invalid, index must be between 0 and |"
1911 Write FIELD in the file specified in the driver given by its index. Use this
1912 method for ASCII drivers (e.g. VTK_DRIVER)
1914 template <class T> inline void FIELD<T>::writeAppend(int index/*=0*/, const string & driverName /*= ""*/)
1916 const char * LOC = "FIELD<T>::write(int index=0, const string & driverName = \"\") : ";
1919 if( _drivers[index] ) {
1920 _drivers[index]->openAppend();
1921 if (driverName != "") _drivers[index]->setFieldName(driverName);
1922 _drivers[index]->writeAppend();
1923 _drivers[index]->close();
1926 throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
1927 << "The index given is invalid, index must be between 0 and |"
1936 Write FIELD with the driver which is equal to the given driver.
1940 template <class T> inline void FIELD<T>::write(const GENDRIVER & genDriver)
1942 const char * LOC = " FIELD<T>::write(const GENDRIVER &) : ";
1945 for (unsigned int index=0; index < _drivers.size(); index++ )
1946 if ( *_drivers[index] == genDriver ) {
1947 _drivers[index]->open();
1948 _drivers[index]->write();
1949 _drivers[index]->close();
1958 Write FIELD with the driver which is equal to the given driver.
1960 Use by MED object. Use this method for ASCII drivers (e.g. VTK_DRIVER).
1962 template <class T> inline void FIELD<T>::writeAppend(const GENDRIVER & genDriver)
1964 const char * LOC = " FIELD<T>::write(const GENDRIVER &) : ";
1967 for (unsigned int index=0; index < _drivers.size(); index++ )
1968 if ( *_drivers[index] == genDriver ) {
1969 _drivers[index]->openAppend();
1970 _drivers[index]->writeAppend();
1971 _drivers[index]->close();
1980 Read FIELD with the driver which is equal to the given driver.
1984 template <class T> inline void FIELD<T>::read(const GENDRIVER & genDriver)
1986 const char * LOC = " FIELD<T>::read(const GENDRIVER &) : ";
1989 for (unsigned int index=0; index < _drivers.size(); index++ )
1990 if ( *_drivers[index] == genDriver ) {
1991 _drivers[index]->open();
1992 _drivers[index]->read();
1993 _drivers[index]->close();
2002 Destroy the MEDARRAY<T> in FIELD and put the new one without copy.
2005 template <class T> inline void FIELD<T>::setValue(MEDARRAY<T> *Value)
2007 if (NULL != _value) delete _value ;
2013 Return a reference to the MEDARRAY<T> in FIELD.
2016 template <class T> inline MEDARRAY<T>* FIELD<T>::getvalue() const
2022 Return the actual length of the reference to values array returned by getValue.
2024 template <class T> inline int FIELD<T>::getValueLength(MED_EN::medModeSwitch Mode) const{
2025 return _numberOfComponents*_numberOfValues;
2029 Return a reference to values array to read them.
2031 template <class T> inline const T* FIELD<T>::getValue(MED_EN::medModeSwitch Mode) const
2033 return _value->get(Mode) ;
2037 Return a reference to i^{th} row or column - component - (depend on Mode value)
2038 of FIELD values array.
2040 template <class T> inline const T* FIELD<T>::getValueI(MED_EN::medModeSwitch Mode,int i) const
2042 if ( Mode == MED_EN::MED_FULL_INTERLACE )
2044 return _value->getRow(i) ;
2046 ASSERT ( Mode == MED_EN::MED_NO_INTERLACE);
2047 return _value->getColumn(i);
2051 Return the value of i^{th} element and j^{th} component.
2053 template <class T> inline T FIELD<T>::getValueIJ(int i,int j) const
2055 return _value->getIJ(i,j) ;
2059 Copy new values array in FIELD according to the given mode.
2061 Array must have right size. If not results are unpredicable.
2063 template <class T> inline void FIELD<T>::setValue(MED_EN::medModeSwitch mode, T* value)
2065 _value->set(mode,value);
2069 Update values array in FIELD with the given ones according to specified mode.
2071 template <class T> inline void FIELD<T>::setValueI(MED_EN::medModeSwitch mode, int i, T* value)
2074 if (MED_EN::MED_FULL_INTERLACE == mode)
2075 _value->setI(i,value);
2076 else if (MED_EN::MED_NO_INTERLACE == mode)
2077 _value->setJ(i,value);
2079 throw MEDEXCEPTION(LOCALIZED("FIELD<T>::setValueI : bad medModeSwitch")) ;
2083 Set the value of i^{th} element and j^{th} component with the given one.
2085 template <class T> inline void FIELD<T>::setValueIJ(int i, int j, T value)
2087 _value->setIJ(i,j,value);
2095 Fill values array with volume values.
2097 template <class T> void FIELD<T>::getVolume() const throw (MEDEXCEPTION)
2099 const char * LOC = "FIELD<double>::getVolume() const : ";
2102 // The field has to be initilised by a non empty support and a
2103 // number of components = 1 and its value type has to be set to MED_REEL64
2104 // (ie a FIELD<double>)
2106 if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2107 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to 1 and a value type set to MED_REEL64"));
2113 Fill values array with area values.
2115 template <class T> void FIELD<T>::getArea() const throw (MEDEXCEPTION)
2117 const char * LOC = "FIELD<double>::getArea() const : ";
2120 // The field has to be initilised by a non empty support and a
2121 // number of components = 1 and its value type has to be set to MED_REEL64
2122 // (ie a FIELD<double>)
2124 if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2125 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to 1 and a value type set to MED_REEL64"));
2131 Fill values array with length values.
2133 template <class T> void FIELD<T>::getLength() const throw (MEDEXCEPTION)
2135 const char * LOC = "FIELD<double>::getLength() const : ";
2138 // The field has to be initilised by a non empty support and a
2139 // number of components = 1 and its value type has to be set to MED_REEL64
2140 // (ie a FIELD<double>)
2142 if ((_support == (SUPPORT *) NULL) || (_numberOfComponents != 1) || (_valueType != MED_REEL64))
2143 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to 1 and a value type set to MED_REEL64"));
2149 Fill values array with normal values.
2151 template <class T> void FIELD<T>::getNormal() const throw (MEDEXCEPTION)
2153 const char * LOC = "FIELD<double>::getNormal() const : ";
2156 // The field has to be initilised by a non empty support and a
2157 // number of components = 1 and its value type has to be set to MED_REEL64
2158 // (ie a FIELD<double>)
2160 if (_support == (SUPPORT *) NULL)
2161 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2163 int dim_space = _support->getMesh()->getSpaceDimension();
2165 if ((_numberOfComponents != dim_space) || (_valueType != MED_REEL64))
2166 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2172 Fill values array with barycenter values.
2174 template <class T> void FIELD<T>::getBarycenter() const throw (MEDEXCEPTION)
2176 const char * LOC = "FIELD<double>::getBarycenter() const : ";
2179 // The field has to be initilised by a non empty support and a number of
2180 //components = space dimension and its value type has to be set to MED_REEL64
2181 // (ie a FIELD<double>)
2183 if (_support == (SUPPORT *) NULL)
2184 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2186 int dim_space = _support->getMesh()->getSpaceDimension();
2188 if ((_numberOfComponents != dim_space) || (_valueType != MED_REEL64))
2189 throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"The field has to be initialised with a non empty support, a number of components set to the space dimension and a value type set to MED_REEL64"));
2194 }//End namespace MEDMEM
2196 #endif /* FIELD_HXX */