Skip to content
Snippets Groups Projects
Commit 503534d3 authored by Oliver Sander's avatar Oliver Sander
Browse files

moving the solver stuff into my application. I seem to be the only one who uses it

[[Imported from SVN: r4462]]
parent 01468954
No related branches found
No related tags found
No related merge requests found
# $Id$
commondir = $(includedir)/dune/solver/common
common_HEADERS = iterativesolver.cc \
iterationstep.hh numproc.hh solver.hh
include $(top_srcdir)/am/global-rules
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_ITERATIONSTEP_HH
#define DUNE_ITERATIONSTEP_HH
#include <dune/solver/common/numproc.hh>
#include <vector>
namespace Dune {
//! Base class for iteration steps being called by a linear solver
template<class OperatorType, class DiscFuncType>
class IterationStep : public NumProc
{
public:
//! Default constructor
IterationStep() {}
/** \brief Destructor */
virtual ~IterationStep() {}
//! Constructor being given linear operator, solution and right hand side
IterationStep(const OperatorType& mat, DiscFuncType& x, DiscFuncType& rhs) {
mat_ = &mat;
x_ = &x;
rhs_ = &rhs;
}
//! Set linear operator, solution and right hand side
virtual void setProblem(const OperatorType& mat, DiscFuncType& x, DiscFuncType& rhs) {
x_ = &x;
rhs_ = &rhs;
mat_ = &mat;
}
//! Do the actual iteration
virtual void iterate() = 0;
//! Return solution object
virtual DiscFuncType getSol() = 0;
//! Return linear operator
virtual const OperatorType* getMatrix() {return mat_;}
/** \brief Checks whether all relevant member variables are set
* \exception SolverError if the iteration step is not set up properly
*/
virtual void check() const {
#if 0
if (!x_)
DUNE_THROW(SolverError, "Iteration step has no solution vector");
if (!rhs_)
DUNE_THROW(SolverError, "Iteration step has no right hand side");
if (!mat_)
DUNE_THROW(SolverError, "Iteration step has no matrix");
#endif
}
//! The solution container
DiscFuncType* x_;
//! The container for the right hand side
DiscFuncType* rhs_;
//! The linear operator
const OperatorType* mat_;
/** \brief A flag for each degree of freedom stating whether the
* dof is dirichlet or not */
const std::vector<bool>* dirichletNodes_;
};
}
#endif
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
template <class OperatorType, class DiscFuncType>
void Dune::IterativeSolver<OperatorType, DiscFuncType>::check() const
{
if (!errorNorm_)
DUNE_THROW(SolverError, "You need to supply a norm-computing routine to an iterative solver!");
if (!iterationStep)
DUNE_THROW(SolverError, "You need to supply an iteration step to an iterative solver!");
iterationStep->check();
}
template <class OperatorType, class DiscFuncType>
inline
void Dune::IterativeSolver<OperatorType, DiscFuncType>::solve()
{
int i;
// Check whether the solver is set up properly
check();
if (verbosity_ != QUIET)
std::cout << "--- LinearSolver ---\n";
double error = std::numeric_limits<double>::max();
double normOfOldCorrection = 0;
// Loop until desired tolerance or maximum number of iterations is reached
for (i=0; i<numIt && (error>tolerance_ || isnan(error)); i++) {
// Backup of the current solution for the error computation later on
DiscFuncType oldSolution = iterationStep->getSol();
// Perform one iteration step
iterationStep->iterate();
// Compute error
double oldNorm = errorNorm_->compute(oldSolution);
oldSolution -= iterationStep->getSol();
double normOfCorrection = errorNorm_->compute(oldSolution);
error = normOfCorrection / oldNorm;
double convRate = normOfCorrection / normOfOldCorrection;
normOfOldCorrection = normOfCorrection;
// Output
if (verbosity_ != QUIET) {
std::cout << i << " -- ||u^{n+1} - u^n||: " << error << ", "
<< "convrate " << convRate << "\n";
}
}
if (verbosity_ != QUIET) {
std::cout << i << " iterations performed\n";
std::cout << "--------------------\n";
}
}
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_NUMPROC_HH
#define DUNE_NUMPROC_HH
namespace Dune {
/** \brief Exception thrown by solvers */
class SolverError : public Exception {};
/** \brief Base class for numerical procedures */
class NumProc
{
public:
NumProc() : verbosity_(FULL) {}
/** \brief Different levels of verbosity */
enum VerbosityMode {QUIET, REDUCED, FULL};
/** \brief Controls the amount of screen output of a numproc */
VerbosityMode verbosity_;
};
} // namespace Dune
#endif
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_SOLVER_HH
#define DUNE_SOLVER_HH
#include <dune/solver/common/numproc.hh>
namespace Dune {
/** \brief The base class for all sorts of solver algorithms */
class Solver : public NumProc
{
public:
/** \brief Virtual destructor */
virtual ~Solver() {}
/** \brief Do the necessary preprocessing */
virtual void preprocess();
/** \brief Derived classes overload this with the actual
* solution algorithm */
virtual void solve() = 0;
/** \brief The requested tolerance of the solver */
double tolerance_;
};
}
void Dune::Solver::preprocess()
{
// Default: Do nothing
}
#endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment