00001 #ifndef DUNE_BVECTOR_HH
00002 #define DUNE_BVECTOR_HH
00003
00004 #include<cmath>
00005 #include<complex>
00006
00007 #include "istlexception.hh"
00008 #include "allocator.hh"
00009 #include "basearray.hh"
00010
00011 #ifdef DUNE_EXPRESSIONTEMPLATES
00012 #include <dune/common/exprtmpl.hh>
00013 #endif
00014
00022 namespace Dune {
00023
00024
00036 template<class B, class A=ISTLAllocator>
00037 class block_vector_unmanaged : public base_array_unmanaged<B,A>
00038 {
00039 public:
00040
00041
00042
00044 typedef typename B::field_type field_type;
00045
00047 typedef B block_type;
00048
00050 typedef A allocator_type;
00051
00053 typedef typename A::size_type size_type;
00054
00056 typedef typename base_array_unmanaged<B,A>::iterator Iterator;
00057
00059 typedef typename base_array_unmanaged<B,A>::const_iterator ConstIterator;
00060
00062 typedef B value_type;
00063
00064
00066
00067 block_vector_unmanaged& operator= (const field_type& k)
00068 {
00069 for (size_type i=0; i<this->n; i++)
00070 (*this)[i] = k;
00071 return *this;
00072 }
00073
00074
00075 #ifndef DUNE_EXPRESSIONTEMPLATES
00077 block_vector_unmanaged& operator+= (const block_vector_unmanaged& y)
00078 {
00079 #ifdef DUNE_ISTL_WITH_CHECKING
00080 if (this->n!=y.N()) DUNE_THROW(ISTLError,"vector size mismatch");
00081 #endif
00082 for (size_type i=0; i<this->n; ++i) (*this)[i] += y[i];
00083 return *this;
00084 }
00085
00087 block_vector_unmanaged& operator-= (const block_vector_unmanaged& y)
00088 {
00089 #ifdef DUNE_ISTL_WITH_CHECKING
00090 if (this->n!=y.N()) DUNE_THROW(ISTLError,"vector size mismatch");
00091 #endif
00092 for (size_type i=0; i<this->n; ++i) (*this)[i] -= y[i];
00093 return *this;
00094 }
00095
00097 block_vector_unmanaged& operator*= (const field_type& k)
00098 {
00099 for (size_type i=0; i<this->n; ++i) (*this)[i] *= k;
00100 return *this;
00101 }
00102
00104 block_vector_unmanaged& operator/= (const field_type& k)
00105 {
00106 for (int i=0; i<this->n; ++i) (*this)[i] /= k;
00107 return *this;
00108 }
00109 #endif
00111 block_vector_unmanaged& axpy (const field_type& a, const block_vector_unmanaged& y)
00112 {
00113 #ifdef DUNE_ISTL_WITH_CHECKING
00114 if (this->n!=y.N()) DUNE_THROW(ISTLError,"vector size mismatch");
00115 #endif
00116 for (size_type i=0; i<this->n; ++i) (*this)[i].axpy(a,y[i]);
00117 return *this;
00118 }
00119
00120
00121
00122
00124 field_type operator* (const block_vector_unmanaged& y) const
00125 {
00126 #ifdef DUNE_ISTL_WITH_CHECKING
00127 if (this->n!=y.N()) DUNE_THROW(ISTLError,"vector size mismatch");
00128 #endif
00129 field_type sum=0;
00130 for (size_type i=0; i<this->n; ++i) sum += (*this)[i]*y[i];
00131 return sum;
00132 }
00133
00134
00135
00136 #ifndef DUNE_EXPRESSIONTEMPLATES
00138 double one_norm () const
00139 {
00140 double sum=0;
00141 for (size_type i=0; i<this->n; ++i) sum += (*this)[i].one_norm();
00142 return sum;
00143 }
00144
00146 double one_norm_real () const
00147 {
00148 double sum=0;
00149 for (size_type i=0; i<this->n; ++i) sum += (*this)[i].one_norm_real();
00150 return sum;
00151 }
00152
00154 double two_norm () const
00155 {
00156 double sum=0;
00157 for (size_type i=0; i<this->n; ++i) sum += (*this)[i].two_norm2();
00158 return sqrt(sum);
00159 }
00160
00162 double two_norm2 () const
00163 {
00164 double sum=0;
00165 for (size_type i=0; i<this->n; ++i) sum += (*this)[i].two_norm2();
00166 return sum;
00167 }
00168
00170 double infinity_norm () const
00171 {
00172 double max=0;
00173 for (size_type i=0; i<this->n; ++i) max = std::max(max,(*this)[i].infinity_norm());
00174 return max;
00175 }
00176
00178 double infinity_norm_real () const
00179 {
00180 double max=0;
00181 for (size_type i=0; i<this->n; ++i) max = std::max(max,(*this)[i].infinity_norm_real());
00182 return max;
00183 }
00184 #endif
00185
00186
00187
00189 size_type N () const
00190 {
00191 return this->n;
00192 }
00193
00195 size_type dim () const
00196 {
00197 size_type d=0;
00198 for (size_type i=0; i<this->n; i++)
00199 d += (*this)[i].dim();
00200 return d;
00201 }
00202
00203 protected:
00205 block_vector_unmanaged () : base_array_unmanaged<B,A>()
00206 { }
00207 };
00208
00224 #ifdef DUNE_EXPRESSIONTEMPLATES
00225 template<class B, class A>
00226 class BlockVector : public block_vector_unmanaged<B,A> ,
00227 public Dune::ExprTmpl::Vector< Dune::BlockVector<B,A> >
00228 #else
00229 template<class B, class A=ISTLAllocator>
00230 class BlockVector : public block_vector_unmanaged<B,A>
00231 #endif
00232 {
00233 public:
00234
00235
00236
00238 typedef typename B::field_type field_type;
00239
00241 typedef B block_type;
00242
00244 typedef A allocator_type;
00245
00247 typedef typename A::size_type size_type;
00248
00250 enum {
00252 blocklevel = B::blocklevel+1};
00253
00255 typedef typename block_vector_unmanaged<B,A>::Iterator Iterator;
00256
00258 typedef typename block_vector_unmanaged<B,A>::ConstIterator ConstIterator;
00259
00260
00261
00263 BlockVector () : block_vector_unmanaged<B,A>(),
00264 capacity_(0)
00265 {
00266 this->n = 0;
00267 }
00268
00270 BlockVector (size_type _n)
00271 {
00272 this->n = _n;
00273 capacity_ = _n;
00274 if (capacity_>0)
00275 this->p = A::template malloc<B>(capacity_);
00276 else
00277 {
00278 this->p = 0;
00279 this->n = 0;
00280 capacity_ = 0;
00281 }
00282 }
00283
00284 BlockVector (size_type _n, size_type capacity)
00285 {
00286 this->n = _n;
00287 if(this->n > capacity)
00288 capacity_ = _n;
00289 else
00290 capacity_ = capacity;
00291
00292 if (capacity_>0)
00293 this->p = A::template malloc<B>(capacity_);
00294 else
00295 {
00296 this->p = 0;
00297 this->n = 0;
00298 capacity_ = 0;
00299 }
00300 }
00301
00302
00319 void reserve(size_type capacity, bool copyOldValues=true)
00320 {
00321 if(capacity >= block_vector_unmanaged<B,A>::N() && capacity != capacity_){
00322
00323 B* pold = this->p;
00324
00325 if(capacity>0){
00326
00327 this->p = A::template malloc<B>(capacity);
00328
00329 if(copyOldValues){
00330
00331 B* to = this->p;
00332 B* from = pold;
00333
00334 for(size_type i=0; i < block_vector_unmanaged<B,A>::N(); ++i, ++from, ++to)
00335 *to = *from;
00336
00337 if(capacity_ > 0)
00338
00339 A::template free<B>(pold);
00340 }
00341 }else{
00342 if(capacity_ > 0)
00343
00344 this->p = 0;
00345 capacity_ = 0;
00346 }
00347
00348 capacity_ = capacity;
00349 }
00350 }
00351
00358 size_type capacity() const
00359 {
00360 return capacity_;
00361 }
00362
00377 void resize(size_type size, bool copyOldValues=true)
00378 {
00379 if(size > block_vector_unmanaged<B,A>::N())
00380 if(capacity_ < size)
00381 this->reserve(size, copyOldValues);
00382
00383 if(size >=0)
00384 this->n=size;
00385 }
00386
00387
00388
00389
00390 #ifdef DUNE_EXPRESSIONTEMPLATES
00392 B& operator[] (size_type i)
00393 {
00394 return block_vector_unmanaged<B,A>::operator[](i);
00395 }
00396
00398 const B& operator[] (size_type i) const
00399 {
00400 return block_vector_unmanaged<B,A>::operator[](i);
00401 }
00402
00404 size_type N () const
00405 {
00406 return block_vector_unmanaged<B,A>::N();
00407 }
00408
00409 BlockVector (const BlockVector& a) {
00410 #ifdef DUNE_VVERBOSE
00411 Dune::dvverb << INDENT << "BlockVector Copy Constructor BlockVector\n";
00412 #endif
00413 capacity_ = a.capacity_;
00414 if (capacity_>0)
00415 this->p = A::template malloc<B>(capacity_);
00416 else
00417 {
00418 this->p = 0;
00419 capacity_ = 0;
00420 }
00421 this->n = a.N();
00422
00423 assignFrom(a);
00424 }
00425 template <class V>
00426 BlockVector (Dune::ExprTmpl::Expression<V> op) {
00427 #ifdef DUNE_VVERBOSE
00428 Dune::dvverb << INDENT << "BlockVector Copy Constructor Expression\n";
00429 #endif
00430 capacity_ = op.N();
00431 if (capacity_>0)
00432 this->p = A::template malloc<B>(capacity_);
00433 else
00434 {
00435 this->p = 0;
00436 capacity_ = 0;
00437 }
00438 this->n = op.N();
00439 assignFrom(op);
00440 }
00441 template <class V>
00442 BlockVector (const Dune::ExprTmpl::Vector<V> & op) {
00443 #ifdef DUNE_VVERBOSE
00444 Dune::dvverb << INDENT << "BlockVector Copy Constructor Vector\n";
00445 #endif
00446 capacity_ = op.N();
00447 if (capacity_>0)
00448 this->p = A::template malloc<B>(capacity_);
00449 else
00450 {
00451 this->p = 0;
00452 capacity_ = 0;
00453 }
00454 this->n = op.N();
00455 assignFrom(op);
00456 }
00457 BlockVector& operator = (const BlockVector& a) {
00458 if (&a!=this)
00459 {
00460 #ifdef DUNE_VVERBOSE
00461 Dune::dvverb << INDENT
00462 << "BlockVector Assignment Operator BlockVector\n";
00463 #endif
00464
00465 if (capacity_!=a.capacity_)
00466 {
00467 #ifdef DUNE_VVERBOSE
00468 Dune::dvverb << INDENT
00469 << "BlockVector Resize to size(BlockVector)\n";
00470 #endif
00471 if (capacity_>0) A::template free<B>(this->p);
00472 capacity_ = a.capacity_;
00473 if (capacity_>0)
00474 this->p = A::template malloc<B>(capacity_);
00475 else
00476 {
00477 this->p = 0;
00478 capacity_ = 0;
00479 }
00480 }
00481 this->n = a.N();
00482
00483 assignFrom(a);
00484 }
00485 return assignFrom(a);
00486 }
00487 template <class E>
00488 BlockVector& operator = (Dune::ExprTmpl::Expression<E> op) {
00489 #ifdef DUNE_VVERBOSE
00490 Dune::dvverb << INDENT << "BlockVector Assignment Operator Expression\n";
00491 #endif
00492
00493 if (capacity_ < op.N())
00494 {
00495 #ifdef DUNE_VVERBOSE
00496 Dune::dvverb << INDENT << "BlockVector Resize to size(Expression)\n";
00497 #endif
00498 if (capacity_>0) A::template free<B>(this->p);
00499 capacity_ = op.N();
00500 if (capacity_>0)
00501 this->p = A::template malloc<B>(capacity_);
00502 else
00503 {
00504 this->p = 0;
00505 capacity_ = 0;
00506 }
00507 }
00508 this->n = op.N();
00509 return assignFrom(op);
00510 }
00511 template <class V>
00512 BlockVector& operator = (const Dune::ExprTmpl::Vector<V> & op) {
00513 #ifdef DUNE_VVERBOSE
00514 Dune::dvverb << INDENT << "BlockVector Assignment Operator Vector\n";
00515 #endif
00516
00517 if (capacity_ < op.N())
00518 {
00519 #ifdef DUNE_VVERBOSE
00520 Dune::dvverb << INDENT << "BlockVector Resize to size(Vector)\n";
00521 #endif
00522 if (capacity_>0) A::template free<B>(this->p);
00523 capacity_ = op.N();
00524 if (capacity_>0)
00525 this->p = A::template malloc<B>(capacity_);
00526 else
00527 {
00528 this->p = 0;
00529 capacity_ = 0;
00530 }
00531 }
00532 this->n = op.N();
00533 return assignFrom(op);
00534 }
00535 #endif
00536
00537 #ifndef DUNE_EXPRESSIONTEMPLATES
00539 BlockVector (const BlockVector& a) :
00540 block_vector_unmanaged<B,A>(a)
00541 {
00542
00543 this->n = a.n;
00544 capacity_ = a.capacity_;
00545
00546 if (capacity_>0)
00547 this->p = A::template malloc<B>(capacity_);
00548 else
00549 {
00550 this->n = 0;
00551 this->p = 0;
00552 }
00553
00554
00555 for (size_type i=0; i<this->n; i++) this->p[i]=a.p[i];
00556 }
00557
00559 BlockVector (const block_vector_unmanaged<B,A>& _a)
00560 {
00561
00562 const BlockVector& a = static_cast<const BlockVector&>(_a);
00563
00564
00565 this->n = a.n;
00566 capacity_ = a.capacity_;
00567
00568 if (capacity_>0)
00569 this->p = A::template malloc<B>(capacity_);
00570 else
00571 {
00572 this->n = 0;
00573 this->p = 0;
00574 capacity_ = 0;
00575 }
00576
00577
00578 for (size_type i=0; i<this->n; i++) this->p[i]=a.p[i];
00579 }
00580 #endif
00581
00583 ~BlockVector ()
00584 {
00585 if (capacity_>0) A::template free<B>(this->p);
00586 }
00587
00589 #ifndef DUNE_EXPRESSIONTEMPLATES
00590 BlockVector& operator= (const BlockVector& a)
00591 {
00592 if (&a!=this)
00593 {
00594
00595 if (capacity_!=a.capacity_)
00596 {
00597 if (capacity_>0) A::template free<B>(this->p);
00598 capacity_ = a.capacity_;
00599 if (capacity_>0)
00600 this->p = A::template malloc<B>(capacity_);
00601 else
00602 {
00603 this->p = 0;
00604 capacity_ = 0;
00605 }
00606 }
00607 this->n = a.n;
00608
00609 for (size_type i=0; i<this->n; i++)
00610 this->p[i]=a.p[i];
00611 }
00612 return *this;
00613 }
00614
00616 BlockVector& operator= (const block_vector_unmanaged<B,A>& a)
00617 {
00618
00619 return this->operator=(static_cast<const BlockVector&>(a));
00620 }
00621 #endif
00622
00624 BlockVector& operator= (const field_type& k)
00625 {
00626
00627 (static_cast<block_vector_unmanaged<B,A>&>(*this)) = k;
00628 return *this;
00629 }
00630 protected:
00631 size_type capacity_;
00632 };
00633
00635
00636 template<class K, class A>
00637 std::ostream& operator<< (std::ostream& s, const BlockVector<K, A>& v)
00638 {
00639 typedef typename BlockVector<K, A>::size_type size_type;
00640
00641 for (size_type i=0; i<v.size(); i++)
00642 s << v[i] << std::endl;
00643
00644 return s;
00645 }
00646
00663 template<class B, class A=ISTLAllocator>
00664 class BlockVectorWindow : public block_vector_unmanaged<B,A>
00665 {
00666 public:
00667
00668
00669
00671 typedef typename B::field_type field_type;
00672
00674 typedef B block_type;
00675
00677 typedef A allocator_type;
00678
00680 typedef typename A::size_type size_type;
00681
00683 enum {
00685 blocklevel = B::blocklevel+1
00686 };
00687
00689 typedef typename block_vector_unmanaged<B,A>::Iterator Iterator;
00690
00692 typedef typename block_vector_unmanaged<B,A>::ConstIterator ConstIterator;
00693
00694
00695
00697 BlockVectorWindow () : block_vector_unmanaged<B,A>()
00698 { }
00699
00701 BlockVectorWindow (B* _p, size_type _n)
00702 {
00703 this->n = _n;
00704 this->p = _p;
00705 }
00706
00708 BlockVectorWindow (const BlockVectorWindow& a)
00709 {
00710 this->n = a.n;
00711 this->p = a.p;
00712 }
00713
00715 BlockVectorWindow (const block_vector_unmanaged<B,A>& _a)
00716 {
00717
00718 const BlockVectorWindow& a = static_cast<const BlockVectorWindow&>(_a);
00719
00720
00721 this->n = a.n;
00722 this->p = a.p;
00723 }
00724
00725
00727 BlockVectorWindow& operator= (const BlockVectorWindow& a)
00728 {
00729
00730 #ifdef DUNE_ISTL_WITH_CHECKING
00731 if (this->n!=a.N()) DUNE_THROW(ISTLError,"vector size mismatch");
00732 #endif
00733
00734 if (&a!=this)
00735 {
00736
00737 for (size_type i=0; i<this->n; i++) this->p[i]=a.p[i];
00738 }
00739 return *this;
00740 }
00741
00743 BlockVectorWindow& operator= (const block_vector_unmanaged<B,A>& a)
00744 {
00745
00746 return this->operator=(static_cast<const BlockVectorWindow&>(a));
00747 }
00748
00750 BlockVectorWindow& operator= (const field_type& k)
00751 {
00752 (static_cast<block_vector_unmanaged<B,A>&>(*this)) = k;
00753 return *this;
00754 }
00755
00756
00757
00758
00760 void set (size_type _n, B* _p)
00761 {
00762 this->n = _n;
00763 this->p = _p;
00764 }
00765
00767 void setsize (size_type _n)
00768 {
00769 this->n = _n;
00770 }
00771
00773 void setptr (B* _p)
00774 {
00775 this->p = _p;
00776 }
00777
00779 B* getptr ()
00780 {
00781 return this->p;
00782 }
00783
00785 size_type getsize ()
00786 {
00787 return this->n;
00788 }
00789 };
00790
00791
00792
00801 template<class B, class A=ISTLAllocator>
00802 class compressed_block_vector_unmanaged : public compressed_base_array_unmanaged<B,A>
00803 {
00804 public:
00805
00806
00807
00809 typedef typename B::field_type field_type;
00810
00812 typedef B block_type;
00813
00815 typedef A allocator_type;
00816
00818 typedef typename compressed_base_array_unmanaged<B,A>::iterator Iterator;
00819
00821 typedef typename compressed_base_array_unmanaged<B,A>::const_iterator ConstIterator;
00822
00824 typedef typename A::size_type size_type;
00825
00826
00827
00828 compressed_block_vector_unmanaged& operator= (const field_type& k)
00829 {
00830 for (size_type i=0; i<this->n; i++)
00831 (this->p)[i] = k;
00832 return *this;
00833 }
00834
00835
00836
00837
00839 template<class V>
00840 compressed_block_vector_unmanaged& operator+= (const V& y)
00841 {
00842 #ifdef DUNE_ISTL_WITH_CHECKING
00843 if (!includesindexset(y)) DUNE_THROW(ISTLError,"index set mismatch");
00844 #endif
00845 for (size_type i=0; i<y.n; ++i) this->operator[](y.j[i]) += y.p[i];
00846 return *this;
00847 }
00848
00850 template<class V>
00851 compressed_block_vector_unmanaged& operator-= (const V& y)
00852 {
00853 #ifdef DUNE_ISTL_WITH_CHECKING
00854 if (!includesindexset(y)) DUNE_THROW(ISTLError,"index set mismatch");
00855 #endif
00856 for (size_type i=0; i<y.n; ++i) this->operator[](y.j[i]) -= y.p[i];
00857 return *this;
00858 }
00859
00861 template<class V>
00862 compressed_block_vector_unmanaged& axpy (const field_type& a, const V& y)
00863 {
00864 #ifdef DUNE_ISTL_WITH_CHECKING
00865 if (!includesindexset(y)) DUNE_THROW(ISTLError,"index set mismatch");
00866 #endif
00867 for (size_type i=0; i<y.n; ++i) (this->operator[](y.j[i])).axpy(a,y.p[i]);
00868 return *this;
00869 }
00870
00872 compressed_block_vector_unmanaged& operator*= (const field_type& k)
00873 {
00874 for (size_type i=0; i<this->n; ++i) (this->p)[i] *= k;
00875 return *this;
00876 }
00877
00879 compressed_block_vector_unmanaged& operator/= (const field_type& k)
00880 {
00881 for (size_type i=0; i<this->n; ++i) (this->p)[i] /= k;
00882 return *this;
00883 }
00884
00885
00886
00887
00889 field_type operator* (const compressed_block_vector_unmanaged& y) const
00890 {
00891 #ifdef DUNE_ISTL_WITH_CHECKING
00892 if (!includesindexset(y) || !y.includesindexset(*this) )
00893 DUNE_THROW(ISTLError,"index set mismatch");
00894 #endif
00895 field_type sum=0;
00896 for (size_type i=0; i<this->n; ++i)
00897 sum += (this->p)[i] * y[(this->j)[i]];
00898 return sum;
00899 }
00900
00901
00902
00903
00905 double one_norm () const
00906 {
00907 double sum=0;
00908 for (size_type i=0; i<this->n; ++i) sum += (this->p)[i].one_norm();
00909 return sum;
00910 }
00911
00913 double one_norm_real () const
00914 {
00915 double sum=0;
00916 for (size_type i=0; i<this->n; ++i) sum += (this->p)[i].one_norm_real();
00917 return sum;
00918 }
00919
00921 double two_norm () const
00922 {
00923 double sum=0;
00924 for (size_type i=0; i<this->n; ++i) sum += (this->p)[i].two_norm2();
00925 return sqrt(sum);
00926 }
00927
00929 double two_norm2 () const
00930 {
00931 double sum=0;
00932 for (size_type i=0; i<this->n; ++i) sum += (this->p)[i].two_norm2();
00933 return sum;
00934 }
00935
00937 double infinity_norm () const
00938 {
00939 double max=0;
00940 for (size_type i=0; i<this->n; ++i) max = std::max(max,(this->p)[i].infinity_norm());
00941 return max;
00942 }
00943
00945 double infinity_norm_real () const
00946 {
00947 double max=0;
00948 for (size_type i=0; i<this->n; ++i) max = std::max(max,(this->p)[i].infinity_norm_real());
00949 return max;
00950 }
00951
00952
00953
00954
00956 size_type N () const
00957 {
00958 return this->n;
00959 }
00960
00962 size_type dim () const
00963 {
00964 size_type d=0;
00965 for (size_type i=0; i<this->n; i++)
00966 d += (this->p)[i].dim();
00967 return d;
00968 }
00969
00970 protected:
00972 compressed_block_vector_unmanaged () : compressed_base_array_unmanaged<B,A>()
00973 { }
00974
00976 template<class V>
00977 bool includesindexset (const V& y)
00978 {
00979 typename V::ConstIterator e=this->end();
00980 for (size_type i=0; i<y.n; i++)
00981 if (find(y.j[i])==e)
00982 return false;
00983 return true;
00984 }
00985 };
00986
00987
01004 template<class B, class A=ISTLAllocator>
01005 class CompressedBlockVectorWindow : public compressed_block_vector_unmanaged<B,A>
01006 {
01007 public:
01008
01009
01010
01012 typedef typename B::field_type field_type;
01013
01015 typedef B block_type;
01016
01018 typedef A allocator_type;
01019
01021 typedef typename A::size_type size_type;
01022
01024 enum {
01026 blocklevel = B::blocklevel+1};
01027
01029 typedef typename compressed_block_vector_unmanaged<B,A>::Iterator Iterator;
01030
01032 typedef typename compressed_block_vector_unmanaged<B,A>::ConstIterator ConstIterator;
01033
01034
01035
01037 CompressedBlockVectorWindow () : compressed_block_vector_unmanaged<B,A>()
01038 { }
01039
01041 CompressedBlockVectorWindow (B* _p, size_type* _j, size_type _n)
01042 {
01043 this->n = _n;
01044 this->p = _p;
01045 this->j = _j;
01046 }
01047
01049 CompressedBlockVectorWindow (const CompressedBlockVectorWindow& a)
01050 {
01051 this->n = a.n;
01052 this->p = a.p;
01053 this->j = a.j;
01054 }
01055
01057 CompressedBlockVectorWindow (const compressed_block_vector_unmanaged<B,A>& _a)
01058 {
01059
01060 const CompressedBlockVectorWindow& a = static_cast<const CompressedBlockVectorWindow&>(_a);
01061
01062
01063 this->n = a.n;
01064 this->p = a.p;
01065 this->j = a.j;
01066 }
01067
01068
01070 CompressedBlockVectorWindow& operator= (const CompressedBlockVectorWindow& a)
01071 {
01072
01073 #ifdef DUNE_ISTL_WITH_CHECKING
01074 if (this->n!=a.N()) DUNE_THROW(ISTLError,"vector size mismatch");
01075 #endif
01076
01077 if (&a!=this)
01078 {
01079
01080 for (size_type i=0; i<this->n; i++) this->p[i]=a.p[i];
01081 for (size_type i=0; i<this->n; i++) this->j[i]=a.j[i];
01082 }
01083 return *this;
01084 }
01085
01087 CompressedBlockVectorWindow& operator= (const compressed_block_vector_unmanaged<B,A>& a)
01088 {
01089
01090 return this->operator=(static_cast<const CompressedBlockVectorWindow&>(a));
01091 }
01092
01094 CompressedBlockVectorWindow& operator= (const field_type& k)
01095 {
01096 (static_cast<compressed_block_vector_unmanaged<B,A>&>(*this)) = k;
01097 return *this;
01098 }
01099
01100
01101
01102
01104 void set (size_type _n, B* _p, size_type* _j)
01105 {
01106 this->n = _n;
01107 this->p = _p;
01108 this->j = _j;
01109 }
01110
01112 void setsize (size_type _n)
01113 {
01114 this->n = _n;
01115 }
01116
01118 void setptr (B* _p)
01119 {
01120 this->p = _p;
01121 }
01122
01124 void setindexptr (size_type* _j)
01125 {
01126 this->j = _j;
01127 }
01128
01130 B* getptr ()
01131 {
01132 return this->p;
01133 }
01134
01136 size_type* getindexptr ()
01137 {
01138 return this->j;
01139 }
01140
01142 const B* getptr () const
01143 {
01144 return this->p;
01145 }
01146
01148 const size_type* getindexptr () const
01149 {
01150 return this->j;
01151 }
01153 size_type getsize () const
01154 {
01155 return this->n;
01156 }
01157 };
01158
01159 #ifdef DUNE_EXPRESSIONTEMPLATES
01160 template <class B, class A>
01161 struct FieldType< BlockVector<B,A> >
01162 {
01163 typedef typename FieldType<B>::type type;
01164 };
01165 template <class B, class A>
01166 struct BlockType< BlockVector<B,A> >
01167 {
01168 typedef B type;
01169 };
01170 #endif
01171
01172 }
01173
01174 #endif