From 50ee728a2f99ff9aae63c8a089cec60a97c62648 Mon Sep 17 00:00:00 2001
From: Uli Sack <usack@math.fu-berlin.de>
Date: Mon, 22 Oct 2012 09:08:55 +0000
Subject: [PATCH] added support for coefficients in the sumoperator, e.g. now
 implementing aA+bB

[[Imported from SVN: r7056]]
---
 dune/solvers/operators/sumoperator.hh | 32 ++++++++++++++++++++-------
 1 file changed, 24 insertions(+), 8 deletions(-)

diff --git a/dune/solvers/operators/sumoperator.hh b/dune/solvers/operators/sumoperator.hh
index ddead966..6e3d185c 100644
--- a/dune/solvers/operators/sumoperator.hh
+++ b/dune/solvers/operators/sumoperator.hh
@@ -22,15 +22,28 @@ class SumOperator
 
         //! default constructor - allocates memory for summand operators internally
         SumOperator():
+            alpha_(1.0),
+            beta_(1.0),
             summands_allocated_internally_(true)
         {
             sparse_matrix_ = new SparseMatrixType;
             lowrank_matrix_ = new LowRankMatrixType;
         }
 
+        //! construct from summands
+        SumOperator(double a, SparseMatrixType& A, double b, LowRankMatrixType& M):
+            alpha_(a),
+            sparse_matrix_(&A),
+            beta_(b),
+            lowrank_matrix_(&M),
+            summands_allocated_internally_(false)
+        {}
+
         //! construct from summands
         SumOperator(SparseMatrixType& A, LowRankMatrixType& M):
+            alpha_(1.0),
             sparse_matrix_(&A),
+            beta_(1.0),
             lowrank_matrix_(&M),
             summands_allocated_internally_(false)
         {}
@@ -48,32 +61,33 @@ class SumOperator
         template <class LVectorType, class RVectorType>
         void umv(const LVectorType& x, RVectorType& b) const
         {
-            sparse_matrix_->umv(x,b);
-            lowrank_matrix_->umv(x,b);
+            sparse_matrix_->usmv(alpha_,x,b);
+            lowrank_matrix_->usmv(beta_,x,b);
         }
 
         //! b -= (A+M)x
         template <class LVectorType, class RVectorType>
         void mmv(const LVectorType& x, RVectorType& b) const
         {
-            sparse_matrix_->mmv(x,b);
-            lowrank_matrix_->mmv(x,b);
+            sparse_matrix_->usmv(-alpha_,x,b);
+            lowrank_matrix_->usmv(-beta_,x,b);
         }
 
         //! b = (A+M)x
         template <class LVectorType, class RVectorType>
         void mv(const LVectorType& x, RVectorType& b) const
         {
-            sparse_matrix_->mv(x,b);
-            lowrank_matrix_->umv(x,b);
+            b = 0.0;
+            sparse_matrix_->usmv(alpha_,x,b);
+            lowrank_matrix_->usmv(beta_,x,b);
         }
 
         //! b += a*(A+M)x
         template <class LVectorType, class RVectorType>
         void usmv(const field_type a, const LVectorType& x, RVectorType& b) const
         {
-            sparse_matrix_->usmv(a,x,b);
-            lowrank_matrix_->usmv(a,x,b);
+            sparse_matrix_->usmv(a*alpha_,x,b);
+            lowrank_matrix_->usmv(a*beta_,x,b);
         }
 
         //! return the number of rows
@@ -112,6 +126,8 @@ class SumOperator
             return *lowrank_matrix_;
         }
 
+        double alpha_,
+               beta_;
     private:
         //! one summand
         SparseMatrixType*  sparse_matrix_;
-- 
GitLab