fvector.hh

Go to the documentation of this file.
00001 // $Id: fvector.hh 5972 2010-04-12 08:13:13Z mdroh_01 $
00002 #ifndef DUNE_FVECTOR_HH
00003 #define DUNE_FVECTOR_HH
00004 
00005 #include<cmath>
00006 #include<cstddef>
00007 #include<cstdlib>
00008 #include<complex>
00009 #include<cstring>
00010 #include<limits>
00011 
00012 #include "exceptions.hh"
00013 #include "genericiterator.hh"
00014 
00015 #include "ftraits.hh"
00016 
00017 #ifdef DUNE_EXPRESSIONTEMPLATES
00018 #include "exprtmpl.hh"
00019 #endif
00020 
00021 namespace Dune {
00022 
00023   // forward declaration of template
00024   template<class K, int SIZE> class FieldVector;
00025 
00026   template<class K, int SIZE>
00027   struct FieldTraits< FieldVector<K,SIZE> >
00028   {
00029     typedef const typename FieldTraits<K>::field_type field_type;
00030     typedef const typename FieldTraits<K>::real_type real_type;
00031   };
00032 
00043 #ifndef DUNE_EXPRESSIONTEMPLATES
00044 
00045 namespace fvmeta
00046 {
00051     template<class K>
00052     inline typename FieldTraits<K>::real_type absreal (const K& k)
00053     {
00054         return std::abs(k);
00055     }
00056 
00061     template<class K>
00062     inline typename FieldTraits<K>::real_type absreal (const std::complex<K>& c)
00063     {
00064         return std::abs(c.real()) + std::abs(c.imag());
00065     }
00066 
00071     template<class K>
00072     inline typename FieldTraits<K>::real_type abs2 (const K& k)
00073     {
00074         return k*k;
00075     }
00076     
00081     template<class K>
00082     inline typename FieldTraits<K>::real_type abs2 (const std::complex<K>& c)
00083     {
00084         return c.real()*c.real() + c.imag()*c.imag();
00085     }
00086     
00091     template<class K, bool isInteger = std::numeric_limits<K>::is_integer>
00092     struct Sqrt
00093     {
00094         static inline typename FieldTraits<K>::real_type sqrt (const K& k)
00095         {
00096             return std::sqrt(k);
00097         }
00098     };
00099     
00104     template<class K>
00105     struct Sqrt<K, true>
00106     {
00107         static inline typename FieldTraits<K>::real_type sqrt (const K& k)
00108         {
00109             return typename FieldTraits<K>::real_type(std::sqrt(double(k)));
00110         }
00111     };
00112     
00117     template<class K>
00118     inline typename FieldTraits<K>::real_type sqrt (const K& k)
00119     {
00120         return Sqrt<K>::sqrt(k);
00121     }
00122 
00123 }
00124 
00125 #endif
00126 
00128   template<class C, class T>
00129   class FieldIterator : 
00130     public Dune::RandomAccessIteratorFacade<FieldIterator<C,T>,T, T&, int>
00131   {
00132     friend class FieldIterator<typename remove_const<C>::type, typename remove_const<T>::type >;
00133     friend class FieldIterator<const typename remove_const<C>::type, const typename remove_const<T>::type >;
00134     
00135   public:
00136     
00140     typedef std::ptrdiff_t DifferenceType;
00141     
00142     // Constructors needed by the base iterators.
00143     FieldIterator()
00144       : container_(0), position_(0)
00145       {}
00146     
00147     FieldIterator(C& cont, DifferenceType pos)
00148       : container_(&cont), position_(pos)
00149       {}
00150     
00151     FieldIterator(const FieldIterator<typename remove_const<C>::type, typename remove_const<T>::type >& other)
00152       : container_(other.container_), position_(other.position_)
00153       {}
00154     
00155 #if 0    
00156     FieldIterator(const FieldIterator<const typename remove_const<C>::type, const typename remove_const<T>::type >& other)
00157       : container_(other.container_), position_(other.position_)
00158       {}
00159 #endif
00160 #if 0
00161     FieldIterator(const FieldIterator<C,T>& other)
00162       : container_(other.container_), position_(other.position_)
00163       {}
00164 #endif
00165     // Methods needed by the forward iterator
00166     bool equals(const FieldIterator<typename remove_const<C>::type,typename remove_const<T>::type>& other) const
00167       {
00168         return position_ == other.position_ && container_ == other.container_;
00169       }
00170     
00171     
00172     bool equals(const FieldIterator<const typename remove_const<C>::type,const typename remove_const<T>::type>& other) const
00173       {
00174         return position_ == other.position_ && container_ == other.container_;
00175       }
00176     
00177     T& dereference() const{
00178       return container_->operator[](position_);
00179     }
00180     
00181     void increment(){
00182       ++position_;
00183     }
00184     
00185     // Additional function needed by BidirectionalIterator
00186     void decrement(){
00187       --position_;
00188     }
00189     
00190     // Additional function needed by RandomAccessIterator
00191     T& elementAt(DifferenceType i)const{
00192       return container_->operator[](position_+i);
00193     }
00194     
00195     void advance(DifferenceType n){
00196       position_=position_+n;
00197     }
00198     
00199     std::ptrdiff_t distanceTo(FieldIterator<const typename remove_const<C>::type,const typename remove_const<T>::type> other)const
00200       {
00201         assert(other.container_==container_);
00202         return other.position_ - position_;
00203       }
00204     
00205     std::ptrdiff_t distanceTo(FieldIterator<typename remove_const<C>::type, typename remove_const<T>::type> other)const
00206       {
00207         assert(other.container_==container_);
00208         return other.position_ - position_;
00209       }
00210     
00212     DifferenceType index () const
00213       {
00214         return this->position_;
00215       }
00216     
00217   private:
00218     C *container_;
00219     DifferenceType position_;
00220   };
00221 
00223   template<class T>
00224   struct IteratorType
00225   {
00226     typedef typename T::Iterator type;
00227   };
00228   
00229   template<class T>
00230   struct IteratorType<const T>
00231   {
00232     typedef typename T::ConstIterator type;
00233   };
00234   
00235 #ifdef DUNE_EXPRESSIONTEMPLATES
00236 
00237   template<class V>
00238   class FlatIterator :
00239     public ForwardIteratorFacade<FlatIterator<V>,
00240                                  typename Dune::FieldType<V>::type,
00241                                  typename Dune::FieldType<V>::type&,
00242                                  int>
00243   {
00244   public:
00245     typedef typename IteratorType<V>::type BlockIterator;
00246     typedef std::ptrdiff_t DifferenceType;
00247 //    typedef typename BlockIterator::DifferenceType DifferenceType;
00248     typedef typename BlockType<V>::type block_type;
00249     typedef typename FieldType<V>::type field_type;
00250     typedef FlatIterator<block_type> SubBlockIterator;
00251     FlatIterator(const BlockIterator & i) :
00252       it(i), bit(i->begin()), bend(i->end()) {};
00253     void increment ()
00254       {
00255         ++bit;
00256         if (bit == bend)
00257         {
00258           ++it;
00259           bit = it->begin();
00260           bend = it->end();
00261         }
00262       }
00263     bool equals (const FlatIterator & fit) const
00264       {
00265         return fit.it == it && fit.bit == bit;
00266       }
00267     field_type& dereference() const
00268       {
00269         return *bit;
00270       }
00272     DifferenceType index () const
00273       {
00274         return bit.index();
00275       }
00276   private:
00277     BlockIterator it;
00278     SubBlockIterator bit;
00279     SubBlockIterator bend;
00280   };
00281 
00284   template<class K, int N>
00285   class FlatIterator< FieldVector<K,N> > :
00286     public ForwardIteratorFacade<FlatIterator< FieldVector<K,N> >,
00287                                  K, K&, int>
00288   {
00289   public:
00290     typedef typename FieldVector<K,N>::Iterator BlockIterator;
00291     typedef std::ptrdiff_t DifferenceType;
00292 //    typedef typename BlockIterator::DifferenceType DifferenceType;
00293     typedef typename FieldVector<K,N>::field_type field_type;
00294     FlatIterator(const BlockIterator & i) : it(i) {};
00295     void increment ()
00296       {
00297         ++it;
00298       }
00299     bool equals (const FlatIterator & fit) const
00300       {
00301         return fit.it == it;
00302       }
00303     field_type& dereference() const
00304       {
00305         return *it;
00306       }
00308     DifferenceType index () const
00309       {
00310         return it.index();
00311       }
00312   private:
00313     BlockIterator it;
00314   };
00315 
00318   template<class K, int N>
00319   class FlatIterator< const FieldVector<K,N> > :
00320     public ForwardIteratorFacade<FlatIterator< const FieldVector<K,N> >,
00321                                  const K, const K&, int>
00322   {
00323   public:
00324     typedef typename FieldVector<K,N>::ConstIterator BlockIterator;
00325     typedef std::ptrdiff_t DifferenceType;
00326 //    typedef typename BlockIterator::DifferenceType DifferenceType;
00327     typedef typename FieldVector<K,N>::field_type field_type;
00328     FlatIterator(const BlockIterator & i) : it(i) {};
00329     void increment ()
00330       {
00331         ++it;
00332       }
00333     bool equals (const FlatIterator & fit) const
00334       {
00335         return fit.it == it;
00336       }
00337     const field_type& dereference() const
00338       {
00339         return *it;
00340       }
00342     DifferenceType index () const
00343       {
00344         return it.index();
00345       }
00346   private:
00347     BlockIterator it;
00348   };
00349 #endif
00350 
00351 
00360   template< class K, int SIZE >
00361   class FieldVector
00362 #ifdef DUNE_EXPRESSIONTEMPLATES
00363   : public Dune :: ExprTmpl :: Vector< FieldVector< K, SIZE > >
00364 #endif
00365   {
00366   public:
00367         // remember size of vector 
00368         enum { dimension = SIZE };
00369 
00370         // standard constructor and everything is sufficient ...
00371 
00372         //===== type definitions and constants
00373 
00375         typedef K field_type;
00376 
00378         typedef K block_type;
00379 
00381     typedef std::size_t size_type;
00382     
00384         enum {
00386           blocklevel = 1
00387         };
00388 
00390         enum {
00392           size = SIZE
00393         };
00394 
00396         FieldVector() {}
00397 
00398 #ifndef DUNE_EXPRESSIONTEMPLATES
00399 
00400         explicit FieldVector (const K& t)
00401         {
00402         for (size_type i=0; i<SIZE; i++) p[i] = t;
00403         }
00404 
00405         //===== assignment from scalar
00407         FieldVector& operator= (const K& k)
00408         {
00409         for (size_type i=0; i<SIZE; i++)
00410             p[i] = k;
00411         return *this;   
00412         }
00413 
00414 #else
00415 
00416         explicit FieldVector (const K& t)
00417         {
00418 #ifdef DUNE_VVERBOSE
00419       Dune::dvverb << INDENT << "FieldVector Copy Constructor Scalar\n";
00420 #endif
00421       assignFrom(t);
00422         }
00424         FieldVector& operator= (const K& k)
00425         {
00426 #ifdef DUNE_VVERBOSE
00427       Dune::dvverb << INDENT << "FieldVector Assignment Operator Scalar\n";
00428 #endif
00429       return assignFrom(k);
00430         }
00431     template <class E>
00432     FieldVector (Dune::ExprTmpl::Expression<E> op) {
00433 #ifdef DUNE_VVERBOSE
00434       Dune::dvverb << INDENT << "FieldVector Copy Constructor Expression\n";
00435 #endif
00436       assignFrom(op);
00437     }
00438     template <class V>
00439     FieldVector (const Dune::ExprTmpl::Vector<V> & op) {
00440 #ifdef DUNE_VVERBOSE
00441       Dune::dvverb << INDENT << "FieldVector Copy Operator Vector\n";
00442 #endif
00443       assignFrom(op);
00444     }
00445     template <class E>
00446     FieldVector&  operator = (Dune::ExprTmpl::Expression<E> op) {
00447 #ifdef DUNE_VVERBOSE
00448       Dune::dvverb << INDENT << "FieldVector Assignment Operator Expression\n";
00449 #endif
00450       return assignFrom(op);
00451     }
00452     template <class V>
00453     FieldVector& operator = (const Dune::ExprTmpl::Vector<V> & op) {
00454 #ifdef DUNE_VVERBOSE
00455       Dune::dvverb << INDENT << "FieldVector Assignment Operator Vector\n";
00456 #endif
00457       return assignFrom(op);
00458     }
00459 #endif
00460 
00462     FieldVector ( const FieldVector &other )
00463     {
00464       for( size_type i = 0; i < SIZE; ++i )
00465         p[ i ] = other[ i ];
00466     }
00467 
00469       FieldVector& operator= (const FieldVector& other) {
00470           for (size_type i=0; i<SIZE; i++)
00471               p[i] = other[i];
00472           return *this;
00473       }
00474       
00475 
00476         //===== access to components
00477 
00479         K& operator[] (size_type i)
00480         {
00481 #ifdef DUNE_ISTL_WITH_CHECKING
00482           if (i<0 || i>=SIZE) DUNE_THROW(MathError,"index out of range");
00483 #endif
00484           return p[i];
00485         }
00486 
00488         const K& operator[] (size_type i) const
00489         {
00490 #ifdef DUNE_ISTL_WITH_CHECKING
00491           if (i<0 || i>=SIZE) DUNE_THROW(MathError,"index out of range");
00492 #endif
00493           return p[i];
00494         }
00495 
00497     typedef FieldIterator<FieldVector<K,SIZE>,K> Iterator;
00499     typedef Iterator iterator;
00500     
00502         Iterator begin ()
00503         {
00504           return Iterator(*this,0);
00505         }
00506 
00508         Iterator end ()
00509         {
00510           return Iterator(*this,SIZE);
00511         }
00512 
00514         Iterator rbegin ()
00515         {
00516           return Iterator(*this,SIZE-1);
00517         }
00518 
00520         Iterator rend ()
00521         {
00522           return Iterator(*this,-1);
00523         }
00524 
00526         Iterator find (size_type i)
00527         {
00528           if (i<SIZE)
00529                 return Iterator(*this,i);
00530           else
00531                 return Iterator(*this,SIZE);
00532         }
00533 
00535     typedef FieldIterator<const FieldVector<K,SIZE>,const K> ConstIterator;
00537     typedef ConstIterator const_iterator;
00538 
00540         ConstIterator begin () const
00541         {
00542           return ConstIterator(*this,0);
00543         }
00544 
00546         ConstIterator end () const
00547         {
00548           return ConstIterator(*this,SIZE);
00549         }
00550 
00552         ConstIterator rbegin () const
00553         {
00554           return ConstIterator(*this,SIZE-1);
00555         }
00556 
00558         ConstIterator rend () const
00559         {
00560           return ConstIterator(*this,-1);
00561         }
00562 
00564         ConstIterator find (size_type i) const
00565         {
00566           if (i<SIZE)
00567                 return ConstIterator(*this,i);
00568           else
00569                 return ConstIterator(*this,SIZE);
00570         }
00571     
00572 #ifndef DUNE_EXPRESSIONTEMPLATES
00573         //===== vector space arithmetic
00574 
00576         FieldVector& operator+= (const FieldVector& y)
00577         {
00578             for (size_type i=0; i<SIZE; i++)
00579                 p[i] += y.p[i];
00580           return *this;
00581         }
00582 
00584         FieldVector& operator-= (const FieldVector& y)
00585         {
00586             for (size_type i=0; i<SIZE; i++)
00587                 p[i] -= y.p[i];
00588           return *this;
00589         }
00590 
00592         FieldVector<K, size> operator+ (const FieldVector<K, size>& b) const
00593         {
00594           FieldVector<K, size> z = *this;
00595           return (z+=b);
00596         }
00597 
00599         FieldVector<K, size> operator- (const FieldVector<K, size>& b) const
00600         {
00601           FieldVector<K, size> z = *this;
00602           return (z-=b);
00603         }
00604 
00606         FieldVector& operator+= (const K& k)
00607         {
00608             for (size_type i=0; i<SIZE; i++)
00609                 p[i] += k;
00610           return *this;
00611         }
00612 
00614         FieldVector& operator-= (const K& k)
00615         {
00616             for (size_type i=0; i<SIZE; i++)
00617                 p[i] -= k;
00618           return *this;
00619         }
00620 
00622         FieldVector& operator*= (const K& k)
00623         {
00624             for (size_type i=0; i<SIZE; i++)
00625                 p[i] *= k;
00626           return *this;
00627         }
00628 
00630         FieldVector& operator/= (const K& k)
00631         {
00632           for (size_type i=0; i<SIZE; i++)
00633               p[i] /= k;
00634           return *this;
00635         }
00636 
00637 #endif
00638 
00640         bool operator== (const FieldVector& y) const
00641         {
00642             for (size_type i=0; i<SIZE; i++)
00643                 if (p[i]!=y.p[i])
00644                     return false;
00645 
00646             return true;
00647         }
00648 
00650         bool operator!= (const FieldVector& y) const
00651         {
00652             return !operator==(y);
00653         }
00654     
00655 
00657         FieldVector& axpy (const K& a, const FieldVector& y)
00658         {
00659 #ifndef DUNE_EXPRESSIONTEMPLATES
00660             for (size_type i=0; i<SIZE; i++)
00661                 p[i] += a*y.p[i];
00662 #else
00663       *this += a*y;
00664 #endif
00665           return *this;
00666         }
00667 
00668 #ifndef DUNE_EXPRESSIONTEMPLATES
00669         //===== Euclidean scalar product
00670 
00672     K operator* (const FieldVector& y) const
00673         {
00674             K result = 0;
00675             for (int i=0; i<size; i++)
00676                 result += p[i]*y[i];
00677             return result;
00678         }
00679 
00680         //===== norms
00681 
00683     typename FieldTraits<K>::real_type one_norm() const {
00684         typename FieldTraits<K>::real_type result = 0;
00685         for (int i=0; i<size; i++)
00686             result += std::abs(p[i]);
00687         return result;
00688     }
00689 
00690 
00692     typename FieldTraits<K>::real_type one_norm_real () const
00693         {
00694         typename FieldTraits<K>::real_type result = 0;
00695         for (int i=0; i<size; i++)
00696             result += fvmeta::absreal(p[i]);
00697         return result;
00698         }
00699 
00701     typename FieldTraits<K>::real_type two_norm () const
00702         {
00703         typename FieldTraits<K>::real_type result = 0;
00704         for (int i=0; i<size; i++)
00705             result += fvmeta::abs2(p[i]);
00706         return fvmeta::sqrt(result);
00707         }
00708 
00710     typename FieldTraits<K>::real_type two_norm2 () const
00711         {
00712         typename FieldTraits<K>::real_type result = 0;
00713         for (int i=0; i<size; i++)
00714             result += fvmeta::abs2(p[i]);
00715         return result;
00716         }
00717 
00719     typename FieldTraits<K>::real_type infinity_norm () const
00720         {
00721         typename FieldTraits<K>::real_type result = 0;
00722         for (int i=0; i<size; i++)
00723             result = std::max(result, std::abs(p[i]));
00724         return result;
00725         }
00726 
00728         typename FieldTraits<K>::real_type infinity_norm_real () const
00729         {
00730         typename FieldTraits<K>::real_type result = 0;
00731         for (int i=0; i<size; i++)
00732             result = std::max(result, fvmeta::absreal(p[i]));
00733         return result;
00734         }
00735 #endif
00736 
00737         //===== sizes
00738 
00740         size_type N () const
00741         {
00742           return SIZE;
00743         }
00744 
00746         size_type dim () const
00747         {
00748           return SIZE;
00749         }
00750 
00751   private:
00752         // the data, very simply a built in array
00753       K p[(SIZE > 0) ? SIZE : 1]; 
00754   };
00755 
00764   template<typename K, int n>
00765   std::ostream& operator<< (std::ostream& s, const FieldVector<K,n>& v)
00766   {
00767       for (typename FieldVector<K,n>::size_type i=0; i<n; i++)
00768               s << ((i>0) ? " " : "") << v[i];
00769       return s;
00770   }
00771 
00783   template< class K, int SIZE >
00784   inline std :: istream &operator>> ( std :: istream &in,
00785                                       FieldVector< K, SIZE > &v )
00786   {
00787     FieldVector< K, SIZE > w;
00788     for( typename FieldVector< K, SIZE > :: size_type i = 0; i < SIZE; ++i )
00789       in >> w[ i ];
00790     if( in )
00791       v = w;
00792     return in;
00793   }
00794 
00795   
00796   // forward declarations
00797   template<class K, int n, int m> class FieldMatrix;
00798 
00799 
00800 #ifndef DOXYGEN
00801 
00803   template< class K >
00804   class FieldVector< K, 1 >
00805 #ifdef DUNE_EXPRESSIONTEMPLATES
00806   : public Dune :: ExprTmpl :: Vector< FieldVector< K, 1 > >
00807 #endif
00808   {
00809     enum { n = 1 };
00810   public:
00811         friend class FieldMatrix<K,1,1>;
00812 
00813         //===== type definitions and constants
00814 
00816         typedef K field_type;
00817 
00819         typedef K block_type;
00820 
00822     typedef std::size_t size_type;
00823     
00825         enum {blocklevel = 1};
00826 
00828         enum {size = 1};
00829 
00831         enum {dimension = 1};
00832 
00833         //===== construction
00834 
00836         FieldVector ()
00837         {       }
00838 
00840         FieldVector (const K& k)
00841         {
00842           p = k;
00843         }
00844 
00846         FieldVector& operator= (const K& k)
00847         {
00848           p = k;
00849           return *this;   
00850         }
00851 
00852 #ifdef DUNE_EXPRESSIONTEMPLATES
00853     template <class E>
00854     FieldVector (Dune::ExprTmpl::Expression<E> op) {
00855 #ifdef DUNE_VVERBOSE
00856       Dune::dvverb << INDENT << "FieldVector<1> Copy Constructor Expression\n";
00857 #endif
00858       assignFrom(op);
00859     }
00860     template <class V>
00861     FieldVector (const Dune::ExprTmpl::Vector<V> & op) {
00862 #ifdef DUNE_VVERBOSE
00863       Dune::dvverb << INDENT << "FieldVector<1> Copy Operator Vector\n";
00864 #endif
00865       assignFrom(op);
00866     }
00867     template <class E>
00868     FieldVector&  operator = (Dune::ExprTmpl::Expression<E> op) {
00869 #ifdef DUNE_VVERBOSE
00870       Dune::dvverb << INDENT << "FieldVector<1> Assignment Operator Expression\n";
00871 #endif
00872       return assignFrom(op);
00873     }
00874     template <class V>
00875     FieldVector& operator = (const Dune::ExprTmpl::Vector<V> & op) {
00876 #ifdef DUNE_VVERBOSE
00877       Dune::dvverb << INDENT << "FieldVector<1> Assignment Operator Vector\n";
00878 #endif
00879       return assignFrom(op);
00880     }
00881 #endif
00882 
00883         //===== access to components
00884 
00886         K& operator[] (size_type i)
00887         {
00888 #ifdef DUNE_ISTL_WITH_CHECKING
00889           if (i != 0) DUNE_THROW(MathError,"index out of range");
00890 #endif
00891           return p;
00892         }
00893 
00895         const K& operator[] (size_type i) const
00896         {
00897 #ifdef DUNE_ISTL_WITH_CHECKING
00898           if (i != 0) DUNE_THROW(MathError,"index out of range");
00899 #endif
00900           return p;
00901         }
00902 
00904     typedef FieldIterator<FieldVector<K,n>,K> Iterator;
00906     typedef Iterator iterator;
00907     
00909         Iterator begin ()
00910         {
00911           return Iterator(*this,0);
00912         }
00913 
00915         Iterator end ()
00916         {
00917           return Iterator(*this,n);
00918         }
00919 
00921         Iterator rbegin ()
00922         {
00923           return Iterator(*this,n-1);
00924         }
00925 
00927         Iterator rend ()
00928         {
00929           return Iterator(*this,-1);
00930         }
00931 
00933         Iterator find (size_type i)
00934         {
00935           if (i<n)
00936                 return Iterator(*this,i);
00937           else
00938                 return Iterator(*this,n);
00939         }
00940 
00942     typedef FieldIterator<const FieldVector<K,n>,const K> ConstIterator;
00944     typedef ConstIterator const_iterator;
00945 
00947         ConstIterator begin () const
00948         {
00949           return ConstIterator(*this,0);
00950         }
00951 
00953         ConstIterator end () const
00954         {
00955           return ConstIterator(*this,n);
00956         }
00957 
00959         ConstIterator rbegin () const
00960         {
00961           return ConstIterator(*this,n-1);
00962         }
00963 
00965         ConstIterator rend () const
00966         {
00967           return ConstIterator(*this,-1);
00968         }
00969 
00971         ConstIterator find (size_type i) const
00972         {
00973           if (i<n)
00974                 return ConstIterator(*this,i);
00975           else
00976                 return ConstIterator(*this,n);
00977         }
00978         //===== vector space arithmetic
00979 
00981         FieldVector& operator+= (const K& k)
00982         {
00983           p += k;
00984           return *this;
00985         }
00986 
00988         FieldVector& operator-= (const K& k)
00989         {
00990           p -= k;
00991           return *this;
00992         }
00993 
00995         FieldVector& operator*= (const K& k)
00996         {
00997           p *= k;
00998           return *this;
00999         }
01000 
01002         FieldVector& operator/= (const K& k)
01003         {
01004           p /= k;
01005           return *this;
01006         }
01007 
01009         FieldVector& axpy (const K& a, const FieldVector& y)
01010         {
01011           p += a*y.p;
01012           return *this;
01013         }
01014         
01016         inline K operator* ( const K & k ) const
01017         {
01018           return p * k;
01019         }
01020 
01021         //===== norms
01022 
01024     typename FieldTraits<K>::real_type one_norm () const
01025         {
01026         return std::abs(p);
01027         }
01028 
01030     typename FieldTraits<K>::real_type one_norm_real () const
01031         {
01032         return fvmeta::absreal(p);
01033         }
01034 
01036     typename FieldTraits<K>::real_type two_norm () const
01037         {
01038         return fvmeta::sqrt(fvmeta::abs2(p));
01039         }
01040 
01042     typename FieldTraits<K>::real_type two_norm2 () const
01043         {
01044         return fvmeta::abs2(p);
01045         }
01046 
01048     typename FieldTraits<K>::field_type infinity_norm () const
01049         {
01050         return std::abs(p);
01051         }
01052 
01054         typename FieldTraits<K>::real_type infinity_norm_real () const
01055         {
01056         return fvmeta::absreal(p);
01057         }
01058 
01059         //===== sizes
01060 
01062         size_type N () const
01063         {
01064           return 1;
01065         }
01066 
01068         size_type dim () const
01069         {
01070           return 1;
01071         }
01072 
01073         //===== conversion operator
01074 
01076         operator K () {return p;}
01077 
01079         operator K () const {return p;}
01080 
01081   private:
01082         // the data
01083         K p; 
01084   };
01085 
01086 #ifndef DUNE_EXPRESSIONTEMPLATES
01087 
01088   template<class K>
01089   inline FieldVector<K,1> operator+ (const FieldVector<K,1>& a, const FieldVector<K,1>& b)
01090   {
01091     FieldVector<K,1> z = a;
01092     return (z+=b);
01093   }
01094   
01096   template<class K>
01097   inline FieldVector<K,1> operator- (const FieldVector<K,1>& a, const FieldVector<K,1>& b)
01098   {
01099     FieldVector<K,1> z = a;
01100     return (z-=b);
01101   }
01102   
01104   template<class K>
01105   inline FieldVector<K,1> operator+ (const FieldVector<K,1>& a, const K b)
01106   {
01107     FieldVector<K,1> z = a;
01108     return (z[0]+=b);
01109   }
01110   
01112   template<class K>
01113   inline FieldVector<K,1> operator- (const FieldVector<K,1>& a, const K b)
01114   {
01115     FieldVector<K,1> z = a;
01116     return (z[0]-=b);
01117   }
01118 
01120   template<class K>
01121   inline FieldVector<K,1> operator+ (const K a, const FieldVector<K,1>& b)
01122   {
01123     FieldVector<K,1> z = a;
01124     return (z[0]+=b);
01125   }
01126   
01128   template<class K>
01129   inline FieldVector<K,1> operator- (const K a, const FieldVector<K,1>& b)
01130   {
01131     FieldVector<K,1> z = a;
01132     return (z[0]-=b);
01133   }
01134 #endif
01135 #endif
01136 
01139 } // end namespace
01140 
01141 #endif
Generated on Mon Apr 26 10:45:21 2010 for dune-common by  doxygen 1.6.3