From 2c730022b44a3a183faaea1429cd69bfcfc20af1 Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Tue, 12 Jul 2016 20:27:22 +0200
Subject: [PATCH] Remove ::CGSolver
Please use Dune::Solvers::CGStep in conjunction with ::LoopSolver
instead.
---
dune/solvers/solvers/CMakeLists.txt | 2 -
dune/solvers/solvers/Makefile.am | 2 +-
dune/solvers/solvers/cgsolver.cc | 166 ----------------------------
dune/solvers/solvers/cgsolver.hh | 61 ----------
dune/solvers/test/cgsteptest.cc | 39 -------
5 files changed, 1 insertion(+), 269 deletions(-)
delete mode 100644 dune/solvers/solvers/cgsolver.cc
delete mode 100644 dune/solvers/solvers/cgsolver.hh
diff --git a/dune/solvers/solvers/CMakeLists.txt b/dune/solvers/solvers/CMakeLists.txt
index 2f2d709b..652164b1 100644
--- a/dune/solvers/solvers/CMakeLists.txt
+++ b/dune/solvers/solvers/CMakeLists.txt
@@ -1,6 +1,4 @@
install(FILES
- cgsolver.cc
- cgsolver.hh
criterion.hh
iterativesolver.cc
iterativesolver.hh
diff --git a/dune/solvers/solvers/Makefile.am b/dune/solvers/solvers/Makefile.am
index a96e9e53..5c9d0e02 100644
--- a/dune/solvers/solvers/Makefile.am
+++ b/dune/solvers/solvers/Makefile.am
@@ -1,7 +1,7 @@
SUBDIRS =
solversdir = $(includedir)/dune/solvers/solvers
-solvers_HEADERS = cgsolver.cc cgsolver.hh criterion.hh iterativesolver.cc iterativesolver.hh \
+solvers_HEADERS = criterion.hh iterativesolver.cc iterativesolver.hh \
loopsolver.cc loopsolver.hh quadraticipopt.hh solver.hh \
tcgsolver.cc tcgsolver.hh trustregionsolver.cc trustregionsolver.hh
diff --git a/dune/solvers/solvers/cgsolver.cc b/dune/solvers/solvers/cgsolver.cc
deleted file mode 100644
index d6cc5f5c..00000000
--- a/dune/solvers/solvers/cgsolver.cc
+++ /dev/null
@@ -1,166 +0,0 @@
-// -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=8 sw=4 sts=4:
-
-#include <cmath>
-#include <limits>
-#include <iomanip>
-
-template <class MatrixType, class VectorType>
-void CGSolver<MatrixType, VectorType>::check() const
-{
- if (!errorNorm_)
- DUNE_THROW(SolverError, "You need to supply a norm to a CG solver!");
-
-}
-
-template <class MatrixType, class VectorType>
-void CGSolver<MatrixType, VectorType>::solve()
-{
-
- int i;
-
- VectorType& b = *rhs_;
-
- // Check whether the solver is set up properly
- check();
-
- // set up preconditioner
- if (preconditioner_)
- preconditioner_->preprocess();
-
- if (this->verbosity_ != NumProc::QUIET)
- std::cout << "--- CGSolver ---\n";
-
- if (this->verbosity_ == NumProc::FULL)
- {
- std::cout << " iter";
- if (this->useRelativeError_)
- std::cout << " correction ";
- else
- std::cout << " abs correction ";
- std::cout << " rate";
- std::string header;
- if (preconditioner_) {
- header = preconditioner_->getOutput();
- std::cout << header;
- }
- std::cout << std::endl;
- std::cout << "-----";
- std::cout << "---------------";
- std::cout << "---------";
- for(size_t i=0; i<header.size(); ++i)
- std::cout << "-";
- std::cout << std::endl;
- }
-
- double error = std::numeric_limits<double>::max();
-
- double normOfOldCorrection = 0;
- double totalConvRate = 1;
- this->maxTotalConvRate_ = 0;
- int convRateCounter = 0;
-
- // overwrite b with residual
- matrix_->mmv(*x_,b);
-
- VectorType p(*x_); // the search direction
- VectorType q(*x_); // a temporary vector
-
- // some local variables
- double rho,rholast,lambda,alpha,beta;
-
- // determine initial search direction
- p = 0; // clear correction
-
- if (preconditioner_) {
- preconditioner_->setProblem(*matrix_,p,b);
- preconditioner_->iterate();
- p = preconditioner_->getSol();
- } else
- p = b;
-
- rholast = p*b; // orthogonalization
-
- // Loop until desired tolerance or maximum number of iterations is reached
- for (i=0; i<this->maxIterations_ && (error>this->tolerance_ || std::isnan(error)); i++) {
-
- // Backup of the current solution for the error computation later on
- VectorType oldSolution = *x_;
-
- // Perform one iteration step
- // minimize in given search direction p
- matrix_->mv(p,q); // q=Ap
- alpha = p*q; // scalar product
- lambda = rholast/alpha; // minimization
- x_->axpy(lambda,p); // update solution
- b.axpy(-lambda,q); // update residual
-
- // determine new search direction
- q = 0; // clear correction
-
- // apply preconditioner
- if (preconditioner_) {
- preconditioner_->setProblem(*matrix_,q,b);
- preconditioner_->iterate();
- q = preconditioner_->getSol();
- } else
- q = b;
-
- rho = q*b; // orthogonalization
- beta = rho/rholast; // scaling factor
- p *= beta; // scale old search direction
- p += q; // orthogonalization with correction
- rholast = rho; // remember rho for recurrence
-
- // write iteration to file, if requested
- if (this->historyBuffer_!="")
- this->writeIterate(*x_, i);
-
- // Compute error
- double oldNorm = errorNorm_->operator()(oldSolution);
- oldSolution -= *x_;
-
- double normOfCorrection = errorNorm_->operator()(oldSolution);
-
- if (this->useRelativeError_)
- error = normOfCorrection / oldNorm;
- else
- error = normOfCorrection;
-
- double convRate = normOfCorrection / normOfOldCorrection;
-
- normOfOldCorrection = normOfCorrection;
-
- if (!std::isinf(convRate) && !std::isnan(convRate)) {
- totalConvRate *= convRate;
- this->maxTotalConvRate_ = std::max(this->maxTotalConvRate_, std::pow(totalConvRate, 1/((double)convRateCounter+1)));
- convRateCounter++;
- }
-
- // Output
- if (this->verbosity_ == NumProc::FULL) {
- std::streamsize const oldPrecision = std::cout.precision();
- std::ios_base::fmtflags const oldFormatFlags = std::cout.flags();
-
- std::cout << std::setw(5) << i
-
- << std::scientific
- << std::setw(15) << std::setprecision(7) << error
-
- << std::fixed
- << std::setw(9) << std::setprecision(5) << convRate
-
- << std::endl;
- std::cout << std::setprecision(oldPrecision);
- std::cout.flags(oldFormatFlags);
- }
- }
-
-
- if (this->verbosity_ != NumProc::QUIET) {
- std::cout << "maxTotalConvRate: " << this->maxTotalConvRate_ << ", "
- << i << " iterations performed\n";
- std::cout << "--------------------\n";
- }
-
-}
diff --git a/dune/solvers/solvers/cgsolver.hh b/dune/solvers/solvers/cgsolver.hh
deleted file mode 100644
index 99561566..00000000
--- a/dune/solvers/solvers/cgsolver.hh
+++ /dev/null
@@ -1,61 +0,0 @@
-// -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=8 sw=4 sts=4:
-#ifndef CG_SOLVER_HH
-#define CG_SOLVER_HH
-
-#include <dune/solvers/solvers/iterativesolver.hh>
-#include <dune/solvers/iterationsteps/lineariterationstep.hh>
-#include <dune/solvers/norms/norm.hh>
-
-/** \brief A conjugate gradient solver
- *
- */
-template <class MatrixType, class VectorType>
-class CGSolver : public IterativeSolver<VectorType>
-{
- static const int blocksize = VectorType::block_type::dimension;
-
-public:
-
- /** \brief Constructor taking all relevant data */
-
- CGSolver(const MatrixType* matrix,
- VectorType* x,
- VectorType* rhs,
- LinearIterationStep<MatrixType,VectorType>* preconditioner,
- int maxIterations,
- double tolerance,
- Norm<VectorType>* errorNorm,
- NumProc::VerbosityMode verbosity,
- bool useRelativeError=true)
- : IterativeSolver<VectorType>(tolerance,maxIterations,verbosity,useRelativeError),
- matrix_(matrix), x_(x), rhs_(rhs),
- preconditioner_(preconditioner), errorNorm_(errorNorm)
- {}
-
- /** \brief Loop, call the iteration procedure
- * and monitor convergence */
- virtual void solve();
-
- /** \brief Checks whether all relevant member variables are set
- * \exception SolverError if the iteration step is not set up properly
- */
- virtual void check() const;
-
- const MatrixType* matrix_;
-
- VectorType* x_;
-
- VectorType* rhs_;
-
- //! The iteration step used by the algorithm
- LinearIterationStep<MatrixType,VectorType>* preconditioner_;
-
- //! The norm used to measure convergence
- Norm<VectorType>* errorNorm_;
-
-};
-
-#include "cgsolver.cc"
-
-#endif
diff --git a/dune/solvers/test/cgsteptest.cc b/dune/solvers/test/cgsteptest.cc
index 75688f53..8afc233b 100644
--- a/dune/solvers/test/cgsteptest.cc
+++ b/dune/solvers/test/cgsteptest.cc
@@ -11,8 +11,6 @@
#include <dune/solvers/solvers/loopsolver.hh>
#include <dune/solvers/iterationsteps/blockgssteps.hh>
#include <dune/solvers/iterationsteps/cgstep.hh>
-#include <dune/solvers/solvers/cgsolver.hh>
-
#include "common.hh"
@@ -66,43 +64,6 @@ struct CGTestSuite
const bool relativeErrors = false;
Analyser<Problem> analyser(p, solveTol);
- {
- std::string const testCase = "CGSolver ";
-
- typename Problem::Vector u_copy = p.u;
- typename Problem::Vector rhs_copy = p.rhs;
-
- CGSolver<typename Problem::Matrix, typename Problem::Vector> solver(
- &p.A, &u_copy, &rhs_copy, nullptr, maxIterations, stepTol,
- &p.energyNorm, verbosity, relativeErrors);
- solver.check();
- solver.preprocess();
- solver.solve();
-
- passed &= analyser.analyse(u_copy, testCase);
- }
- {
- std::string const testCase = "CGSolver, preconditioned ";
-
- typename Problem::Vector u_copy = p.u;
- typename Problem::Vector rhs_copy = p.rhs;
-
- auto blockgs =
- Dune::Solvers::BlockGSStepFactory<typename Problem::Matrix,
- typename Problem::Vector>::
- create(Dune::Solvers::BlockGS::LocalSolvers::gs(),
- Dune::Solvers::BlockGS::Direction::SYMMETRIC);
- blockgs.setIgnore(p.ignore);
-
- CGSolver<typename Problem::Matrix, typename Problem::Vector> solver(
- &p.A, &u_copy, &rhs_copy, &blockgs, maxIterations, stepTol,
- &p.energyNorm, verbosity, relativeErrors);
- solver.check();
- solver.preprocess();
- solver.solve();
-
- passed &= analyser.analyse(u_copy, testCase);
- }
{
std::string const testCase = "Dune::Solvers::CGStep ";
--
GitLab