Merge branch 'feature/resize-in-matrix-vector' into 'master'

Feature/resize in matrix vector

See merge request !19

dune/solvers/common/genericvectortools.hh

@@ -22,6 +22,7 @@
 #include <dune/istl/bvector.hh>
 
 #include <dune/matrix-vector/genericvectortools.hh>
+#include <dune/matrix-vector/resize.hh>
 
 #include <dune/solvers/operators/sumoperator.hh>
 
@@ -49,44 +50,12 @@ struct GenericVector
 
     //! Resize vector recursivly to size of given vector/matrix
     template <class VectorTypeA, class VectorTypeB>
+    DUNE_DEPRECATED_MSG("Please use Dune::MatrixVector::resize instead.")
     static void resize(VectorTypeA& a, const VectorTypeB& b)
     {
-        a.resize(b.size());
-        typename VectorTypeB::const_iterator it = b.begin();
-        typename VectorTypeB::const_iterator end = b.end();
-        for(; it!=end; ++it)
-            GenericVector::resize(a[it.index()], *it);
+        Dune::MatrixVector::resize(a, b);
     }
 
-    template <class VectorTypeA, class T, class A>
-    static void resize(VectorTypeA& a, const Dune::BCRSMatrix<T,A>& b)
-    {
-        a.resize(b.N());
-        typedef typename Dune::BCRSMatrix<T,A> VectorTypeB;
-        typename VectorTypeB::const_iterator it = b.begin();
-        typename VectorTypeB::const_iterator end = b.end();
-        for(; it!=end; ++it)
-            GenericVector::resize(a[it.index()], *(it->begin()));
-    }
-
-    template <class VectorTypeA, class SparseMatrixType, class LowRankMatrixType>
-    static void resize(VectorTypeA& a, const SumOperator<SparseMatrixType, LowRankMatrixType>& b)
-    {
-        a.resize(b.N());
-        typename SparseMatrixType::const_iterator it  = b.sparseMatrix().begin();
-        typename SparseMatrixType::const_iterator end = b.sparseMatrix().end();
-        for(; it!=end; ++it)
-            GenericVector::resize(a[it.index()], *(it->begin()));
-    }
-
-    template <class field_type, int n, class VectorTypeB>
-    static void resize(Dune::FieldVector<field_type,n>& a DUNE_UNUSED, const VectorTypeB& b DUNE_UNUSED)
-    {}
-
-    template <int n, class VectorTypeB>
-    static void resize(std::bitset<n>& a DUNE_UNUSED, const VectorTypeB& b DUNE_UNUSED)
-    {}
-
     //! Set vector to zero at indices that are true in bitvector recursivly
     template <class VectorType, class BitVectorType>
     DUNE_DEPRECATED_MSG("Please use Dune::MatrixVector::Generic::truncate instead.")

dune/solvers/operators/sumoperator.hh

@@ -4,6 +4,7 @@
 #define SUMOPERATOR_HH
 
 #include <cstddef>
+#include <dune/matrix-vector/resize.hh>
 
 /** \brief represents the sum of two linear operators
  *
@@ -139,5 +140,14 @@ class SumOperator
         const bool summands_allocated_internally_;
 };
 
+//! inject the resize from SumOperator to generic vector tools
+namespace Dune { namespace MatrixVector {
+template <class Vector, class SparseMatrix, class LowRankMatrix>
+//! Resize vector from SumOperator. Note that we only consider the sparse part.
+void resize(Vector& v, const SumOperator<SparseMatrix, LowRankMatrix>& sumOp) {
+  resize(v, sumOp.sparseMatrix());
+}
+}}
+
 #endif
 

dune/solvers/test/sumoperatortest.cc

@@ -170,6 +170,30 @@ bool check()
     return passed;
 }
 
+bool checkResize() try {
+  std::cout << "Checking resize from SumOperator ... ";
+  using namespace Dune;
+  using Block = FieldMatrix<double, 2, 6>;
+  using SparseMatrix = BCRSMatrix<Block>;
+  SparseMatrix sparseMatrix;
+  MatrixIndexSet indices(3,4);
+  indices.exportIdx(sparseMatrix);
+  using LowRankFactor = Matrix<Block>;
+  SumOperator<SparseMatrix, LowRankOperator<LowRankFactor>> sumOp;
+  sumOp.sparseMatrix() = sparseMatrix;
+  sumOp.lowRankMatrix().lowRankFactor() = LowRankFactor(5,6); // "weird" sizes should not matter
+  BlockVector<Dune::FieldVector<double, 2>> vector;
+  MatrixVector::resize(vector, sumOp);
+  if(vector.size() != 3)
+    DUNE_THROW(Exception, "Resize yields unexpected vector size.");
+  std::cout << "passed" << std::endl;
+  return true;
+} catch(Dune::Exception e) {
+  std::cout << " FAILED." << std::endl;
+  std::cout << "Exception was: " << e << std::endl;
+  return false;
+}
+
 int main(int argc, char** argv)
 {
     Dune::MPIHelper::instance(argc, argv);
@@ -194,5 +218,7 @@ int main(int argc, char** argv)
     passed = passed and check<Dune::BCRSMatrix<Block3>, Dune::Matrix<Block2> >();
     passed = passed and check<Dune::BCRSMatrix<Block3>, Dune::Matrix<Block3> >();
 
+    passed = passed and checkResize();
+
     return passed ? 0 : 1;
 }