ElementVect.h

Go to the documentation of this file.
00001 /*==============================================================================
00002 
00003                                     O  F  E  L  I
00004 
00005                            Object  Finite  Element  Library
00006 
00007   ==============================================================================
00008 
00009    Copyright (C) 1998 - 2008 Rachid Touzani
00010 
00011    This program is free software; you can redistribute it and/or modify it under
00012    the terms of the GNU General Public License as published by the Free
00013    Software Foundation; Version 2 of the License.
00014 
00015    This program is distributed in the hope that it will be useful, but WITHOUT
00016    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00017    FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
00018    details.
00019 
00020    You should have received a copy of the GNU General Public License
00021    along with this program; if not, write to the :
00022 
00023    Free Software Foundation
00024    Inc., 59 Temple Place - Suite 330
00025    Boston, MA  02111-1307, USA
00026 
00027   ==============================================================================
00028 
00029              Definition of class 'ElementVect' for ElementWise Vectors
00030 
00031   ==============================================================================*/
00032 
00033 
00034 #ifndef __ELEMENT_VECT_H
00035 #define __ELEMENT_VECT_H
00036 
00037 #include "AbsVect.h"
00038 #include "Vect.h"
00039 
00040 namespace OFELI {
00041 
00046 template <class T_> class Vect;
00047 template <class T_> class AbsVect;
00048 
00049 enum {
00050    EL_CENTER_DOF  =  1,
00051    EL_NODE_DOF    =  2,
00052    EL_SIDE_DOF    =  3
00053 };
00054 
00055 
00056 //-----------------------------------------------------------------------------
00057 // Class ErrorInElementVect
00058 //-----------------------------------------------------------------------------
00059 
00060 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00061 template<class T_>
00062 struct ErrorInElementVect {
00063    void Message(const char *file, int line, int code, int p1=0) {
00064      cerr << "\n*** Fatal Error in OFELI ***\n";
00065      cerr << "File : " << file << ", line : " << line << "\nIn ElementVect<>::";
00066      switch (code) {
00067 
00068         case 21:
00069            cerr << "ElementVect(mesh,nb_dof) : Illegal value of nb_dof : " << p1 << endl;
00070            break;
00071 
00072         case 31:
00073            cerr << "ElementVect(mesh,_nb_dof,name,time) : Illegal value of nb_dof : " << p1 << endl;
00074            break;
00075 
00076         case 121:
00077            cerr << "Get(ff) : Nb. of elements in mesh and input file do not match\n";
00078            break;
00079 
00080         case 131:
00081            cerr << "Get(ff,time) : Nb. of elements in mesh and input file do not match\n";
00082            break;
00083 
00084         case 141:
00085            cerr << "View(ff) : Nb. of elements in mesh and input file do not match\n";
00086            exit(1);
00087      }
00088      cerr << "Program stops" << endl;
00089      exit(1);
00090    }
00091 };
00092 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00093 
00094 
00108 template <class T_>
00109 class ElementVect : virtual public AbsVect<T_>
00110 {
00111 
00112 public:
00113 
00114    using AbsVect<T_>::Init;
00115    using AbsVect<T_>::_kw;
00116    using AbsVect<T_>::getName;
00117    using AbsVect<T_>::setName;
00118 
00120     ElementVect();
00121 
00131     ElementVect(const class Mesh &mesh, size_t nb_dof=1, size_t type=EL_CENTER_DOF);
00132 
00144     ElementVect(const class Mesh &mesh, size_t nb_dof, const char *name, double time=0, size_t type=EL_CENTER_DOF);
00145 
00148     ElementVect(const ElementVect<T_> &v);
00149 
00151     ~ElementVect();
00152 
00156     void setMesh(const class Mesh &ms, size_t nb_dof=1);
00157 
00167     void setMesh(const class Mesh &mesh, size_t nb_dof, size_t type);
00168 
00171     int get(const string &file);
00172 
00177     int get(const string &file, double time);
00178 
00179 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00180     int Get(const string &file) { return get(file); }
00181     int Get(const string &file, double time) { return get(file,time); }
00182 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00183 
00185     size_t getNbElements() const;
00186 
00189     ElementVect<T_> & operator=(T_ a);
00190 
00195     ElementVect<T_> & operator+=(const ElementVect<T_> &x);
00196 
00201     ElementVect<T_> & operator+=(const T_ &a);
00202 
00205     ElementVect<T_> & operator-=(const ElementVect<T_> &x);
00206 
00211     ElementVect<T_> & operator-=(const T_ &a);
00212 
00215     ElementVect<T_> &operator*=(const T_ &a);
00216 
00219     ElementVect<T_> &operator/=(const T_ &a);
00220 
00221    using AbsVect<T_>::_nbn;
00222    using AbsVect<T_>::_nb;
00223    using AbsVect<T_>::_nb_dof;
00224    using AbsVect<T_>::_first_dof;
00225    using AbsVect<T_>::_length;
00226    using AbsVect<T_>::_time;
00227    using AbsVect<T_>::_name;
00228    using AbsVect<T_>::_theMesh;
00229    using AbsVect<T_>::_v;
00230 
00231 protected:
00232 
00233 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00234    ErrorInElementVect<T_> _e;
00235    size_t                 _nb_elements, _type;
00236 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00237 };
00238 
00239 
00240 
00242 //                         I M P L E M E N T A T I O N                       //
00244 
00245 
00246 template<class T_>
00247 ElementVect<T_>::ElementVect()
00248 {
00249 #ifdef _OFELI_DEBUG
00250    std::clog << "An instance of class ElementVect is constructed.\n";
00251    std::clog << "File: " << __FILE__ << ", Line: " << __LINE__ << endl;
00252 #endif
00253 }
00254 
00255 
00256 template<class T_>
00257 ElementVect<T_>::ElementVect(const class Mesh &mesh, size_t nb_dof, size_t type)
00258 {
00259    setMesh(mesh,nb_dof,type);
00260 #ifdef _OFELI_DEBUG
00261    std::clog << "An instance of class ElementVect is constructed.\n";
00262    std::clog << "File: " << __FILE__ << ", Line: " << __LINE__ << endl;
00263 #endif
00264 }
00265 
00266 
00267 template<class T_>
00268 ElementVect<T_>::ElementVect(const class Mesh &mesh, size_t nb_dof, const char *name, double time, size_t type)
00269                 : AbsVect<T_>(mesh,nb_dof,name,time)
00270 {
00271    setMesh(mesh,nb_dof,type);
00272    Init();
00273 #ifdef _OFELI_DEBUG
00274    std::clog << "An instance of class ElementVect is constructed.\n";
00275    std::clog << "File: " << __FILE__ << ", Line: " << __LINE__ << endl;
00276 #endif
00277 }
00278 
00279 
00280 template<class T_>
00281 ElementVect<T_>::ElementVect(const ElementVect<T_> &v) : AbsVect<T_>(v)
00282 {
00283    _time = v._time;
00284    _name = v._name;
00285    _length = v.getLength();
00286    _nbn = v._nbn;
00287    _theMesh = v._theMesh;
00288    _first_dof = 1;
00289    _nb_dof = v._nb_dof;
00290    Init();
00291    _v = v._v;
00292 #ifdef _OFELI_DEBUG
00293    std::clog << "An instance of class ElementVect is constructed.\n";
00294    std::clog << "File: " << __FILE__ << ", Line: " << __LINE__ << endl;
00295 #endif
00296 }
00297 
00298 
00299 template<class T_>
00300 ElementVect<T_>::~ElementVect()
00301 {
00302 #ifdef _OFELI_DEBUG
00303    std::clog << "An instance of class ElementVect is destructed.\n";
00304    std::clog << "File: " << __FILE__ << ", Line: " << __LINE__ << endl;
00305 #endif
00306 }
00307 
00308 
00309 template<class T_>
00310 void ElementVect<T_>::setMesh(const class Mesh &ms, size_t nb_dof)
00311 {
00312    if (nb_dof == 0)
00313       _e.Message(__FILE__,__LINE__,21,int(nb_dof));
00314    _theMesh = &ms;
00315    _first_dof = 1;
00316    _nbn = 1;
00317    _nb_dof = nb_dof;
00318    _nb_elements = ms.getNbElements();
00319    _nb = _nb_elements*_nbn;
00320    _length = _nb_dof * _nb;
00321    _type = EL_CENTER_DOF;
00322    Init();
00323 }
00324 
00325 
00326 template<class T_>
00327 void ElementVect<T_>::setMesh(const class Mesh &ms, size_t nb_dof, size_t type)
00328 {
00329    if (nb_dof == 0)
00330       _e.Message(__FILE__,__LINE__,21,int(nb_dof));
00331    _theMesh = &ms;
00332    _first_dof = 1;
00333    _nb_dof = nb_dof;
00334    if (type==EL_CENTER_DOF)
00335       _nbn = 1;
00336    else if (type==EL_NODE_DOF)
00337       _nbn = ms.getPtrElement(1)->getNbNodes();
00338    else
00339       _nbn = ms.getPtrElement(1)->getNbSides();
00340    _nb_elements = ms.getNbElements();
00341    _nb = _nb_elements*_nbn;
00342    _length = _nb_dof * _nb;
00343    _type = type;
00344    Init();
00345 }
00346 
00347 
00348 template<class T_>
00349 int ElementVect<T_>::get(const string &file)
00350 {
00351    _kw[3] = "NbElements$";
00352    size_t nb_elements, i, j, n;
00353    int key;
00354    FFI ff(file,"#ELEMENT_FIELD!");
00355    ff.setKeywords(_kw);
00356    do {
00357       switch (key=ff.getKW()) {
00358 
00359          case -1:
00360             exit(1);
00361 
00362          case 1:
00363             setName(ff.getS());
00364             break;
00365 
00366          case 3:
00367             nb_elements = ff.getI();
00368             if (nb_elements != _nb_elements)
00369                _e.Message(__FILE__,__LINE__,121);
00370             break;
00371 
00372          case 2:
00373             _nb_dof = ff.getI();
00374             _length = _nb_dof * _nb_elements;
00375             break;
00376 
00377          case 4:
00378             _time = ff.getD();
00379             break;
00380 
00381          case 5:
00382             _type = ff.getI();
00383             break;
00384 
00385          case 6:
00386             Init();
00387             for (j=0; j<_nb_elements; ++j) {
00388                n = ff.getI();
00389                for (i=0; i<_nb_dof; i++)
00390                   for (size_t k=0; k<1; ++k)
00391                      _v[_nb_dof*(n-1)+_nbn*i+k] = ff.getD();
00392             }
00393             break;
00394       }
00395    } while (key);
00396    ff.Close();
00397    return 0;
00398 }
00399 
00400 
00401 template<class T_>
00402 int ElementVect<T_>::get(const string &file, double time)
00403 {
00404    _kw[3] = "NbElements$";
00405    size_t nb_elements, i, j, n, k;
00406    int key, time_found=0;
00407    FFI ff(file,"#ELEMENT_FIELD!");
00408    ff.setKeywords(_kw);
00409    do {
00410       switch (key=ff.getKW()) {
00411 
00412          case -1:
00413             exit(1);
00414 
00415          case 1:
00416             setName(ff.getS());
00417             break;
00418 
00419          case 3:
00420             nb_elements = ff.getI();
00421             if (nb_elements != _nb_elements)
00422                _e.Message(__FILE__,__LINE__,121);
00423             break;
00424 
00425          case 2:
00426             _nb_dof = ff.getI();
00427             _length = _nb_dof * _nb_elements;
00428             break;
00429 
00430          case 4:
00431             _time = ff.getD();
00432             if (_time == time)
00433                time_found = 1;
00434             break;
00435 
00436          case 5:
00437             _type = ff.getI();
00438             break;
00439 
00440          case 6:
00441             Init();
00442             for (j=0; j<_nb_elements; ++j) {
00443                n = ff.getI();
00444                for (i=0; i<_nb_dof; ++i)
00445                   for (k=0; k<1; ++k)
00446                      _v[_nb_dof*(n-1)+_nbn*i+k] = ff.getD();
00447             }
00448             if (time_found==1) {
00449                ff.Close();
00450                return 0;
00451             }
00452             break;
00453       }
00454    } while (key);
00455    ff.Close();
00456    return 1;
00457 }
00458 
00459 
00460 template<class T_>
00461 size_t ElementVect<T_>::getNbElements() const
00462 {
00463    return _nb_elements;
00464 }
00465 
00466 
00467 template<class T_>
00468 ElementVect<T_> & ElementVect<T_>::operator=(T_ a)
00469 {
00470   for (size_t i=0; i<_length; ++i)
00471      _v[i] = a;
00472   return *this;
00473 }
00474 
00475 
00476 template <class T_>
00477 ElementVect<T_> & ElementVect<T_>::operator+=(const ElementVect<T_> &x)
00478 {
00479    for (size_t i=0; i<_length; ++i)
00480       _v[i] += x._v[i];
00481    return *this;
00482 }
00483 
00484 
00485 template <class T_>
00486 ElementVect<T_> & ElementVect<T_>::operator+=(const T_ &a)
00487 {
00488    for (size_t i=0; i<_length; ++i)
00489       _v[i] += a;
00490    return *this;
00491 }
00492 
00493 
00494 template <class T_>
00495 ElementVect<T_> & ElementVect<T_>::operator-=(const ElementVect<T_> &x)
00496 {
00497    for (size_t i=0; i<_length; ++i)
00498       _v[i] -= x._v[i];
00499    return *this;
00500 }
00501 
00502 
00503 template <class T_>
00504 ElementVect<T_> & ElementVect<T_>::operator-=(const T_ &a)
00505 {
00506    for (size_t i=0; i<_length; ++i)
00507       _v[i] -= a;
00508    return *this;
00509 }
00510 
00511 
00512 template <class T_>
00513 ElementVect<T_> & ElementVect<T_>::operator*=(const T_ &a)
00514 {
00515    for (size_t i=0; i<_length; ++i)
00516       _v[i] *= a;
00517    return *this;
00518 }
00519 
00520 
00521 template <class T_>
00522 ElementVect<T_> & ElementVect<T_>::operator/=(const T_ &a)
00523 {
00524    for (size_t i=0; i<_length; ++i)
00525       _v[i] /= a;
00526    return *this;
00527 }
00528 
00529 
00531 //                             Associated Functions                          //
00533 
00537 template <class T_>
00538 istream& operator>>(istream& s, ElementVect<T_> &v)
00539 {
00540    size_t j, n, type;
00541    double time;
00542    static char name[40];
00543    s >> name;
00544    v.setName(name);
00545    s >> n >> time >> type;
00546    v.Time(time);
00547 
00548    for (size_t i=1; i<=n; i++) {
00549       s >> j;
00550       for (size_t k=1; k<=v.NbDOF(); k++)
00551          s >> v(i,j) >> "    ";
00552    }
00553    return s;
00554 }
00555 
00556 
00560 template <class T_>
00561 ostream& operator<<(ostream& s, const ElementVect<T_> &v)
00562 {
00563    const Element *el;
00564    s << v.getName() << " at time = " << v.getTime() << endl << endl;
00565    s.setf(ios::scientific);
00566    for (v.getMesh().topElement(); (el=v.getMesh().getElement());) {
00567       s << setw(6) << el->getLabel() << "   ";
00568       for (size_t j=1; j<=v.getNbDOF(); j++)
00569          s << setprecision(8) << setw(18) << v(el->getLabel(),j);
00570       s << endl;
00571    }
00572    s << endl;
00573    return s;
00574 }
00575 
00576 } /* namespace OFELI */
00577 
00578 #endif