LocalVect.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 Template Class LocalVect for small size vectors
00030 
00031   ==============================================================================*/
00032 
00033 
00034 #ifndef __LOCAL_VECT_H
00035 #define __LOCAL_VECT_H
00036 
00037 #include "OFELI_Config.h"
00038 #include "Element.h"
00039 #include "Side.h"
00040 
00041 namespace OFELI {
00042 
00075 template<class T_> class Vect;
00076 template<class T_> class NodeVect;
00077 
00078 template<class T_,size_t N_>
00079 class LocalVect
00080 {
00081 
00082 public:
00083 
00086     LocalVect();
00087 
00089     LocalVect(const T_ *a);
00090 
00092     LocalVect(const Element *el);
00093 
00095     LocalVect(const Side *sd);
00096 
00098     LocalVect(const LocalVect<T_,N_> &v);
00099 
00110     LocalVect(const Element *el, const Vect<T_> &v, int opt=0);
00111 
00122     LocalVect(const Element *el, const NodeVect<T_> &v, int opt=0);
00123 
00134     LocalVect(const Side *sd, const Vect<T_> &v, int opt=0);
00135 
00137     ~LocalVect();
00138 
00146     void Localize(const Element *el, const Vect<T_> &v);
00147 
00155     void Localize(const Side *sd, const Vect<T_> &v);
00156 
00159     T_ &operator[](size_t i) { return _v[i]; }
00160 
00163     T_ operator[](size_t i) const { return _v[i]; }
00164 
00168     T_ &operator()(size_t i) { return _v[i-1]; }
00169 
00173     T_ operator()(size_t i) const { return _v[i-1]; }
00174 
00177     Element *El() { return _el; }
00178 
00181     Side *Sd() { return _sd; }
00182 
00185     LocalVect<T_,N_> & operator=(const LocalVect<T_,N_> &v);
00186 
00189     LocalVect<T_,N_> & operator=(const T_ &x);
00190 
00193     LocalVect<T_,N_> & operator+=(const LocalVect<T_,N_> &v);
00194 
00197     LocalVect<T_,N_> & operator+=(const T_ &a);
00198 
00201     LocalVect<T_,N_> & operator-=(const LocalVect<T_,N_> &v);
00202 
00205     LocalVect<T_,N_> & operator-=(const T_ &a);
00206 
00209     LocalVect<T_,N_> & operator*=(const T_ &a);
00210 
00213     LocalVect<T_,N_> & operator/=(const T_ &a);
00214 
00216     T_ *getArray() { return _v; }
00217 
00218  private:
00219 
00220    T_      _v[N_];
00221    Element *_el;
00222    Side    *_sd;
00223 };
00224 
00226 //                         I M P L E M E N T A T I O N                       //
00228 
00229 
00230 template<class T_,size_t N_>
00231 LocalVect<T_,N_>::LocalVect()
00232 {
00233    for (size_t i=0; i<N_; ++i)
00234       _v[i] = 0;
00235    _el = NULL;
00236    _sd = NULL;
00237 #ifdef _OFELI_DEBUG
00238    std::clog << "An instance of class LocalVect is constructed.\n";
00239    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00240 #endif
00241 }
00242 
00243 
00244 template<class T_,size_t N_>
00245 LocalVect<T_,N_>::LocalVect(const T_ *a)
00246 {
00247    for (size_t i=0; i<N_; ++i)
00248       _v[i] = a[i];
00249    _el = NULL;
00250    _sd = NULL;
00251 #ifdef _OFELI_DEBUG
00252    std::clog << "An instance of class LocalVect is constructed.\n";
00253    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00254 #endif
00255 }
00256 
00257 
00258 template<class T_,size_t N_>
00259 LocalVect<T_,N_>::LocalVect(const Element *el)
00260 {
00261    for (size_t i=0; i<N_; ++i)
00262       _v[i] = 0;
00263    _el = el;
00264    _sd = NULL;
00265 #ifdef _OFELI_DEBUG
00266    std::clog << "An instance of class LocalVect is constructed.\n";
00267    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00268 #endif
00269 }
00270 
00271 
00272 template<class T_,size_t N_>
00273 LocalVect<T_,N_>::LocalVect(const Side *sd)
00274 {
00275    for (size_t i=0; i<N_; ++i)
00276       _v[i] = 0;
00277    _el = NULL;
00278    _sd = sd;
00279 #ifdef _OFELI_DEBUG
00280    std::clog << "An instance of class LocalVect is constructed.\n";
00281    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00282 #endif
00283 }
00284 
00285 
00286 template<class T_,size_t N_>
00287 LocalVect<T_,N_>::LocalVect(const LocalVect<T_,N_> &v)
00288 {
00289    for (size_t i=0; i<N_; ++i)
00290       _v[i] = v[i];
00291    _el = v._el;
00292    _sd = v._sd;
00293 #ifdef _OFELI_DEBUG
00294    std::clog << "An instance of class LocalVect is constructed.\n";
00295    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00296 #endif
00297 }
00298 
00299 
00300 template<class T_,size_t N_>
00301 LocalVect<T_,N_>::LocalVect(const Element *el, const Vect<T_> &v, int opt)
00302 {
00303    if (opt==0)
00304       Localize(el,v);
00305    else {
00306       size_t i = 0;
00307       for (size_t n=1; n<=el->getNbNodes(); ++n)
00308          _v[i++] = v(el->getNodeLabel(n));
00309    }
00310 #ifdef _OFELI_DEBUG
00311    std::clog << "An instance of class LocalVect is constructed.\n";
00312    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00313 #endif
00314 }
00315 
00316 
00317 template<class T_,size_t N_>
00318 LocalVect<T_,N_>::LocalVect(const Element *el, const NodeVect<T_> &v, int opt)
00319 {
00320    size_t i;
00321    if (opt==0)
00322       for (i=0; i<el->getNbNodes(); ++i)
00323          _v[i] = v(el->getNodeLabel(i+1),1);
00324    else
00325       for (i=0; i<el->getNbNodes(); ++i)
00326          _v[i] = v(el->getNodeLabel(i+1),opt);
00327 #ifdef _OFELI_DEBUG
00328    std::clog << "An instance of class LocalVect is constructed.\n";
00329    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00330 #endif
00331 }
00332 
00333 
00334 template<class T_,size_t N_>
00335 LocalVect<T_,N_>::LocalVect(const Side *sd, const Vect<T_> &v, int opt)
00336 {
00337    if (opt==0)
00338       Localize(sd,v);
00339    else {
00340       size_t i = 0;
00341       for (size_t n=1; n<=sd->getNbNodes(); ++n)
00342          _v[i++] = v(sd->getNodeLabel(n));
00343    }
00344 #ifdef _OFELI_DEBUG
00345    std::clog << "An instance of class LocalVect is constructed.\n";
00346    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00347 #endif
00348 }
00349 
00350 
00351 template<class T_,size_t N_>
00352 LocalVect<T_,N_>::~LocalVect()
00353 {
00354 #ifdef _OFELI_DEBUG
00355    std::clog << "An instance of class LocalVect is destructed.\n";
00356    std::clog << "File : " << __FILE__ << ", Line : " << __LINE__ << endl;
00357 #endif
00358 }
00359 
00360 
00361 template<class T_,size_t N_>
00362 void LocalVect<T_,N_>::Localize(const Element *el, const Vect<T_> &v)
00363 {
00364    size_t i=0, k;
00365    for (size_t n=1; n<=el->getNbNodes(); ++n) {
00366       Node *nd = el->getPtrNode(n);
00367       k = nd->getFirstDOF();
00368       for (size_t j=1; j<=nd->getNbDOF(); ++j)
00369          _v[i++] = v(k++);
00370    }
00371 }
00372 
00373 
00374 template<class T_,size_t N_>
00375 void LocalVect<T_,N_>::Localize(const Side *sd, const Vect<T_> &v)
00376 {
00377    size_t i = 0;
00378    for (size_t n=1; n<=sd->getNbNodes(); ++n) {
00379       Node *nd = sd->getPtrNode(n);
00380       for (size_t j=1; j<=nd->getNbDOF(); ++j)
00381          _v[i++] = v(nd->getDOF(j));
00382    }
00383 }
00384 
00385 
00386 template<class T_,size_t N_>
00387 LocalVect<T_,N_> & LocalVect<T_,N_>::operator=(const LocalVect<T_,N_> &v)
00388 {
00389    for (size_t i=0; i<N_; ++i)
00390       _v[i] = v._v[i];
00391    return *this;
00392 }
00393 
00394 
00395 template<class T_,size_t N_>
00396 LocalVect<T_,N_> & LocalVect<T_,N_>::operator=(const T_ &x)
00397 {
00398    for (size_t i=0; i<N_; ++i)
00399       _v[i] = x;
00400    return *this;
00401 }
00402 
00403 
00404 template<class T_,size_t N_>
00405 LocalVect<T_,N_> & LocalVect<T_,N_>::operator+=(const LocalVect<T_,N_> &v)
00406 {
00407    for (size_t i=0; i<N_; ++i)
00408       _v[i] += v[i];
00409    return *this;
00410 }
00411 
00412 
00413 template<class T_,size_t N_>
00414 LocalVect<T_,N_> & LocalVect<T_,N_>::operator+=(const T_ &a)
00415 {
00416    for (size_t i=0; i<N_; ++i)
00417       _v[i] += a;
00418    return *this;
00419 }
00420 
00421 
00422 template<class T_,size_t N_>
00423 LocalVect<T_,N_> & LocalVect<T_,N_>::operator-=(const LocalVect<T_,N_> &v)
00424 {
00425    for (size_t i=0; i<N_; ++i)
00426       _v[i] -= v[i];
00427    return *this;
00428 }
00429 
00430 
00431 template<class T_,size_t N_>
00432 LocalVect<T_,N_> & LocalVect<T_,N_>::operator-=(const T_ &a)
00433 {
00434    for (size_t i=0; i<N_; ++i)
00435       _v[i] -= a;
00436    return *this;
00437 }
00438 
00439 
00440 template<class T_,size_t N_>
00441 LocalVect<T_,N_> & LocalVect<T_,N_>::operator*=(const T_ &a)
00442 {
00443    for (size_t i=0; i<N_; ++i)
00444       _v[i] *= a;
00445    return *this;
00446 }
00447 
00448 
00449 template<class T_,size_t N_>
00450 LocalVect<T_,N_> & LocalVect<T_,N_>::operator/=(const T_ &a)
00451 {
00452    for (size_t i=0; i<N_; ++i)
00453       _v[i] /= a;
00454    return *this;
00455 }
00456 
00457             
00459 //                A S S O C I A T E D    F U N C T I O N S                   //
00461 
00467 template<class T_,size_t N_>
00468 LocalVect<T_,N_> operator+(const LocalVect<T_,N_> &x, const LocalVect<T_,N_> &y)
00469 {
00470    LocalVect<T_,N_> z(x);
00471    for (size_t i=0; i<N_; ++i)
00472       z[i] += y[i];
00473    return z;
00474 }
00475 
00481 template<class T_,size_t N_>
00482 LocalVect<T_,N_> operator-(const LocalVect<T_,N_> &x, const LocalVect<T_,N_> &y)
00483 {
00484    LocalVect<T_,N_> z(x);
00485    for (size_t i=0; i<N_; ++i)
00486       z[i] -= y[i];
00487    return z;
00488 }
00489 
00490 
00496 template<class T_,size_t N_>
00497 LocalVect<T_,N_> operator*(T_ a, const LocalVect<T_,N_> &x)
00498 {
00499    LocalVect<T_,N_> v(x);
00500    for (size_t i=0; i<N_; ++i)
00501       v[i] *= a;
00502    return v;
00503 }
00504 
00505 
00511 template<class T_,size_t N_>
00512 LocalVect<T_,N_> operator/(T_ a, const LocalVect<T_,N_> &x)
00513 {
00514    LocalVect<T_,N_> v(x);
00515    for (size_t i=0; i<N_; ++i)
00516       v[i] /= a;
00517    return v;
00518 }
00519 
00520 
00522 //                       T H E   M E T A P R O G R A M S                     //
00524 
00525 // Dot Product
00526 
00527 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00528 template<size_t I_>
00529 struct _meta_dot
00530 {
00531   template<class T_, size_t N_>
00532   static T_ _dot(const LocalVect<T_,N_> &a, const LocalVect<T_,N_> &b)
00533   { return a[I_]*b[I_] + _meta_dot<I_-1>::_dot(a,b); }
00534 };
00535 
00536 
00537 template<>
00538 struct _meta_dot<0>
00539 {
00540   template<class T_, size_t N_>
00541   static T_ _dot(const LocalVect<T_,N_> &a, const LocalVect<T_,N_> &b)
00542   { return a[0]*b[0]; }
00543 };
00544 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00545 
00551 template<class T_, size_t N_>
00552 inline double Dot(const LocalVect<T_,N_> &a, const LocalVect<T_,N_> &b)
00553 {
00554    return _meta_dot<N_-1>::_dot(a,b);
00555 }
00556 
00557 
00558 // scale
00559 
00560 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00561 template<size_t I_>
00562 struct _meta_scale
00563 {
00564   template<class T_, size_t N_>
00565   static void _scale(T_ a, const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00566   {
00567      y[I_] = a*x[I_];
00568      _meta_scale<I_-1>::_scale(a,x,y);
00569   }
00570 };
00571 
00572 template<>
00573 struct _meta_scale<0>
00574 {
00575   template<class T_, size_t N_>
00576   static void _scale(T_ a, const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00577   {
00578      y[0] = a*x[0];
00579   }
00580 };
00581 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00582 
00583 
00587 template<class T_, size_t N_>
00588 inline void Scale(T_ a, const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00589 {
00590    _meta_scale<N_-1>::_scale(a,x,y);
00591 }
00592 
00593 
00594 // scale (in place)
00595 
00596 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00597 template<size_t I_>
00598 struct _meta_scale1
00599 {
00600   template<class T_, size_t N_>
00601   static void _scale(T_ a, LocalVect<T_,N_> &x)
00602   { x[I_] *= a; _meta_scale1<I_-1>::_scale(a,x); }
00603 };
00604 
00605 template<>
00606 struct _meta_scale1<0>
00607 {
00608   template<class T_, size_t N_>
00609   static void scale_(T_ a, LocalVect<T_,N_> &x)
00610   { x[0] *= a; }
00611 };
00612 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00613 
00617 template<class T_, size_t N_>
00618 inline void Scale(T_ a, LocalVect<T_,N_> &x)
00619 {
00620    _meta_scale1<N_-1>::_scale(a,x);
00621 }
00622 
00623 // axpy
00624 
00625 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00626 template<size_t I_>
00627 struct _meta_axpy
00628 {
00629   template<class T_, size_t N_>
00630   static void _axpy(T_ a, const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00631   { y[I_] += a*x[I_]; _meta_axpy<I_-1>::_axpy(a,x,y); }
00632 };
00633 
00634 
00635 template<>
00636 struct _meta_axpy<0>
00637 {
00638   template<class T_, size_t N_>
00639   static void _axpy(T_ a, const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00640   { y[0] += a*x[0]; }
00641 };
00642 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00643 
00647 template<class T_, size_t N_>
00648 inline void Axpy(T_ a, const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00649 {
00650    _meta_axpy<N_-1>::_axpy(a,x,y);
00651 }
00652 
00653 
00654 // xpy
00655 
00656 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00657 template<size_t I_>
00658 struct _meta_xpy
00659 {
00660   template<class T_, size_t N_>
00661   void xpy_(const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00662   { y[I_] += x[I_]; _meta_axpy<I_-1>::_xpy(x,y); }
00663 };
00664 
00665 
00666 template<>
00667 struct _meta_xpy<0>
00668 {
00669   template<class T_, size_t N_>
00670   void _xpy(const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00671   { y[0] += x[0]; }
00672 };
00673 
00674 
00678 template<class T_, size_t N_>
00679 inline void Xpy(const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00680 {
00681    _meta_xpy<N_-1>::_xpy(x,y);
00682 }
00683 
00684 // copy
00685 
00686 template<size_t I_>
00687 struct _meta_copy
00688 {
00689   template<class T_, size_t N_>
00690   void _copy(const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00691   { y[I_] = x[I_]; _meta_axpy<I_-1>::_xpy(x,y); }
00692 };
00693 
00694 template<>
00695 struct _meta_copy<0>
00696 {
00697   template<class T_, size_t N_>
00698   void _copy(const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00699   { y[0] = x[0]; }
00700 };
00701 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
00702 
00703 
00707 template<class T_, size_t N_>
00708 inline void Copy(const LocalVect<T_,N_> &x, LocalVect<T_,N_> &y)
00709 {
00710    _meta_copy<N_-1>::_copy(x,y);
00711 }
00712 
00713 
00717 template<class T_,size_t N_>
00718 ostream& operator<<(ostream& s, const LocalVect<T_,N_> & v)
00719 {
00720    s.setf(ios::scientific);
00721    for (size_t i=1; i<=N_; i++)
00722       s << setw(6) << i << "  " << setprecision(8) << setw(18) << v(i) << endl;
00723    return s;
00724 }
00725 
00726 } /* namespace OFELI */
00727 
00728 #endif