#include <base/logstream.h>
#include <lac/block_sparsity_pattern.h>
#include <lac/block_sparse_matrix.h>
+#include <lac/block_vector.h>
#include <fstream>
#include <algorithm>
for (unsigned int row=0; row<19; ++row)
for (unsigned int i=0; i<10; ++i)
bsm.add (row, (row*5+i*9)%29, 0.5);
+
+ // now allocate two block vectors
+ // and see what we can get after
+ // vmults:
+ BlockVector<2,double> src;
+ vector<unsigned int> src_sizes (2);
+ src_sizes[0] = 10;
+ src_sizes[1] = 19;
+ src.reinit (src_sizes);
+
+ BlockVector<3,double> dst;
+ vector<unsigned int> dst_sizes (3);
+ dst_sizes[0] = 2;
+ dst_sizes[1] = 7;
+ dst_sizes[2] = 10;
+ dst.reinit (dst_sizes);
+
+ for (unsigned int i=0; i<29; ++i)
+ src(i) = i;
+
+ bsm.vmult (dst, src);
+ // now check what came out
+ for (unsigned int row=0; row<19; ++row)
+ {
+ vector<double> t(29, 0);
+ // first check which elements
+ // in this row exist
+ for (unsigned int i=0; i<10; ++i)
+ t[(row*5+i*9)%29] = row*((row*5+i*9)%29);
+
+ for (unsigned int i=0; i<10; ++i)
+ t[(row*5+i*9)%29] += 0.5;
+
+ // compute the exact result
+ double row_sum = 0;
+ for (unsigned int i=0; i<29; ++i)
+ row_sum += t[i]*i;
+
+ // compare to vmult result
+ Assert (row_sum == dst(row), ExcInternalError());
+ deallog << "vmult " << row << ' ' << row_sum << ' ' << dst(row) << endl;
+ };
+
+
+ // test matrix_scalar_product. note that dst=M*src
+ const double msp1 = dst.norm_sqr ();
+ const double msp2 = bsm.matrix_scalar_product (dst, src);
+ Assert (msp1 == msp2, ExcInternalError());
+ deallog << "matrix_scalar_product " << msp1 << ' ' << msp2 << endl;
};
DEAL::Row=18: expected length=11, actual length=11
DEAL::196==196
DEAL::196==196
+DEAL::vmult 0 72.000 72.000
+DEAL::vmult 1 2658.000 2658.000
+DEAL::vmult 2 6332.500 6332.500
+DEAL::vmult 3 6288.000 6288.000
+DEAL::vmult 4 10810.500 10810.500
+DEAL::vmult 5 12221.500 12221.500
+DEAL::vmult 6 18749.000 18749.000
+DEAL::vmult 7 14275.500 14275.500
+DEAL::vmult 8 20949.000 20949.000
+DEAL::vmult 9 20999.000 20999.000
+DEAL::vmult 10 29845.500 29845.500
+DEAL::vmult 11 29848.500 29848.500
+DEAL::vmult 12 31135.500 31135.500
+DEAL::vmult 13 29638.500 29638.500
+DEAL::vmult 14 40618.000 40618.000
+DEAL::vmult 15 39010.000 39010.000
+DEAL::vmult 16 39234.000 39234.000
+DEAL::vmult 17 51127.500 51127.500
+DEAL::vmult 18 51714.500 51714.500
+DEAL::matrix_scalar_product 15416941535.500 15416941535.500