From db0a17a105a216775e922122c917021d2124bcc3 Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Thu, 11 Jul 2013 15:54:02 +0000
Subject: [PATCH] Use the field_type instead of hard-coded 'double'
[[Imported from SVN: r11761]]
---
.../iterationsteps/trustregiongsstep.cc | 24 +++++++++----------
.../iterationsteps/trustregiongsstep.hh | 2 +-
2 files changed, 13 insertions(+), 13 deletions(-)
diff --git a/dune/solvers/iterationsteps/trustregiongsstep.cc b/dune/solvers/iterationsteps/trustregiongsstep.cc
index 6f2c9ce1..cd465133 100644
--- a/dune/solvers/iterationsteps/trustregiongsstep.cc
+++ b/dune/solvers/iterationsteps/trustregiongsstep.cc
@@ -7,7 +7,7 @@ DiscFuncType TrustRegionGSStep<OperatorType, DiscFuncType>::getSol()
template<class OperatorType, class DiscFuncType>
inline
-double TrustRegionGSStep<OperatorType, DiscFuncType>::
+typename DiscFuncType::field_type TrustRegionGSStep<OperatorType, DiscFuncType>::
residual(unsigned int index, int blockIndex) const
{
const OperatorType& mat = *this->mat_;
@@ -18,7 +18,7 @@ residual(unsigned int index, int blockIndex) const
/* The following loop computes
* \f[ sum_i = \sum_{j \ne i} A_{ij}x_j \f]
*/
- double r = (*this->rhs_)[index][blockIndex];
+ field_type r = (*this->rhs_)[index][blockIndex];
ColumnIterator cIt = mat[index].template begin();
ColumnIterator cEndIt = mat[index].template end();
@@ -62,23 +62,23 @@ void TrustRegionGSStep<OperatorType, DiscFuncType>::iterate()
continue;
// Compute residual
- double r = residual(i, j);
- const double& diag = mat[i][i][j][j];
+ field_type r = residual(i, j);
+ const field_type& diag = mat[i][i][j][j];
// Find local minimum
- double& x = (*this->x_)[i][j];
+ field_type& x = (*this->x_)[i][j];
if (diag > 0) {
//printf("%d %d konvex\n", i, j);
//double oldEnergy = computeEnergy(mat, *this->x_, *this->rhs_);
// 1d problem is convex
#ifndef NDEBUG
- double oldX = x;
+ field_type oldX = x;
#endif
x = r / diag;
#ifndef NDEBUG
- double energyX = 0.5*diag*x*x - r*x;
- double energyOldX = 0.5*diag*oldX*oldX - r*oldX;
+ field_type energyX = 0.5*diag*x*x - r*x;
+ field_type energyOldX = 0.5*diag*oldX*oldX - r*oldX;
assert(energyX <= energyOldX + 1e-6);
#endif
if ( x < obstacles[i].lower(j)) {
@@ -104,11 +104,11 @@ void TrustRegionGSStep<OperatorType, DiscFuncType>::iterate()
//printf("%d %d konkav\n", i, j);
// 1d problem is concave or linear.
// Minimum is attained at one of the boundaries
- double lBound = obstacles[i].lower(j);
- double uBound = obstacles[i].upper(j);
+ field_type lBound = obstacles[i].lower(j);
+ field_type uBound = obstacles[i].upper(j);
- double lValue = 0.5*diag*lBound*lBound - r*lBound;
- double uValue = 0.5*diag*uBound*uBound - r*uBound;
+ field_type lValue = 0.5*diag*lBound*lBound - r*lBound;
+ field_type uValue = 0.5*diag*uBound*uBound - r*uBound;
#if 0
double max = 0.5*diag*(r/diag)*(r/diag) - r*r/diag;
diff --git a/dune/solvers/iterationsteps/trustregiongsstep.hh b/dune/solvers/iterationsteps/trustregiongsstep.hh
index 79d142ee..b43054b0 100644
--- a/dune/solvers/iterationsteps/trustregiongsstep.hh
+++ b/dune/solvers/iterationsteps/trustregiongsstep.hh
@@ -37,7 +37,7 @@
virtual DiscFuncType getSol();
- double residual(unsigned int index, int blockIndex) const;
+ field_type residual(unsigned int index, int blockIndex) const;
};
--
GitLab