From c81a840382ac396821ac344f4e6dd8e9c6785111 Mon Sep 17 00:00:00 2001 From: Xavier Raynaud Date: Wed, 6 Jun 2012 15:14:46 +0200 Subject: [PATCH] Added possibility to set tolerance for linear solvers. --- opm/core/linalg/LinearSolverFactory.cpp | 9 +++++++++ opm/core/linalg/LinearSolverFactory.hpp | 11 +++++++++++ opm/core/linalg/LinearSolverInterface.hpp | 8 ++++++++ opm/core/linalg/LinearSolverIstl.cpp | 11 +++++++++++ opm/core/linalg/LinearSolverIstl.hpp | 9 +++++++++ opm/core/linalg/LinearSolverUmfpack.cpp | 10 ++++++++++ opm/core/linalg/LinearSolverUmfpack.hpp | 12 ++++++++++++ opm/core/pressure/IncompTpfa.cpp | 16 +++++++++++++--- opm/core/pressure/IncompTpfa.hpp | 13 +++++++++++-- 9 files changed, 94 insertions(+), 5 deletions(-) diff --git a/opm/core/linalg/LinearSolverFactory.cpp b/opm/core/linalg/LinearSolverFactory.cpp index 9ba1e1388..6943b88cd 100644 --- a/opm/core/linalg/LinearSolverFactory.cpp +++ b/opm/core/linalg/LinearSolverFactory.cpp @@ -93,6 +93,15 @@ namespace Opm return solver_->solve(size, nonzeros, ia, ja, sa, rhs, solution); } + void LinearSolverFactory::setTolerance(const double tol) + { + solver_->setTolerance(tol); + } + + double LinearSolverFactory::getTolerance() const + { + return solver_->getTolerance(); + } diff --git a/opm/core/linalg/LinearSolverFactory.hpp b/opm/core/linalg/LinearSolverFactory.hpp index ad5377cf3..e6349c9e2 100644 --- a/opm/core/linalg/LinearSolverFactory.hpp +++ b/opm/core/linalg/LinearSolverFactory.hpp @@ -75,6 +75,17 @@ namespace Opm const double* sa, const double* rhs, double* solution) const; + + /// Set tolerance for the linear solver. + /// \param[in] tol tolerance value + /// Not used for LinearSolverFactory + virtual void setTolerance(const double tol); + + /// Get tolerance for the linear solver. + /// \param[out] tolerance value + /// Not used for LinearSolverFactory. Returns -1. + virtual double getTolerance() const; + private: std::tr1::shared_ptr solver_; }; diff --git a/opm/core/linalg/LinearSolverInterface.hpp b/opm/core/linalg/LinearSolverInterface.hpp index 46e4cc5ff..b7c8da40d 100644 --- a/opm/core/linalg/LinearSolverInterface.hpp +++ b/opm/core/linalg/LinearSolverInterface.hpp @@ -71,6 +71,14 @@ namespace Opm const double* rhs, double* solution) const = 0; + /// Set tolerance for the linear solver. + /// \param[in] tol tolerance value + virtual void setTolerance(const double tol) = 0; + + /// Get tolerance for the linear solver. + /// \param[out] tolerance value + virtual double getTolerance() const = 0; + }; diff --git a/opm/core/linalg/LinearSolverIstl.cpp b/opm/core/linalg/LinearSolverIstl.cpp index 33f06b190..791fcb179 100644 --- a/opm/core/linalg/LinearSolverIstl.cpp +++ b/opm/core/linalg/LinearSolverIstl.cpp @@ -172,6 +172,17 @@ namespace Opm return res; } + void LinearSolverIstl::setTolerance(const double tol) + { + linsolver_residual_tolerance_ = tol; + } + + double LinearSolverIstl::getTolerance() const + { + return linsolver_residual_tolerance_; + } + + namespace { diff --git a/opm/core/linalg/LinearSolverIstl.hpp b/opm/core/linalg/LinearSolverIstl.hpp index 8e4e223ff..6de42e725 100644 --- a/opm/core/linalg/LinearSolverIstl.hpp +++ b/opm/core/linalg/LinearSolverIstl.hpp @@ -71,6 +71,15 @@ namespace Opm const double* sa, const double* rhs, double* solution) const; + + /// Set tolerance for the residual in dune istl linear solver. + /// \param[in] tol tolerance value + virtual void setTolerance(const double tol); + + /// Get tolerance ofthe linear solver. + /// \param[out] tolerance value + virtual double getTolerance() const; + private: double linsolver_residual_tolerance_; int linsolver_verbosity_; diff --git a/opm/core/linalg/LinearSolverUmfpack.cpp b/opm/core/linalg/LinearSolverUmfpack.cpp index aec76f40a..c9054d938 100644 --- a/opm/core/linalg/LinearSolverUmfpack.cpp +++ b/opm/core/linalg/LinearSolverUmfpack.cpp @@ -60,5 +60,15 @@ namespace Opm return rep; } + void LinearSolverUmfpack::setTolerance(const double /*tol*/) + { + } + + double LinearSolverUmfpack::getTolerance() const + { + return -1.; + } + + } // namespace Opm diff --git a/opm/core/linalg/LinearSolverUmfpack.hpp b/opm/core/linalg/LinearSolverUmfpack.hpp index 666825c02..f2562a156 100644 --- a/opm/core/linalg/LinearSolverUmfpack.hpp +++ b/opm/core/linalg/LinearSolverUmfpack.hpp @@ -56,6 +56,18 @@ namespace Opm const double* sa, const double* rhs, double* solution) const; + + /// Set tolerance for the linear solver. + /// \param[in] tol tolerance value + /// Not used for UMFPACK solver. + virtual void setTolerance(const double /*tol*/); + + /// Get tolerance for the linear solver. + /// \param[out] tolerance value + /// Not used for UMFPACK solver. Returns -1. + virtual double getTolerance() const; + + }; diff --git a/opm/core/pressure/IncompTpfa.cpp b/opm/core/pressure/IncompTpfa.cpp index 2e551224b..b8ea2f32f 100644 --- a/opm/core/pressure/IncompTpfa.cpp +++ b/opm/core/pressure/IncompTpfa.cpp @@ -40,12 +40,12 @@ namespace Opm /// and N == g.number_of_cells. /// \param[in] gravity Gravity vector. If nonzero, the array should /// have D elements. - /// \param[in] wells The wells argument. Will be used in solution, + /// \param[in] wells The wells argument. Will be used in solution, /// is ignored if NULL IncompTpfa::IncompTpfa(const UnstructuredGrid& g, const double* permeability, const double* gravity, - const LinearSolverInterface& linsolver, + LinearSolverInterface& linsolver, const struct Wells* wells) : grid_(g), linsolver_(linsolver), @@ -294,6 +294,16 @@ namespace Opm } + void IncompTpfa::setTolerance(const double tol) + { + linsolver_.setTolerance(tol); + } + + double IncompTpfa::getTolerance() const { + return linsolver_.getTolerance(); + } + + void IncompTpfa::computeFaceFlux(const std::vector& totmob, const std::vector& omega, const std::vector& src, @@ -327,7 +337,7 @@ namespace Opm if (! wdp.empty()) { F.wdp = &wdp[0]; } ifs_tpfa_assemble(gg, &F, &trans_[0], &gpress_omegaweighted_[0], h_); - + faceflux.resize(grid_.number_of_faces); ifs_tpfa_solution soln = { NULL, NULL, NULL, NULL }; diff --git a/opm/core/pressure/IncompTpfa.hpp b/opm/core/pressure/IncompTpfa.hpp index cbba4fea6..d66535e22 100644 --- a/opm/core/pressure/IncompTpfa.hpp +++ b/opm/core/pressure/IncompTpfa.hpp @@ -56,7 +56,7 @@ namespace Opm IncompTpfa(const UnstructuredGrid& g, const double* permeability, const double* gravity, - const LinearSolverInterface& linsolver, + LinearSolverInterface& linsolver, const struct Wells* wells = 0); /// Destructor. @@ -155,9 +155,18 @@ namespace Opm /// Expose read-only reference to internal half-transmissibility. const ::std::vector& getHalfTrans() const { return htrans_; } + /// Set tolerance for the linear solver. + /// \param[in] tol tolerance value + void setTolerance(const double tol); + + /// Get tolerance of the linear solver. + /// \param[out] tolerance value + double getTolerance() const; + + private: const UnstructuredGrid& grid_; - const LinearSolverInterface& linsolver_; + LinearSolverInterface& linsolver_; ::std::vector htrans_; ::std::vector trans_ ; ::std::vector gpress_;