00001
00002 #ifndef DUNE_FVECTOR_HH
00003 #define DUNE_FVECTOR_HH
00004
00005 #include<cmath>
00006 #include<cstddef>
00007 #include<cstdlib>
00008 #include<complex>
00009
00010 #include "exceptions.hh"
00011 #include "genericiterator.hh"
00012
00013 #ifdef DUNE_EXPRESSIONTEMPLATES
00014 #include "exprtmpl.hh"
00015 #endif
00016
00017 namespace Dune {
00018
00019 #ifndef DUNE_EXPRESSIONTEMPLATES
00020
00031
00032 template<class K, int SIZE> class FieldVector;
00033
00034 #endif
00035
00036 #ifndef DUNE_EXPRESSIONTEMPLATES
00037
00038 template<class K>
00039 inline double fvmeta_absreal (const K& k)
00040 {
00041 return std::abs(k);
00042 }
00043
00044 template<class K>
00045 inline double fvmeta_absreal (const std::complex<K>& c)
00046 {
00047 return fvmeta_abs(c.real()) + fvmeta_abs(c.imag());
00048 }
00049
00050 template<class K>
00051 inline double fvmeta_abs2 (const K& k)
00052 {
00053 return k*k;
00054 }
00055
00056 template<class K>
00057 inline double fvmeta_abs2 (const std::complex<K>& c)
00058 {
00059 return c.real()*c.real() + c.imag()*c.imag();
00060 }
00061
00062 #endif
00063
00065 template<class C, class T>
00066 class FieldIterator :
00067 public Dune::RandomAccessIteratorFacade<FieldIterator<C,T>,T, T&, int>
00068 {
00069 friend class FieldIterator<typename remove_const<C>::type, typename remove_const<T>::type >;
00070 friend class FieldIterator<const typename remove_const<C>::type, const typename remove_const<T>::type >;
00071
00072 public:
00073
00077 typedef std::ptrdiff_t DifferenceType;
00078
00079
00080 FieldIterator()
00081 : container_(0), position_(0)
00082 {}
00083
00084 FieldIterator(C& cont, DifferenceType pos)
00085 : container_(&cont), position_(pos)
00086 {}
00087
00088 FieldIterator(const FieldIterator<typename remove_const<C>::type, typename remove_const<T>::type >& other)
00089 : container_(other.container_), position_(other.position_)
00090 {}
00091
00092 #if 0
00093 FieldIterator(const FieldIterator<const typename remove_const<C>::type, const typename remove_const<T>::type >& other)
00094 : container_(other.container_), position_(other.position_)
00095 {}
00096 #endif
00097 #if 0
00098 FieldIterator(const FieldIterator<C,T>& other)
00099 : container_(other.container_), position_(other.position_)
00100 {}
00101 #endif
00102
00103 bool equals(const FieldIterator<typename remove_const<C>::type,typename remove_const<T>::type>& other) const
00104 {
00105 return position_ == other.position_ && container_ == other.container_;
00106 }
00107
00108
00109 bool equals(const FieldIterator<const typename remove_const<C>::type,const typename remove_const<T>::type>& other) const
00110 {
00111 return position_ == other.position_ && container_ == other.container_;
00112 }
00113
00114 T& dereference() const{
00115 return container_->operator[](position_);
00116 }
00117
00118 void increment(){
00119 ++position_;
00120 }
00121
00122
00123 void decrement(){
00124 --position_;
00125 }
00126
00127
00128 T& elementAt(DifferenceType i)const{
00129 return container_->operator[](position_+i);
00130 }
00131
00132 void advance(DifferenceType n){
00133 position_=position_+n;
00134 }
00135
00136 std::ptrdiff_t distanceTo(FieldIterator<const typename remove_const<C>::type,const typename remove_const<T>::type> other)const
00137 {
00138 assert(other.container_==container_);
00139 return other.position_ - position_;
00140 }
00141
00142 std::ptrdiff_t distanceTo(FieldIterator<typename remove_const<C>::type, typename remove_const<T>::type> other)const
00143 {
00144 assert(other.container_==container_);
00145 return other.position_ - position_;
00146 }
00147
00149 DifferenceType index () const
00150 {
00151 return this->position_;
00152 }
00153
00154 private:
00155 C *container_;
00156 DifferenceType position_;
00157 };
00158
00160 template<class T>
00161 struct IteratorType
00162 {
00163 typedef typename T::Iterator type;
00164 };
00165
00166 template<class T>
00167 struct IteratorType<const T>
00168 {
00169 typedef typename T::ConstIterator type;
00170 };
00171
00172 #ifdef DUNE_EXPRESSIONTEMPLATES
00174 template<class V>
00175 class FlatIterator :
00176 public ForwardIteratorFacade<FlatIterator<V>,
00177 typename Dune::FieldType<V>::type,
00178 typename Dune::FieldType<V>::type&,
00179 int>
00180 {
00181 public:
00182 typedef typename IteratorType<V>::type BlockIterator;
00183 typedef std::ptrdiff_t DifferenceType;
00184
00185 typedef typename BlockType<V>::type block_type;
00186 typedef typename FieldType<V>::type field_type;
00187 typedef FlatIterator<block_type> SubBlockIterator;
00188 FlatIterator(const BlockIterator & i) :
00189 it(i), bit(i->begin()), bend(i->end()) {};
00190 void increment ()
00191 {
00192 ++bit;
00193 if (bit == bend)
00194 {
00195 ++it;
00196 bit = it->begin();
00197 bend = it->end();
00198 }
00199 }
00200 bool equals (const FlatIterator & fit) const
00201 {
00202 return fit.it == it && fit.bit == bit;
00203 }
00204 field_type& dereference() const
00205 {
00206 return *bit;
00207 }
00209 DifferenceType index () const
00210 {
00211 return bit.index();
00212 }
00213 private:
00214 BlockIterator it;
00215 SubBlockIterator bit;
00216 SubBlockIterator bend;
00217 };
00218
00221 template<class K, int N>
00222 class FlatIterator< FieldVector<K,N> > :
00223 public ForwardIteratorFacade<FlatIterator< FieldVector<K,N> >,
00224 K, K&, int>
00225 {
00226 public:
00227 typedef typename FieldVector<K,N>::Iterator BlockIterator;
00228 typedef std::ptrdiff_t DifferenceType;
00229
00230 typedef typename FieldVector<K,N>::field_type field_type;
00231 FlatIterator(const BlockIterator & i) : it(i) {};
00232 void increment ()
00233 {
00234 ++it;
00235 }
00236 bool equals (const FlatIterator & fit) const
00237 {
00238 return fit.it == it;
00239 }
00240 field_type& dereference() const
00241 {
00242 return *it;
00243 }
00245 DifferenceType index () const
00246 {
00247 return it.index();
00248 }
00249 private:
00250 BlockIterator it;
00251 };
00252
00255 template<class K, int N>
00256 class FlatIterator< const FieldVector<K,N> > :
00257 public ForwardIteratorFacade<FlatIterator< const FieldVector<K,N> >,
00258 const K, const K&, int>
00259 {
00260 public:
00261 typedef typename FieldVector<K,N>::ConstIterator BlockIterator;
00262 typedef std::ptrdiff_t DifferenceType;
00263
00264 typedef typename FieldVector<K,N>::field_type field_type;
00265 FlatIterator(const BlockIterator & i) : it(i) {};
00266 void increment ()
00267 {
00268 ++it;
00269 }
00270 bool equals (const FlatIterator & fit) const
00271 {
00272 return fit.it == it;
00273 }
00274 const field_type& dereference() const
00275 {
00276 return *it;
00277 }
00279 DifferenceType index () const
00280 {
00281 return it.index();
00282 }
00283 private:
00284 BlockIterator it;
00285 };
00286 #endif
00287
00288
00297 template< class K, int SIZE >
00298 class FieldVector
00299 #ifdef DUNE_EXPRESSIONTEMPLATES
00300 : public Dune :: ExprTmpl :: Vector< FieldVector< K, SIZE > >
00301 #endif
00302 {
00303 public:
00304
00305 enum { dimension = SIZE };
00306
00307
00308
00309
00310
00312 typedef K field_type;
00313
00315 typedef K block_type;
00316
00318 typedef std::size_t size_type;
00319
00321 enum {
00323 blocklevel = 1
00324 };
00325
00327 enum {
00329 size = SIZE
00330 };
00331
00333 FieldVector() {}
00334
00335 #ifndef DUNE_EXPRESSIONTEMPLATES
00337 explicit FieldVector (const K& t)
00338 {
00339 for (size_type i=0; i<SIZE; i++) p[i] = t;
00340 }
00341
00342
00344 FieldVector& operator= (const K& k)
00345 {
00346
00347 for (size_type i=0; i<SIZE; i++)
00348 p[i] = k;
00349 return *this;
00350 }
00351
00352 #else
00354 explicit FieldVector (const K& t)
00355 {
00356 #ifdef DUNE_VVERBOSE
00357 Dune::dvverb << INDENT << "FieldVector Copy Constructor Scalar\n";
00358 #endif
00359 assignFrom(t);
00360 }
00362 FieldVector& operator= (const K& k)
00363 {
00364 #ifdef DUNE_VVERBOSE
00365 Dune::dvverb << INDENT << "FieldVector Assignment Operator Scalar\n";
00366 #endif
00367 return assignFrom(k);
00368 }
00369 template <class E>
00370 FieldVector (Dune::ExprTmpl::Expression<E> op) {
00371 #ifdef DUNE_VVERBOSE
00372 Dune::dvverb << INDENT << "FieldVector Copy Constructor Expression\n";
00373 #endif
00374 assignFrom(op);
00375 }
00376 template <class V>
00377 FieldVector (const Dune::ExprTmpl::Vector<V> & op) {
00378 #ifdef DUNE_VVERBOSE
00379 Dune::dvverb << INDENT << "FieldVector Copy Operator Vector\n";
00380 #endif
00381 assignFrom(op);
00382 }
00383 template <class E>
00384 FieldVector& operator = (Dune::ExprTmpl::Expression<E> op) {
00385 #ifdef DUNE_VVERBOSE
00386 Dune::dvverb << INDENT << "FieldVector Assignment Operator Expression\n";
00387 #endif
00388 return assignFrom(op);
00389 }
00390 template <class V>
00391 FieldVector& operator = (const Dune::ExprTmpl::Vector<V> & op) {
00392 #ifdef DUNE_VVERBOSE
00393 Dune::dvverb << INDENT << "FieldVector Assignment Operator Vector\n";
00394 #endif
00395 return assignFrom(op);
00396 }
00397 #endif
00398
00400 FieldVector ( const FieldVector &other )
00401 {
00402 for( size_type i = 0; i < SIZE; ++i )
00403 p[ i ] = other[ i ];
00404 }
00405
00407 FieldVector& operator= (const FieldVector& other) {
00408 for (size_type i=0; i<SIZE; i++)
00409 p[i] = other[i];
00410 return *this;
00411 }
00412
00413
00414
00415
00417 K& operator[] (size_type i)
00418 {
00419 #ifdef DUNE_ISTL_WITH_CHECKING
00420 if (i<0 || i>=SIZE) DUNE_THROW(MathError,"index out of range");
00421 #endif
00422 return p[i];
00423 }
00424
00426 const K& operator[] (size_type i) const
00427 {
00428 #ifdef DUNE_ISTL_WITH_CHECKING
00429 if (i<0 || i>=SIZE) DUNE_THROW(MathError,"index out of range");
00430 #endif
00431 return p[i];
00432 }
00433
00435 typedef FieldIterator<FieldVector<K,SIZE>,K> Iterator;
00437 typedef Iterator iterator;
00438
00440 Iterator begin ()
00441 {
00442 return Iterator(*this,0);
00443 }
00444
00446 Iterator end ()
00447 {
00448 return Iterator(*this,SIZE);
00449 }
00450
00452 Iterator rbegin ()
00453 {
00454 return Iterator(*this,SIZE-1);
00455 }
00456
00458 Iterator rend ()
00459 {
00460 return Iterator(*this,-1);
00461 }
00462
00464 Iterator find (size_type i)
00465 {
00466 if (i<SIZE)
00467 return Iterator(*this,i);
00468 else
00469 return Iterator(*this,SIZE);
00470 }
00471
00473 typedef FieldIterator<const FieldVector<K,SIZE>,const K> ConstIterator;
00475 typedef ConstIterator const_iterator;
00476
00478 ConstIterator begin () const
00479 {
00480 return ConstIterator(*this,0);
00481 }
00482
00484 ConstIterator end () const
00485 {
00486 return ConstIterator(*this,SIZE);
00487 }
00488
00490 ConstIterator rbegin () const
00491 {
00492 return ConstIterator(*this,SIZE-1);
00493 }
00494
00496 ConstIterator rend () const
00497 {
00498 return ConstIterator(*this,-1);
00499 }
00500
00502 ConstIterator find (size_type i) const
00503 {
00504 if (i<SIZE)
00505 return ConstIterator(*this,i);
00506 else
00507 return ConstIterator(*this,SIZE);
00508 }
00509
00510 #ifndef DUNE_EXPRESSIONTEMPLATES
00511
00512
00514 FieldVector& operator+= (const FieldVector& y)
00515 {
00516 for (size_type i=0; i<SIZE; i++)
00517 p[i] += y.p[i];
00518 return *this;
00519 }
00520
00522 FieldVector& operator-= (const FieldVector& y)
00523 {
00524 for (size_type i=0; i<SIZE; i++)
00525 p[i] -= y.p[i];
00526 return *this;
00527 }
00528
00530 FieldVector<K, size> operator+ (const FieldVector<K, size>& b) const
00531 {
00532 FieldVector<K, size> z = *this;
00533 return (z+=b);
00534 }
00535
00537 FieldVector<K, size> operator- (const FieldVector<K, size>& b) const
00538 {
00539 FieldVector<K, size> z = *this;
00540 return (z-=b);
00541 }
00542
00544 FieldVector& operator+= (const K& k)
00545 {
00546 for (size_type i=0; i<SIZE; i++)
00547 p[i] += k;
00548 return *this;
00549 }
00550
00552 FieldVector& operator-= (const K& k)
00553 {
00554 for (size_type i=0; i<SIZE; i++)
00555 p[i] -= k;
00556 return *this;
00557 }
00558
00560 FieldVector& operator*= (const K& k)
00561 {
00562 for (size_type i=0; i<SIZE; i++)
00563 p[i] *= k;
00564 return *this;
00565 }
00566
00568 FieldVector& operator/= (const K& k)
00569 {
00570 for (size_type i=0; i<SIZE; i++)
00571 p[i] /= k;
00572 return *this;
00573 }
00574
00575 #endif
00576
00578 bool operator== (const FieldVector& y) const
00579 {
00580 for (size_type i=0; i<SIZE; i++)
00581 if (p[i]!=y.p[i])
00582 return false;
00583
00584 return true;
00585 }
00586
00588 bool operator!= (const FieldVector& y) const
00589 {
00590 return !operator==(y);
00591 }
00592
00593
00595 FieldVector& axpy (const K& a, const FieldVector& y)
00596 {
00597 #ifndef DUNE_EXPRESSIONTEMPLATES
00598 for (size_type i=0; i<SIZE; i++)
00599 p[i] += a*y.p[i];
00600 #else
00601 *this += a*y;
00602 #endif
00603 return *this;
00604 }
00605
00606 #ifndef DUNE_EXPRESSIONTEMPLATES
00607
00608
00610 K operator* (const FieldVector& y) const
00611 {
00612 K result = 0;
00613 for (int i=0; i<size; i++)
00614 result += p[i]*y[i];
00615 return result;
00616 }
00617
00618
00619
00621 double one_norm() const {
00622 double result = 0;
00623 for (int i=0; i<size; i++)
00624 result += std::abs(p[i]);
00625 return result;
00626 }
00627
00628
00630 double one_norm_real () const
00631 {
00632 double result = 0;
00633 for (int i=0; i<size; i++)
00634 result += fvmeta_absreal(p[i]);
00635 return result;
00636 }
00637
00639 double two_norm () const
00640 {
00641 double result = 0;
00642 for (int i=0; i<size; i++)
00643 result += fvmeta_abs2(p[i]);
00644 return std::sqrt(result);
00645 }
00646
00648 double two_norm2 () const
00649 {
00650 double result = 0;
00651 for (int i=0; i<size; i++)
00652 result += fvmeta_abs2(p[i]);
00653 return result;
00654 }
00655
00657 double infinity_norm () const
00658 {
00659 double result = 0;
00660 for (int i=0; i<size; i++)
00661 result = std::max(result, std::abs(p[i]));
00662 return result;
00663 }
00664
00666 double infinity_norm_real () const
00667 {
00668 double result = 0;
00669 for (int i=0; i<size; i++)
00670 result = std::max(result, fvmeta_absreal(p[i]));
00671 return result;
00672 }
00673 #endif
00674
00675
00676
00678 size_type N () const
00679 {
00680 return SIZE;
00681 }
00682
00684 size_type dim () const
00685 {
00686 return SIZE;
00687 }
00688
00689 private:
00690
00691 K p[(SIZE > 0) ? SIZE : 1];
00692 };
00693
00702 template<typename K, int n>
00703 std::ostream& operator<< (std::ostream& s, const FieldVector<K,n>& v)
00704 {
00705 for (typename FieldVector<K,n>::size_type i=0; i<n; i++)
00706 s << ((i>0) ? " " : "") << v[i];
00707 return s;
00708 }
00709
00721 template< class K, int SIZE >
00722 inline std :: istream &operator>> ( std :: istream &in,
00723 FieldVector< K, SIZE > &v )
00724 {
00725 FieldVector< K, SIZE > w;
00726 for( typename FieldVector< K, SIZE > :: size_type i = 0; i < SIZE; ++i )
00727 in >> w[ i ];
00728 if( in )
00729 v = w;
00730 return in;
00731 }
00732
00733
00734
00735 template<class K, int n, int m> class FieldMatrix;
00736
00737
00738
00741 template< class K >
00742 class FieldVector< K, 1 >
00743 #ifdef DUNE_EXPRESSIONTEMPLATES
00744 : public Dune :: ExprTmpl :: Vector< FieldVector< K, 1 > >
00745 #endif
00746 {
00747 enum { n = 1 };
00748 public:
00749 friend class FieldMatrix<K,1,1>;
00750
00751
00752
00754 typedef K field_type;
00755
00757 typedef K block_type;
00758
00760 typedef std::size_t size_type;
00761
00763 enum {blocklevel = 1};
00764
00766 enum {size = 1};
00767
00769 enum {dimension = 1};
00770
00771
00772
00774 FieldVector ()
00775 { }
00776
00778 FieldVector (const K& k)
00779 {
00780 p = k;
00781 }
00782
00784 FieldVector& operator= (const K& k)
00785 {
00786 p = k;
00787 return *this;
00788 }
00789
00790 #ifdef DUNE_EXPRESSIONTEMPLATES
00791 template <class E>
00792 FieldVector (Dune::ExprTmpl::Expression<E> op) {
00793 #ifdef DUNE_VVERBOSE
00794 Dune::dvverb << INDENT << "FieldVector<1> Copy Constructor Expression\n";
00795 #endif
00796 assignFrom(op);
00797 }
00798 template <class V>
00799 FieldVector (const Dune::ExprTmpl::Vector<V> & op) {
00800 #ifdef DUNE_VVERBOSE
00801 Dune::dvverb << INDENT << "FieldVector<1> Copy Operator Vector\n";
00802 #endif
00803 assignFrom(op);
00804 }
00805 template <class E>
00806 FieldVector& operator = (Dune::ExprTmpl::Expression<E> op) {
00807 #ifdef DUNE_VVERBOSE
00808 Dune::dvverb << INDENT << "FieldVector<1> Assignment Operator Expression\n";
00809 #endif
00810 return assignFrom(op);
00811 }
00812 template <class V>
00813 FieldVector& operator = (const Dune::ExprTmpl::Vector<V> & op) {
00814 #ifdef DUNE_VVERBOSE
00815 Dune::dvverb << INDENT << "FieldVector<1> Assignment Operator Vector\n";
00816 #endif
00817 return assignFrom(op);
00818 }
00819 #endif
00820
00821
00822
00824 K& operator[] (size_type i)
00825 {
00826 #ifdef DUNE_ISTL_WITH_CHECKING
00827 if (i != 0) DUNE_THROW(MathError,"index out of range");
00828 #endif
00829 return p;
00830 }
00831
00833 const K& operator[] (size_type i) const
00834 {
00835 #ifdef DUNE_ISTL_WITH_CHECKING
00836 if (i != 0) DUNE_THROW(MathError,"index out of range");
00837 #endif
00838 return p;
00839 }
00840
00842 typedef FieldIterator<FieldVector<K,n>,K> Iterator;
00844 typedef Iterator iterator;
00845
00847 Iterator begin ()
00848 {
00849 return Iterator(*this,0);
00850 }
00851
00853 Iterator end ()
00854 {
00855 return Iterator(*this,n);
00856 }
00857
00859 Iterator rbegin ()
00860 {
00861 return Iterator(*this,n-1);
00862 }
00863
00865 Iterator rend ()
00866 {
00867 return Iterator(*this,-1);
00868 }
00869
00871 Iterator find (size_type i)
00872 {
00873 if (i<n)
00874 return Iterator(*this,i);
00875 else
00876 return Iterator(*this,n);
00877 }
00878
00880 typedef FieldIterator<const FieldVector<K,n>,const K> ConstIterator;
00882 typedef ConstIterator const_iterator;
00883
00885 ConstIterator begin () const
00886 {
00887 return ConstIterator(*this,0);
00888 }
00889
00891 ConstIterator end () const
00892 {
00893 return ConstIterator(*this,n);
00894 }
00895
00897 ConstIterator rbegin () const
00898 {
00899 return ConstIterator(*this,n-1);
00900 }
00901
00903 ConstIterator rend () const
00904 {
00905 return ConstIterator(*this,-1);
00906 }
00907
00909 ConstIterator find (size_type i) const
00910 {
00911 if (i<n)
00912 return ConstIterator(*this,i);
00913 else
00914 return ConstIterator(*this,n);
00915 }
00916
00917
00919 FieldVector& operator+= (const K& k)
00920 {
00921 p += k;
00922 return *this;
00923 }
00924
00926 FieldVector& operator-= (const K& k)
00927 {
00928 p -= k;
00929 return *this;
00930 }
00931
00933 FieldVector& operator*= (const K& k)
00934 {
00935 p *= k;
00936 return *this;
00937 }
00938
00940 FieldVector& operator/= (const K& k)
00941 {
00942 p /= k;
00943 return *this;
00944 }
00945
00947 FieldVector& axpy (const K& a, const FieldVector& y)
00948 {
00949 p += a*y.p;
00950 return *this;
00951 }
00952
00954 inline K operator* ( const K & k ) const
00955 {
00956 return p * k;
00957 }
00958
00959
00960
00962 double one_norm () const
00963 {
00964 return std::abs(p);
00965 }
00966
00968 double one_norm_real () const
00969 {
00970 return fvmeta_abs_real(p);
00971 }
00972
00974 double two_norm () const
00975 {
00976 return sqrt(fvmeta_abs2(p));
00977 }
00978
00980 double two_norm2 () const
00981 {
00982 return fvmeta_abs2(p);
00983 }
00984
00986 double infinity_norm () const
00987 {
00988 return std::abs(p);
00989 }
00990
00992 double infinity_norm_real () const
00993 {
00994 return fvmeta_abs_real(p);
00995 }
00996
00997
00998
01000 size_type N () const
01001 {
01002 return 1;
01003 }
01004
01006 size_type dim () const
01007 {
01008 return 1;
01009 }
01010
01011
01012
01014 operator K () {return p;}
01015
01017 operator K () const {return p;}
01018
01019 private:
01020
01021 K p;
01022 };
01023
01024 #ifndef DUNE_EXPRESSIONTEMPLATES
01026 template<class K>
01027 inline FieldVector<K,1> operator+ (const FieldVector<K,1>& a, const FieldVector<K,1>& b)
01028 {
01029 FieldVector<K,1> z = a;
01030 return (z+=b);
01031 }
01032
01034 template<class K>
01035 inline FieldVector<K,1> operator- (const FieldVector<K,1>& a, const FieldVector<K,1>& b)
01036 {
01037 FieldVector<K,1> z = a;
01038 return (z-=b);
01039 }
01040
01042 template<class K>
01043 inline FieldVector<K,1> operator+ (const FieldVector<K,1>& a, const K b)
01044 {
01045 FieldVector<K,1> z = a;
01046 return (z[0]+=b);
01047 }
01048
01050 template<class K>
01051 inline FieldVector<K,1> operator- (const FieldVector<K,1>& a, const K b)
01052 {
01053 FieldVector<K,1> z = a;
01054 return (z[0]-=b);
01055 }
01056
01058 template<class K>
01059 inline FieldVector<K,1> operator+ (const K a, const FieldVector<K,1>& b)
01060 {
01061 FieldVector<K,1> z = a;
01062 return (z[0]+=b);
01063 }
01064
01066 template<class K>
01067 inline FieldVector<K,1> operator- (const K a, const FieldVector<K,1>& b)
01068 {
01069 FieldVector<K,1> z = a;
01070 return (z[0]-=b);
01071 }
01072 #endif
01073
01076 }
01077
01078 #endif