From 62e628c811498fd7f4509640e854aa7f46b3bf30 Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Sun, 10 May 2009 20:38:33 +0000
Subject: [PATCH] minor fixes concerning trust-region with matrix-induced norms

[[Imported from SVN: r2427]]
---
 dune-solvers/solvers/tcgsolver.cc | 6 +++---
 dune-solvers/solvers/tcgsolver.hh | 5 ++++-
 2 files changed, 7 insertions(+), 4 deletions(-)

diff --git a/dune-solvers/solvers/tcgsolver.cc b/dune-solvers/solvers/tcgsolver.cc
index 511d1016..e67c3af7 100644
--- a/dune-solvers/solvers/tcgsolver.cc
+++ b/dune-solvers/solvers/tcgsolver.cc
@@ -69,7 +69,7 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
 
             *x_ = eta;  // return eta
 
-            std::cout << "Abort1, radius = " << std::sqrt((*x_)*(*x_)) << std::endl;
+            std::cout << "Abort1, radius = " << std::sqrt(trustRegionScalarProduct(*x_,*x_)) << std::endl;
 
             if (this->historyBuffer_!="")
                 this->writeIterate(*x_, i);
@@ -92,7 +92,7 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
             *x_ = eta_j;  // return eta + tau*delta_j
             x_->axpy(tau, delta_j);
 
-            std::cout << "Abort2, radius = " << std::sqrt((*x_)*(*x_)) << std::endl;
+            std::cout << "Abort2, radius = " << std::sqrt(trustRegionScalarProduct(*x_,*x_)) << std::endl;
 
             if (this->historyBuffer_!="")
                 this->writeIterate(*x_, i);
@@ -113,7 +113,7 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
         *x_ = eta_jPlus1;
         eta_j = eta_jPlus1;
 
-        std::cout << "Radius = " << sqrt((*x_)*(*x_)) << std::endl;
+        std::cout << "Radius = " << sqrt(trustRegionScalarProduct(*x_,*x_)) << std::endl;
 
         // /////////////////////////////////////////////
         //   write iteration to file, if requested
diff --git a/dune-solvers/solvers/tcgsolver.hh b/dune-solvers/solvers/tcgsolver.hh
index 7144e667..e0fea8e2 100644
--- a/dune-solvers/solvers/tcgsolver.hh
+++ b/dune-solvers/solvers/tcgsolver.hh
@@ -28,6 +28,8 @@ class TruncatedCGSolver : public IterativeSolver<VectorType>
         // There must be one positive and one negative root
         assert(root1*root2 <= 0);
 
+        std::cout << "root1: " << root1 << ",  root2: " << root2 << std::endl;
+
         return std::max(root1, root2);
     }
 
@@ -71,13 +73,14 @@ public:
     TruncatedCGSolver(int maxIterations,
                       double tolerance,
                       Norm<VectorType>* errorNorm,
+                      const MatrixType* trustRegionNormMatrix,
                       NumProc::VerbosityMode verbosity,
                       bool useRelativeError=true)
         : IterativeSolver<VectorType>(tolerance,maxIterations,verbosity,useRelativeError),
           matrix_(NULL), x_(NULL), rhs_(NULL),
           errorNorm_(errorNorm),
           trustRegionRadius_(0),
-          trustRegionNormMatrix_(NULL)
+          trustRegionNormMatrix_(trustRegionNormMatrix)
     {}
     
     void setProblem(const MatrixType& matrix,
-- 
GitLab