From b4807eba25011ba8abf453c1cfc79dc766efb03a Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Mon, 9 Dec 2013 14:13:32 +0000
Subject: [PATCH] Replace hard-coded 'double' by 'field_type'
[[Imported from SVN: r12436]]
---
dune/solvers/solvers/tcgsolver.cc | 24 ++++++++++++------------
dune/solvers/solvers/tcgsolver.hh | 2 +-
2 files changed, 13 insertions(+), 13 deletions(-)
diff --git a/dune/solvers/solvers/tcgsolver.cc b/dune/solvers/solvers/tcgsolver.cc
index c46f0b6c..787dfc0b 100644
--- a/dune/solvers/solvers/tcgsolver.cc
+++ b/dune/solvers/solvers/tcgsolver.cc
@@ -53,10 +53,10 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
std::cout << std::endl;
}
- double error = std::numeric_limits<double>::max();
+ field_type error = std::numeric_limits<field_type>::max();
- double normOfOldCorrection = 0;
- double totalConvRate = 1;
+ field_type normOfOldCorrection = 0;
+ field_type totalConvRate = 1;
this->maxTotalConvRate_ = 0;
int convRateCounter = 0;
@@ -67,7 +67,7 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
VectorType q(*x_); // a temporary vector
// some local variables
- double rho,rholast,beta;
+ field_type rho,rholast,beta;
// determine initial search direction
p = 0; // clear correction
@@ -82,9 +82,9 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
rholast = p*b; // orthogonalization
- double solutionNormSquared = 0;
- double solutionTimesCorrection = 0;
- double correctionNormSquared = rholast;
+ field_type solutionNormSquared = 0;
+ field_type solutionTimesCorrection = 0;
+ field_type correctionNormSquared = rholast;
// Loop until desired tolerance or maximum number of iterations is reached
for (i=0; i<this->maxIterations_ && (error>this->tolerance_ || std::isnan(error)); i++) {
@@ -147,7 +147,7 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
// minimize in given search direction p
matrix_->mv(p,q); // q=Ap
- double alpha = rholast/(p*q); // scalar product
+ field_type alpha = rholast/(p*q); // scalar product
// ////////////////////////////////////////////////////////////////////////
// Compute the new iterate. If it is outside of the admissible set
@@ -159,7 +159,7 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
tentativeNewIterate.axpy(alpha,p); // update solution
// Compute length of this tentative new iterate
- double tentativeLengthSquared = solutionNormSquared
+ field_type tentativeLengthSquared = solutionNormSquared
+ 2*alpha*solutionTimesCorrection
+ alpha*alpha*correctionNormSquared;
@@ -221,17 +221,17 @@ void TruncatedCGSolver<MatrixType, VectorType>::solve()
// //////////////////////////////////////////////////
// Compute error
// //////////////////////////////////////////////////
- double oldNorm = errorNorm_->operator()(oldSolution);
+ field_type oldNorm = errorNorm_->operator()(oldSolution);
oldSolution -= *x_;
- double normOfCorrection = errorNorm_->operator()(oldSolution);
+ field_type normOfCorrection = errorNorm_->operator()(oldSolution);
if (this->useRelativeError_)
error = normOfCorrection / oldNorm;
else
error = normOfCorrection;
- double convRate = normOfCorrection / normOfOldCorrection;
+ field_type convRate = normOfCorrection / normOfOldCorrection;
normOfOldCorrection = normOfCorrection;
diff --git a/dune/solvers/solvers/tcgsolver.hh b/dune/solvers/solvers/tcgsolver.hh
index a050eb78..da48e8e8 100644
--- a/dune/solvers/solvers/tcgsolver.hh
+++ b/dune/solvers/solvers/tcgsolver.hh
@@ -86,7 +86,7 @@ public:
void setProblem(const MatrixType& matrix,
VectorType* x,
VectorType* rhs,
- double trustRegionRadius)
+ field_type trustRegionRadius)
{
matrix_ = &matrix;
x_ = x;
--
GitLab