4 #ifndef DUNE_ISTL_BCRSMATRIX_HH 5 #define DUNE_ISTL_BCRSMATRIX_HH 19 #include <dune/common/stdstreams.hh> 20 #include <dune/common/iteratorfacades.hh> 21 #include <dune/common/typetraits.hh> 22 #include <dune/common/ftraits.hh> 79 template<
typename size_type>
137 return _m.entry(_i,j);
166 if (m.buildMode() != Matrix::implicit)
167 DUNE_THROW(
BCRSMatrixError,
"You can only create an ImplicitBuilder for a matrix in implicit build mode");
168 if (m.buildStage() != Matrix::building)
169 DUNE_THROW(
BCRSMatrixError,
"You can only create an ImplicitBuilder for a matrix with set size that has not been compressed() yet");
187 ImplicitMatrixBuilder(Matrix& m, size_type rows, size_type cols, size_type avg_cols_per_row,
double overflow_fraction)
190 if (m.buildStage() != Matrix::notAllocated)
191 DUNE_THROW(
BCRSMatrixError,
"You can only set up a matrix for this ImplicitBuilder if it has no memory allocated yet");
192 m.setBuildMode(Matrix::implicit);
193 m.setImplicitBuildModeParameters(avg_cols_per_row,overflow_fraction);
194 m.setSize(rows,cols);
421 template<
class B,
class A=std::allocator<B> >
465 blocklevel = B::blocklevel+1
508 row_type& operator[] (size_type i)
510 #ifdef DUNE_ISTL_WITH_CHECKING 511 if (build_mode == implicit && ready != built)
512 DUNE_THROW(
BCRSMatrixError,
"You cannot use operator[] in implicit build mode before calling compress()");
520 const row_type& operator[] (size_type i)
const 522 #ifdef DUNE_ISTL_WITH_CHECKING 523 if (build_mode == implicit && ready != built)
524 DUNE_THROW(
BCRSMatrixError,
"You cannot use operator[] in implicit build mode before calling compress()");
525 if (built!=ready) DUNE_THROW(
BCRSMatrixError,
"row not initialized yet");
537 :
public RandomAccessIteratorFacade<RealRowIterator<T>, T>
544 friend class RandomAccessIteratorFacade<
RealRowIterator<const ValueType>, const ValueType>;
609 void advance(std::ptrdiff_t diff)
614 T& elementAt(std::ptrdiff_t diff)
const 620 row_type& dereference ()
const 636 return Iterator(r,0);
642 return Iterator(r,n);
649 return Iterator(r,n-1);
656 return Iterator(r,-1);
673 return ConstIterator(r,0);
677 ConstIterator
end ()
const 679 return ConstIterator(r,n);
686 return ConstIterator(r,n-1);
693 return ConstIterator(r,-1);
709 : build_mode(unknown), ready(notAllocated), n(0), m(0), nnz_(0),
710 allocationSize_(0), r(0), a(0),
711 avg(0), overflowsize(-1.0)
716 : build_mode(bm), ready(notAllocated), n(0), m(0), nnz_(0),
717 allocationSize_(0), r(0), a(0),
718 avg(0), overflowsize(-1.0)
720 allocate(_n, _m, _nnz,
true,
false);
725 : build_mode(bm), ready(notAllocated), n(0), m(0), nnz_(0),
726 allocationSize_(0), r(0), a(0),
727 avg(0), overflowsize(-1.0)
729 allocate(_n, _m,0,
true,
false);
744 : build_mode(bm), ready(notAllocated), n(0), m(0), nnz_(0),
745 allocationSize_(0), r(0), a(0),
746 avg(_avg), overflowsize(_overflowsize)
749 DUNE_THROW(
BCRSMatrixError,
"Only call this constructor when using the implicit build mode");
754 if (_overflowsize < 0.0)
755 DUNE_THROW(
BCRSMatrixError,
"You cannot set a negative overflow fraction");
756 implicit_allocate(_n,_m);
765 : build_mode(Mat.build_mode), ready(notAllocated), n(0), m(0), nnz_(0),
766 allocationSize_(0), r(0), a(0),
767 avg(Mat.
avg), overflowsize(Mat.overflowsize)
769 if (!(Mat.
ready == notAllocated || Mat.
ready == built))
770 DUNE_THROW(InvalidStateException,
"BCRSMatrix can only be copy-constructed when source matrix is completely empty (size not set) or fully built)");
773 size_type _nnz = Mat.
nnz_;
779 for (size_type i=0; i<Mat.
n; i++)
784 allocate(Mat.
n, Mat.
m, _nnz,
true,
true);
787 copyWindowStructure(Mat);
802 if (ready == notAllocated)
807 if (ready == building && (build_mode == unknown || build_mode == random || build_mode == row_wise) && (bm == row_wise || bm == random))
810 DUNE_THROW(InvalidStateException,
"Matrix structure cannot be changed at this stage anymore (ready == "<<ready<<
").");
828 void setSize(size_type rows, size_type columns, size_type nnz=0)
833 if (build_mode == implicit)
836 DUNE_THROW(
Dune::BCRSMatrixError,
"number of non-zeroes may not be set in implicit mode, use setImplicitBuildModeParameters() instead");
839 implicit_allocate(rows,columns);
844 allocate(rows, columns, nnz,
true,
false);
863 DUNE_THROW(
BCRSMatrixError,
"You cannot set a negative overflow fraction");
866 if (ready != notAllocated)
867 DUNE_THROW(InvalidStateException,
"You cannot modify build mode parameters at this stage anymore");
869 overflowsize = _overflow;
881 if (&Mat==
this)
return *
this;
883 if (!((ready == notAllocated || ready == built) && (Mat.
ready == notAllocated || Mat.
ready == built)))
884 DUNE_THROW(InvalidStateException,
"BCRSMatrix can only be copied when both target and source are empty or fully built)");
888 deallocate(n!=Mat.
n);
891 if (n>0 && n!=Mat.
n) {
893 for(row_type *riter=r+(n-1), *rend=r-1; riter!=rend; --riter)
894 rowAllocator_.destroy(riter);
895 rowAllocator_.deallocate(r,n);
901 for (size_type i=0; i<Mat.
n; i++)
907 allocate(Mat.
n, Mat.
m, nnz_, n!=Mat.
n,
true);
910 copyWindowStructure(Mat);
918 if (!(ready == notAllocated || ready == built))
919 DUNE_THROW(InvalidStateException,
"Scalar assignment only works on fully built BCRSMatrix)");
921 for (size_type i=0; i<n; i++) r[i] = k;
933 : Mat(_Mat), i(_i), nnz(0), current_row(nullptr, Mat.j_.
get(), 0)
935 if (Mat.build_mode == unknown && Mat.ready == building)
937 Mat.build_mode = row_wise;
939 if (i==0 && Mat.ready != building)
940 DUNE_THROW(
BCRSMatrixError,
"creation only allowed for uninitialized matrix");
941 if(Mat.build_mode!=row_wise)
942 DUNE_THROW(
BCRSMatrixError,
"creation only allowed if row wise allocation was requested in the constructor");
949 if (Mat.ready != building)
973 Mat.r[i].set(s,
nullptr,current_row.getindexptr());
974 current_row.setindexptr(current_row.getindexptr()+s);
978 B* b = Mat.allocator_.allocate(s);
982 size_type* j = Mat.sizeAllocator_.allocate(s);
987 Mat.r[i].set(0,
nullptr,
nullptr);
990 std::copy(
pattern.cbegin(),
pattern.cend(), Mat.r[i].getindexptr());
1007 Mat.setDataPointers();
1017 return (i!=it.i) || (&Mat!=&it.Mat);
1023 return (i==it.i) && (&Mat==&it.Mat);
1057 typedef std::set<size_type,std::less<size_type> > PatternType;
1059 row_type current_row;
1088 if (build_mode!=random)
1090 if (ready != building)
1099 #ifdef DUNE_ISTL_WITH_CHECKING 1103 return r[i].getsize();
1109 if (build_mode!=random)
1111 if (ready != building)
1114 r[i].setsize(r[i].getsize()+s);
1120 if (build_mode!=random)
1122 if (ready != building)
1127 for (size_type i=0; i<n; i++)
1129 total += r[i].getsize();
1134 allocate(n,m,total,
false,
false);
1135 else if(nnz_ < total)
1136 DUNE_THROW(
BCRSMatrixError,
"Specified number of nonzeros ("<<nnz_<<
") not " 1137 <<
"sufficient for calculated nonzeros ("<<total<<
"! ");
1140 setColumnPointers(begin());
1144 for (size_type k=0; k<nnz_; k++)
1146 ready = rowSizesBuilt;
1162 if (build_mode!=random)
1166 if (ready==building)
1168 if (ready==notAllocated)
1169 DUNE_THROW(
BCRSMatrixError,
"matrix size not set and no memory allocated yet");
1175 size_type*
const first = r[row].getindexptr();
1176 size_type*
const last = first + r[row].getsize();
1179 size_type* pos = std::lower_bound(first,last,col);
1182 if (pos!=last && *pos == col)
return;
1185 size_type* end = std::lower_bound(pos,last,m);
1190 std::copy_backward(pos,end,end+1);
1202 template<
typename It>
1205 size_type row_size = r[row].size();
1206 size_type* col_begin = r[row].getindexptr();
1209 if ((col_end = std::copy(begin,end,r[row].getindexptr())) != col_begin + row_size)
1212 <<
") does not match number of passed entries (" << (col_end - col_begin) <<
")");
1213 std::sort(col_begin,col_end);
1219 if (build_mode!=random)
1223 if (ready==building)
1225 if (ready==notAllocated)
1226 DUNE_THROW(
BCRSMatrixError,
"matrix size not set and no memory allocated yet");
1229 RowIterator endi=end();
1230 for (RowIterator i=begin(); i!=endi; ++i)
1232 ColIterator endj = (*i).end();
1233 for (ColIterator j=(*i).begin(); j!=endj; ++j) {
1234 if (j.index() >= m) {
1235 dwarn <<
"WARNING: size of row "<< i.index()<<
" is "<<j.offset()<<
". But was specified as being "<< (*i).end().offset()
1236 <<
". This means you are wasting valuable space and creating additional cache misses!"<<std::endl;
1237 r[i.index()].setsize(j.offset());
1266 #ifdef DUNE_ISTL_WITH_CHECKING 1267 if (build_mode!=implicit)
1270 DUNE_THROW(
BCRSMatrixError,
"matrix already built up, use operator[] for entry access now");
1271 if (ready==notAllocated)
1272 DUNE_THROW(
BCRSMatrixError,
"matrix size not set and no memory allocated yet");
1273 if (ready!=building)
1274 DUNE_THROW(InvalidStateException,
"You may only use entry() during the 'building' stage");
1282 size_type* begin = r[row].getindexptr();
1283 size_type* end = begin + r[row].getsize();
1285 size_type* pos = std::find(begin, end, col);
1291 std::ptrdiff_t offset = pos - r[row].getindexptr();
1292 B* aptr = r[row].getptr() + offset;
1298 if (r[row].getsize() ==
avg)
1299 return overflow[std::make_pair(row,col)];
1306 std::ptrdiff_t offset = end - r[row].getindexptr();
1307 B* apos = r[row].getptr() + offset;
1310 r[row].setsize(r[row].getsize()+1);
1330 if (build_mode!=implicit)
1333 DUNE_THROW(
BCRSMatrixError,
"matrix already built up, no more need for compression");
1334 if (ready==notAllocated)
1335 DUNE_THROW(
BCRSMatrixError,
"matrix size not set and no memory allocated yet");
1336 if (ready!=building)
1337 DUNE_THROW(InvalidStateException,
"You may only call compress() at the end of the 'building' stage");
1340 CompressionStatistics stats;
1345 size_type* jiit = j_.get();
1349 typename OverflowType::iterator oit = overflow.begin();
1352 std::vector<size_type*> perm;
1355 for (size_type i=0; i<n; i++)
1358 size_type* begin = r[i].getindexptr();
1360 size_type size = r[i].getsize();
1364 typename std::vector<size_type*>::iterator it = perm.begin();
1365 for (size_type* iit = begin; iit < begin + size; ++iit, ++it)
1372 r[i].setindexptr(jiit);
1375 for (it = perm.begin(); it != perm.end(); ++it)
1378 while ((oit!=overflow.end()) && (oit->first < std::make_pair(i,**it)))
1383 "Allocated memory for BCRSMatrix exhausted during compress()!" 1384 "Please increase either the average number of entries per row or the overflow fraction." 1387 *jiit = oit->first.second;
1389 *aiit = oit->second;
1392 r[i].setsize(r[i].getsize()+1);
1398 "Allocated memory for BCRSMatrix exhausted during compress()!" 1399 "Please increase either the average number of entries per row or the overflow fraction." 1405 B* apos = *it - j_.get() + a;
1411 while ((oit!=overflow.end()) && (oit->first.first == i))
1416 "Allocated memory for BCRSMatrix exhausted during compress()!" 1417 "Please increase either the average number of entries per row or the overflow fraction." 1421 *jiit = oit->first.second;
1423 *aiit = oit->second;
1426 r[i].setsize(r[i].getsize()+1);
1430 if (r[i].getsize()>stats.
maximum)
1431 stats.
maximum = r[i].getsize();
1438 std::ptrdiff_t diff = (r[n-1].getindexptr() + r[n-1].getsize() - j_.get());
1440 stats.
avg = (double) (nnz_) / (double) n;
1441 stats.
mem_ratio = (double) (nnz_) / (double) allocationSize_;
1454 #ifdef DUNE_ISTL_WITH_CHECKING 1456 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1462 for (size_type i=0; i<nnz_; i++)
1467 RowIterator endi=end();
1468 for (RowIterator i=begin(); i!=endi; ++i)
1470 ColIterator endj = (*i).end();
1471 for (ColIterator j=(*i).begin(); j!=endj; ++j)
1482 #ifdef DUNE_ISTL_WITH_CHECKING 1484 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1490 for (size_type i=0; i<nnz_; i++)
1495 RowIterator endi=end();
1496 for (RowIterator i=begin(); i!=endi; ++i)
1498 ColIterator endj = (*i).end();
1499 for (ColIterator j=(*i).begin(); j!=endj; ++j)
1515 #ifdef DUNE_ISTL_WITH_CHECKING 1516 if (ready != built || b.
ready != built)
1517 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1518 if(N()!=b.
N() || M() != b.
M())
1519 DUNE_THROW(RangeError,
"Matrix sizes do not match!");
1521 RowIterator endi=end();
1522 ConstRowIterator j=b.
begin();
1523 for (RowIterator i=begin(); i!=endi; ++i, ++j) {
1537 #ifdef DUNE_ISTL_WITH_CHECKING 1538 if (ready != built || b.
ready != built)
1539 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1540 if(N()!=b.
N() || M() != b.
M())
1541 DUNE_THROW(RangeError,
"Matrix sizes do not match!");
1543 RowIterator endi=end();
1544 ConstRowIterator j=b.
begin();
1545 for (RowIterator i=begin(); i!=endi; ++i, ++j) {
1562 #ifdef DUNE_ISTL_WITH_CHECKING 1563 if (ready != built || b.
ready != built)
1564 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1565 if(N()!=b.
N() || M() != b.
M())
1566 DUNE_THROW(RangeError,
"Matrix sizes do not match!");
1568 RowIterator endi=end();
1569 ConstRowIterator j=b.
begin();
1570 for(RowIterator i=begin(); i!=endi; ++i, ++j)
1579 template<
class X,
class Y>
1580 void mv (
const X& x, Y& y)
const 1582 #ifdef DUNE_ISTL_WITH_CHECKING 1584 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1586 "Size mismatch: M: " << N() <<
"x" << M() <<
" x: " << x.N());
1588 "Size mismatch: M: " << N() <<
"x" << M() <<
" y: " << y.N());
1590 ConstRowIterator endi=end();
1591 for (ConstRowIterator i=begin(); i!=endi; ++i)
1594 ConstColIterator endj = (*i).end();
1595 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1596 (*j).umv(x[j.index()],y[i.index()]);
1601 template<
class X,
class Y>
1602 void umv (
const X& x, Y& y)
const 1604 #ifdef DUNE_ISTL_WITH_CHECKING 1606 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1610 ConstRowIterator endi=end();
1611 for (ConstRowIterator i=begin(); i!=endi; ++i)
1613 ConstColIterator endj = (*i).end();
1614 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1615 (*j).umv(x[j.index()],y[i.index()]);
1620 template<
class X,
class Y>
1621 void mmv (
const X& x, Y& y)
const 1623 #ifdef DUNE_ISTL_WITH_CHECKING 1625 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1629 ConstRowIterator endi=end();
1630 for (ConstRowIterator i=begin(); i!=endi; ++i)
1632 ConstColIterator endj = (*i).end();
1633 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1634 (*j).mmv(x[j.index()],y[i.index()]);
1639 template<
class X,
class Y,
class F>
1640 void usmv (F&& alpha,
const X& x, Y& y)
const 1642 #ifdef DUNE_ISTL_WITH_CHECKING 1644 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1648 ConstRowIterator endi=end();
1649 for (ConstRowIterator i=begin(); i!=endi; ++i)
1651 ConstColIterator endj = (*i).end();
1652 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1653 (*j).usmv(alpha,x[j.index()],y[i.index()]);
1658 template<
class X,
class Y>
1659 void mtv (
const X& x, Y& y)
const 1661 #ifdef DUNE_ISTL_WITH_CHECKING 1663 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1667 for(size_type i=0; i<y.N(); ++i)
1673 template<
class X,
class Y>
1676 #ifdef DUNE_ISTL_WITH_CHECKING 1678 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1682 ConstRowIterator endi=end();
1683 for (ConstRowIterator i=begin(); i!=endi; ++i)
1685 ConstColIterator endj = (*i).end();
1686 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1687 (*j).umtv(x[i.index()],y[j.index()]);
1692 template<
class X,
class Y>
1695 #ifdef DUNE_ISTL_WITH_CHECKING 1699 ConstRowIterator endi=end();
1700 for (ConstRowIterator i=begin(); i!=endi; ++i)
1702 ConstColIterator endj = (*i).end();
1703 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1704 (*j).mmtv(x[i.index()],y[j.index()]);
1709 template<
class X,
class Y>
1710 void usmtv (
const field_type& alpha,
const X& x, Y& y)
const 1712 #ifdef DUNE_ISTL_WITH_CHECKING 1714 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1718 ConstRowIterator endi=end();
1719 for (ConstRowIterator i=begin(); i!=endi; ++i)
1721 ConstColIterator endj = (*i).end();
1722 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1723 (*j).usmtv(alpha,x[i.index()],y[j.index()]);
1728 template<
class X,
class Y>
1731 #ifdef DUNE_ISTL_WITH_CHECKING 1733 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1737 ConstRowIterator endi=end();
1738 for (ConstRowIterator i=begin(); i!=endi; ++i)
1740 ConstColIterator endj = (*i).end();
1741 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1742 (*j).umhv(x[i.index()],y[j.index()]);
1747 template<
class X,
class Y>
1750 #ifdef DUNE_ISTL_WITH_CHECKING 1752 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1756 ConstRowIterator endi=end();
1757 for (ConstRowIterator i=begin(); i!=endi; ++i)
1759 ConstColIterator endj = (*i).end();
1760 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1761 (*j).mmhv(x[i.index()],y[j.index()]);
1766 template<
class X,
class Y>
1767 void usmhv (
const field_type& alpha,
const X& x, Y& y)
const 1769 #ifdef DUNE_ISTL_WITH_CHECKING 1771 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1775 ConstRowIterator endi=end();
1776 for (ConstRowIterator i=begin(); i!=endi; ++i)
1778 ConstColIterator endj = (*i).end();
1779 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1780 (*j).usmhv(alpha,x[i.index()],y[j.index()]);
1790 #ifdef DUNE_ISTL_WITH_CHECKING 1792 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1797 ConstRowIterator endi=end();
1798 for (ConstRowIterator i=begin(); i!=endi; ++i)
1800 ConstColIterator endj = (*i).end();
1801 for (ConstColIterator j=(*i).begin(); j!=endj; ++j)
1802 sum += (*j).frobenius_norm2();
1811 return sqrt(frobenius_norm2());
1815 template <
typename ft = field_type,
1816 typename std::enable_if<!has_nan<ft>::value,
int>::type = 0>
1819 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1821 using real_type =
typename FieldTraits<ft>::real_type;
1825 for (
auto const &x : *
this) {
1827 for (
auto const &y : x)
1828 sum += y.infinity_norm();
1829 norm = max(sum, norm);
1835 template <
typename ft = field_type,
1836 typename std::enable_if<!has_nan<ft>::value,
int>::type = 0>
1839 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1841 using real_type =
typename FieldTraits<ft>::real_type;
1845 for (
auto const &x : *
this) {
1847 for (
auto const &y : x)
1848 sum += y.infinity_norm_real();
1849 norm = max(sum, norm);
1855 template <
typename ft = field_type,
1856 typename std::enable_if<has_nan<ft>::value,
int>::type = 0>
1859 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1861 using real_type =
typename FieldTraits<ft>::real_type;
1865 real_type isNaN = 1;
1866 for (
auto const &x : *
this) {
1868 for (
auto const &y : x)
1869 sum += y.infinity_norm();
1870 norm = max(sum, norm);
1874 return norm * isNaN;
1878 template <
typename ft = field_type,
1879 typename std::enable_if<has_nan<ft>::value,
int>::type = 0>
1882 DUNE_THROW(
BCRSMatrixError,
"You can only call arithmetic operations on fully built BCRSMatrix instances");
1884 using real_type =
typename FieldTraits<ft>::real_type;
1888 real_type isNaN = 1;
1889 for (
auto const &x : *
this) {
1891 for (
auto const &y : x)
1892 sum += y.infinity_norm_real();
1893 norm = max(sum, norm);
1897 return norm * isNaN;
1903 size_type
N ()
const 1909 size_type
M ()
const 1937 #ifdef DUNE_ISTL_WITH_CHECKING 1938 if (i<0 || i>=n) DUNE_THROW(
BCRSMatrixError,
"row index out of range");
1939 if (j<0 || j>=m) DUNE_THROW(
BCRSMatrixError,
"column index out of range");
1941 return (r[i].size() && r[i].find(j) != r[i].end());
1971 std::shared_ptr<size_type>
j_;
1982 row_type current_row(a,j_.get(),0);
1983 for (size_type i=0; i<n; i++, ++row) {
1985 size_type s = row->getsize();
1989 r[i].
set(s,current_row.getptr(), current_row.getindexptr());
1991 current_row.setptr(current_row.getptr()+s);
1992 current_row.setindexptr(current_row.getindexptr()+s);
1995 r[i].
set(0,
nullptr,
nullptr);
2007 size_type* jptr = j_.get();
2008 for (size_type i=0; i<n; ++i, ++row) {
2010 size_type s = row->getsize();
2018 r[i].
set(0,
nullptr,
nullptr);
2034 for (size_type i=0; i<n; ++i) {
2036 if (r[i].getsize() > 0) {
2041 r[i].
set(0,
nullptr,
nullptr);
2052 setWindowPointers(Mat.
begin());
2055 for (size_type i=0; i<n; i++) r[i] = Mat.
r[i];
2058 build_mode = row_wise;
2073 if (allocationSize_>0)
2079 for(B *aiter=a+(allocationSize_-1), *aend=a-1; aiter!=aend; --aiter)
2080 allocator_.destroy(aiter);
2081 allocator_.deallocate(a,allocationSize_);
2088 for (size_type i=0; i<n; i++)
2089 if (r[i].getsize()>0)
2091 for (B *
col=r[i].getptr()+(r[i].getsize()-1),
2092 *colend = r[i].getptr()-1;
col!=colend; --
col) {
2093 allocator_.destroy(
col);
2095 sizeAllocator_.deallocate(r[i].getindexptr(),1);
2096 allocator_.deallocate(r[i].getptr(),1);
2099 r[i].
set(0,
nullptr,
nullptr);
2104 if (n>0 && deallocateRows && r) {
2105 for(row_type *riter=r+(n-1), *rend=r-1; riter!=rend; --riter)
2106 rowAllocator_.destroy(riter);
2107 rowAllocator_.deallocate(r,n);
2119 typename A::template rebind<size_type>::other& sizeAllocator_;
2122 Deallocator(
typename A::template rebind<size_type>::other& sizeAllocator)
2123 : sizeAllocator_(sizeAllocator)
2126 void operator()(size_type* p) { sizeAllocator_.deallocate(p,1); }
2147 void allocate(size_type rows, size_type columns, size_type allocationSize,
bool allocateRows,
bool allocate_data)
2152 nnz_ = allocationSize;
2153 allocationSize_ = allocationSize;
2159 DUNE_THROW(InvalidStateException,
"Rows have already been allocated, cannot allocate a second time");
2160 r = rowAllocator_.allocate(rows);
2169 if (allocationSize_>0) {
2172 j_.reset(sizeAllocator_.allocate(allocationSize_),
Deallocator(sizeAllocator_));
2175 for(row_type* ri=r; ri!=r+rows; ++ri)
2176 rowAllocator_.construct(ri, row_type());
2186 DUNE_THROW(InvalidStateException,
"Cannot allocate data array (already allocated)");
2187 if (allocationSize_>0) {
2188 a = allocator_.allocate(allocationSize_);
2191 new (a) B[allocationSize_];
2204 if (build_mode != implicit)
2205 DUNE_THROW(InvalidStateException,
"implicit_allocate() may only be called in implicit build mode");
2206 if (ready != notAllocated)
2207 DUNE_THROW(InvalidStateException,
"memory has already been allocated");
2210 if (overflowsize < 0)
2211 DUNE_THROW(InvalidStateException,
"You have to set the implicit build mode parameters before starting to build the matrix");
2213 size_type osize = (size_type) (_n*avg)*overflowsize + 4*
avg;
2214 allocationSize_ = _n*avg + osize;
2216 allocate(_n, _m, allocationSize_,
true,
true);
2219 size_type* jptr = j_.get() + osize;
2220 B* aptr = a + osize;
2221 for (size_type i=0; i<n; i++)
2223 r[i].
set(0,aptr,jptr);
void insert(size_type j)
put column index in row
Definition: bcrsmatrix.hh:1033
Class used by shared_ptr to deallocate memory using the proper allocator.
Definition: bcrsmatrix.hh:2117
BuildMode
we support two modes
Definition: bcrsmatrix.hh:469
Col col
Definition: matrixmatrix.hh:349
ConstIterator begin() const
Get const iterator to first row.
Definition: bcrsmatrix.hh:671
OverflowType overflow
Definition: bcrsmatrix.hh:1978
BCRSMatrix(size_type _n, size_type _m, BuildMode bm)
matrix with unknown number of nonzeroes
Definition: bcrsmatrix.hh:724
::Dune::CompressionStatistics< size_type > CompressionStatistics
The type for the statistics object returned by compress()
Definition: bcrsmatrix.hh:460
BCRSMatrix(size_type _n, size_type _m, size_type _nnz, BuildMode bm)
matrix with known number of nonzeroes
Definition: bcrsmatrix.hh:715
void mtv(const X &x, Y &y) const
y = A^T x
Definition: bcrsmatrix.hh:1659
void setImplicitBuildModeParameters(size_type _avg, double _overflow)
Set parameters needed for creation in implicit build mode.
Definition: bcrsmatrix.hh:856
A::size_type size_type
Type for indices and sizes.
Definition: matrix.hh:571
Definition: bvector.hh:1024
void umhv(const X &x, Y &y) const
y += A^H x
Definition: bcrsmatrix.hh:1729
void setsize(size_type _n)
set size only
Definition: bvector.hh:1112
FieldTraits< field_type >::real_type frobenius_norm() const
frobenius norm: sqrt(sum over squared values of entries)
Definition: bcrsmatrix.hh:1809
void mv(const X &x, Y &y) const
y = A x
Definition: bcrsmatrix.hh:1580
std::ptrdiff_t distanceTo(const RealRowIterator< ValueType > &other) const
Definition: bcrsmatrix.hh:570
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:457
size_type allocationSize_
Definition: bcrsmatrix.hh:1961
size_type M() const
The number of columns in the matrix.
Definition: bcrsmatrix.hh:210
BuildStage
Definition: bcrsmatrix.hh:426
ConstIterator beforeBegin() const
Definition: bcrsmatrix.hh:691
void setrowsize(size_type i, size_type s)
Set number of indices in row i to s.
Definition: bcrsmatrix.hh:1086
~BCRSMatrix()
destructor
Definition: bcrsmatrix.hh:791
size_type m
Definition: bcrsmatrix.hh:1959
ImplicitMatrixBuilder(Matrix &m, size_type rows, size_type cols, size_type avg_cols_per_row, double overflow_fraction)
Sets up matrix m for implicit construction using the given parameters and creates an ImplicitBmatrixu...
Definition: bcrsmatrix.hh:187
size_type nonzeroes() const
number of blocks that are stored (the number of blocks that possibly are nonzero) ...
Definition: bcrsmatrix.hh:1915
void endindices()
indicate that all indices are defined, check consistency
Definition: bcrsmatrix.hh:1217
Iterator end()
Get iterator to one beyond last row.
Definition: bcrsmatrix.hh:640
void set(size_type _n, B *_p, size_type *_j)
set size and pointer
Definition: bvector.hh:1104
Definition: matrixutils.hh:511
Iterator begin()
Get iterator to first row.
Definition: bcrsmatrix.hh:634
BCRSMatrix()
an empty matrix
Definition: bcrsmatrix.hh:708
B::field_type field_type
export the type representing the field
Definition: bcrsmatrix.hh:445
size_type n
Definition: bcrsmatrix.hh:1958
std::remove_const< T >::type ValueType
The unqualified value type.
Definition: bcrsmatrix.hh:542
void setBuildMode(BuildMode bm)
Sets the build mode of the matrix.
Definition: bcrsmatrix.hh:800
double avg
average number of non-zeroes per row.
Definition: bcrsmatrix.hh:83
void setWindowPointers(ConstRowIterator row)
Definition: bcrsmatrix.hh:1980
ConstIterator end() const
Get const iterator to one beyond last row.
Definition: bcrsmatrix.hh:677
bool equals(const RealRowIterator< ValueType > &other) const
equality
Definition: bcrsmatrix.hh:583
Definition: bcrsmatrix.hh:70
Error specific to BCRSMatrix.
Definition: istlexception.hh:19
FieldTraits< ft >::real_type infinity_norm_real() const
simplified infinity norm (uses Manhattan norm for complex values)
Definition: bcrsmatrix.hh:1837
void operator()(size_type *p)
Definition: bcrsmatrix.hh:2126
RealRowIterator< row_type > iterator
The iterator over the (mutable matrix rows.
Definition: bcrsmatrix.hh:630
double mem_ratio
fraction of wasted memory resulting from non-used overflow area.
Definition: bcrsmatrix.hh:92
T block_type
Export the type representing the components.
Definition: matrix.hh:562
B & entry(size_type row, size_type col)
Returns reference to entry (row,col) of the matrix.
Definition: bcrsmatrix.hh:1264
PropertyMapTypeSelector< Amg::VertexVisitedTag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > >::Type get(const Amg::VertexVisitedTag &tag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > &graph)
Definition: dependency.hh:292
Definition: matrixmarket.hh:258
Iterator beforeBegin()
Definition: bcrsmatrix.hh:654
block_type & operator[](size_type j) const
Returns entry in column j.
Definition: bcrsmatrix.hh:135
void copyWindowStructure(const BCRSMatrix &Mat)
Copy the window structure from another matrix.
Definition: bcrsmatrix.hh:2050
size_type avg
Definition: bcrsmatrix.hh:1974
bool contains(size_type j)
return true if column index is in row
Definition: bcrsmatrix.hh:1039
CompressedBlockVectorWindow< B, A > row_type
implement row_type with compressed vector
Definition: bcrsmatrix.hh:454
CreateIterator(BCRSMatrix &_Mat, size_type _i)
constructor
Definition: bcrsmatrix.hh:932
BuildStage ready
Definition: bcrsmatrix.hh:1948
void setindexptr(size_type *_j)
set pointer only
Definition: bvector.hh:1124
Build entries randomly with an educated guess for the number of entries per row.
Definition: bcrsmatrix.hh:498
B * a
Definition: bcrsmatrix.hh:1968
BuildStage buildStage() const
The current build stage of the matrix.
Definition: bcrsmatrix.hh:1921
Iterator RowIterator
rename the iterators for easier access
Definition: bcrsmatrix.hh:660
Matrix::block_type block_type
The block_type of the underlying matrix.
Definition: bcrsmatrix.hh:118
Matrix::size_type size_type
The size_type of the underlying matrix.
Definition: bcrsmatrix.hh:121
void usmtv(const field_type &alpha, const X &x, Y &y) const
y += alpha A^T x
Definition: bcrsmatrix.hh:1710
CreateIterator createend()
get create iterator pointing to one after the last block
Definition: bcrsmatrix.hh:1072
FieldTraits< ft >::real_type infinity_norm() const
infinity norm (row sum norm, how to generalize for blocks?)
Definition: bcrsmatrix.hh:1817
BuildMode build_mode
Definition: bcrsmatrix.hh:1947
size_type getsize() const
get size
Definition: bvector.hh:1153
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:422
void usmv(F &&alpha, const X &x, Y &y) const
y += alpha A x
Definition: bcrsmatrix.hh:1640
std::map< std::pair< size_type, size_type >, B > OverflowType
Definition: bcrsmatrix.hh:1977
void mmv(const X &x, Y &y) const
y -= A x
Definition: bcrsmatrix.hh:1621
BCRSMatrix & axpy(field_type alpha, const BCRSMatrix &b)
Add the scaled entries of another matrix to this one.
Definition: bcrsmatrix.hh:1560
void setptr(B *_p)
set pointer only
Definition: bvector.hh:1118
Statistics about compression achieved in implicit mode.
Definition: bcrsmatrix.hh:80
B block_type
export the type representing the components
Definition: bcrsmatrix.hh:448
void setSize(size_type rows, size_type columns, size_type nnz=0)
Set the size of the matrix.
Definition: bcrsmatrix.hh:828
void usmhv(const field_type &alpha, const X &x, Y &y) const
y += alpha A^H x
Definition: bcrsmatrix.hh:1767
void setIndices(size_type row, It begin, It end)
Set all column indices for row from the given iterator range.
Definition: bcrsmatrix.hh:1203
Definition: basearray.hh:19
void endrowsizes()
indicate that size of all rows is defined
Definition: bcrsmatrix.hh:1118
RealRowIterator()
empty constructor, use with care!
Definition: bcrsmatrix.hh:555
BCRSMatrix(const BCRSMatrix &Mat)
copy constructor
Definition: bcrsmatrix.hh:764
A allocator_type
export the allocator type
Definition: bcrsmatrix.hh:451
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Iterator class for sequential creation of blocks
Definition: bcrsmatrix.hh:928
size_type overflow_total
total number of elements written to the overflow area during construction.
Definition: bcrsmatrix.hh:87
The overflow error used during implicit BCRSMatrix construction was exhausted.
Definition: istlexception.hh:32
void setDataPointers()
Set data pointers for all rows.
Definition: bcrsmatrix.hh:2031
Build in a row-wise manner.
Definition: bcrsmatrix.hh:480
M_ Matrix
The underlying matrix.
Definition: bcrsmatrix.hh:115
ImplicitMatrixBuilder(Matrix &m)
Creates an ImplicitMatrixBuilder for matrix m.
Definition: bcrsmatrix.hh:163
RealRowIterator(const RealRowIterator< ValueType > &it)
Definition: bcrsmatrix.hh:559
size_type N() const
The number of rows in the matrix.
Definition: bcrsmatrix.hh:204
bool exists(size_type i, size_type j) const
return true if (i,j) is in pattern
Definition: bcrsmatrix.hh:1935
A::template rebind< size_type >::other sizeAllocator_
Definition: bcrsmatrix.hh:1955
void mmtv(const X &x, Y &y) const
y -= A^T x
Definition: bcrsmatrix.hh:1693
Iterator beforeEnd()
Definition: bcrsmatrix.hh:647
void implicit_allocate(size_type _n, size_type _m)
organizes allocation implicit mode calculates correct array size to be allocated and sets the the win...
Definition: bcrsmatrix.hh:2202
void umv(const X &x, Y &y) const
y += A x
Definition: bcrsmatrix.hh:1602
BCRSMatrix(size_type _n, size_type _m, size_type _avg, double _overflowsize, BuildMode bm)
construct matrix with a known average number of entries per row
Definition: bcrsmatrix.hh:743
size_type index() const
return index
Definition: bcrsmatrix.hh:565
RealRowIterator(row_type *_p, size_type _i)
constructor
Definition: bcrsmatrix.hh:550
Build entries randomly.
Definition: bcrsmatrix.hh:489
std::shared_ptr< size_type > j_
Definition: bcrsmatrix.hh:1971
size_type N() const
number of rows (counted in blocks)
Definition: bcrsmatrix.hh:1903
row_type * r
Definition: bcrsmatrix.hh:1965
double overflowsize
Definition: bcrsmatrix.hh:1975
size_type index() const
The number of the row that the iterator currently points to.
Definition: bcrsmatrix.hh:1027
size_type maximum
maximum number of non-zeroes per row.
Definition: bcrsmatrix.hh:85
RealRowIterator< const row_type > const_iterator
The const iterator over the matrix rows.
Definition: bcrsmatrix.hh:666
ConstIterator ConstRowIterator
rename the const row iterator for easier access
Definition: bcrsmatrix.hh:697
ConstIterator beforeEnd() const
Definition: bcrsmatrix.hh:684
row_object operator[](size_type i) const
Returns a proxy for entries in row i.
Definition: bcrsmatrix.hh:198
size_type size() const
Get the current row size.
Definition: bcrsmatrix.hh:1048
void setColumnPointers(ConstRowIterator row)
Copy row sizes from iterator range starting at row and set column index pointers for all rows...
Definition: bcrsmatrix.hh:2005
void deallocate(bool deallocateRows=true)
deallocate memory of the matrix.
Definition: bcrsmatrix.hh:2067
row_type::Iterator ColIterator
Iterator for the entries of each row.
Definition: bcrsmatrix.hh:663
RealRowIterator< row_type > Iterator
Definition: bcrsmatrix.hh:631
Proxy row object for entry access.
Definition: bcrsmatrix.hh:129
CreateIterator createbegin()
get initial create iterator
Definition: bcrsmatrix.hh:1066
A::template rebind< row_type >::other rowAllocator_
Definition: bcrsmatrix.hh:1953
Deallocator(typename A::template rebind< size_type >::other &sizeAllocator)
Definition: bcrsmatrix.hh:2122
RealRowIterator< const row_type > ConstIterator
Definition: bcrsmatrix.hh:667
std::ptrdiff_t distanceTo(const RealRowIterator< const ValueType > &other) const
Definition: bcrsmatrix.hh:576
size_type M() const
number of columns (counted in blocks)
Definition: bcrsmatrix.hh:1909
Some handy generic functions for ISTL matrices.
void umtv(const X &x, Y &y) const
y += A^T x
Definition: bcrsmatrix.hh:1674
A::template rebind< B >::other allocator_
Definition: bcrsmatrix.hh:1951
iterator class for sequential access
Definition: basearray.hh:563
size_type getrowsize(size_type i) const
get current number of indices in row i
Definition: bcrsmatrix.hh:1097
FieldTraits< field_type >::real_type frobenius_norm2() const
square of frobenius norm, need for block recursion
Definition: bcrsmatrix.hh:1788
void incrementrowsize(size_type i, size_type s=1)
increment size of row i by s (1 by default)
Definition: bcrsmatrix.hh:1107
Iterator access to matrix rows
Definition: bcrsmatrix.hh:536
A wrapper for uniform access to the BCRSMatrix during and after the build stage in implicit build mod...
Definition: bcrsmatrix.hh:109
void allocateData()
Definition: bcrsmatrix.hh:2183
CreateIterator & operator++()
prefix increment
Definition: bcrsmatrix.hh:946
bool equals(const RealRowIterator< const ValueType > &other) const
equality
Definition: bcrsmatrix.hh:590
void mmhv(const X &x, Y &y) const
y -= A^H x
Definition: bcrsmatrix.hh:1748
size_type nnz_
Definition: bcrsmatrix.hh:1960
BuildMode buildMode() const
The currently selected build mode of the matrix.
Definition: bcrsmatrix.hh:1927
row_type::ConstIterator ConstColIterator
Const iterator to the entries of a row.
Definition: bcrsmatrix.hh:700
void allocate(size_type rows, size_type columns, size_type allocationSize, bool allocateRows, bool allocate_data)
Allocate memory for the matrix structure.
Definition: bcrsmatrix.hh:2147
void addindex(size_type row, size_type col)
add index (row,col) to the matrix
Definition: bcrsmatrix.hh:1160
CompressionStatistics compress()
Finishes the buildstage in implicit mode.
Definition: bcrsmatrix.hh:1328