Skip to content
Snippets Groups Projects
Commit 03254765 authored by Kilian Weishaupt's avatar Kilian Weishaupt
Browse files

[linearSolver] Add version check to avoid deprecation warning

* dune-istl has unified SeqILU0 and SeqILUn to SeqILU
parent f9d21793
No related branches found
No related tags found
1 merge request!777[linearSolver] Add version check to avoid deprecation warning
......@@ -30,6 +30,7 @@
#include <dune/istl/solvers.hh>
#include <dune/istl/superlu.hh>
#include <dune/istl/umfpack.hh>
#include <dune/common/version.hh>
#include <dumux/common/parameters.hh>
#include <dumux/common/properties.hh>
......@@ -177,7 +178,11 @@ public:
template<int precondBlockLevel = 1, class Matrix, class Vector>
bool solve(const Matrix& A, Vector& x, const Vector& b)
{
#if DUNE_VERSION_NEWER(DUNE_ISTL,2,6)
using Preconditioner = Dune::SeqILU<Matrix, Vector, Vector, precondBlockLevel>;
#else
using Preconditioner = Dune::SeqILUn<Matrix, Vector, Vector, precondBlockLevel>;
#endif
using Solver = Dune::BiCGSTABSolver<Vector>;
return IterativePreconditionedSolverImpl::template solve<Preconditioner, Solver>(*this, A, x, b, this->paramGroup());
......@@ -365,7 +370,11 @@ public:
template<int precondBlockLevel = 1, class Matrix, class Vector>
bool solve(const Matrix& A, Vector& x, const Vector& b)
{
#if DUNE_VERSION_NEWER(DUNE_ISTL,2,6)
using Preconditioner = Dune::SeqILU<Matrix, Vector, Vector, precondBlockLevel>;
#else
using Preconditioner = Dune::SeqILUn<Matrix, Vector, Vector, precondBlockLevel>;
#endif
using Solver = Dune::CGSolver<Vector>;
return IterativePreconditionedSolverImpl::template solve<Preconditioner, Solver>(*this, A, x, b, this->paramGroup());
......@@ -587,7 +596,11 @@ public:
template<int precondBlockLevel = 1, class Matrix, class Vector>
bool solve(const Matrix& A, Vector& x, const Vector& b)
{
#if DUNE_VERSION_NEWER(DUNE_ISTL,2,6)
using Preconditioner = Dune::SeqILU<Matrix, Vector, Vector, precondBlockLevel>;
#else
using Preconditioner = Dune::SeqILU0<Matrix, Vector, Vector, precondBlockLevel>;
#endif
using Solver = Dune::BiCGSTABSolver<Vector>;
return IterativePreconditionedSolverImpl::template solveWithILU0Prec<Preconditioner, Solver>(*this, A, x, b, this->paramGroup());
......@@ -623,7 +636,11 @@ public:
template<int precondBlockLevel = 1, class Matrix, class Vector>
bool solve(const Matrix& A, Vector& x, const Vector& b)
{
#if DUNE_VERSION_NEWER(DUNE_ISTL,2,6)
using Preconditioner = Dune::SeqILU<Matrix, Vector, Vector, precondBlockLevel>;
#else
using Preconditioner = Dune::SeqILU0<Matrix, Vector, Vector, precondBlockLevel>;
#endif
using Solver = Dune::CGSolver<Vector>;
return IterativePreconditionedSolverImpl::template solveWithILU0Prec<Preconditioner, Solver>(*this, A, x, b, this->paramGroup());
......@@ -660,7 +677,11 @@ public:
template<int precondBlockLevel = 1, class Matrix, class Vector>
bool solve(const Matrix& A, Vector& x, const Vector& b)
{
#if DUNE_VERSION_NEWER(DUNE_ISTL,2,6)
using Preconditioner = Dune::SeqILU<Matrix, Vector, Vector, precondBlockLevel>;
#else
using Preconditioner = Dune::SeqILU0<Matrix, Vector, Vector, precondBlockLevel>;
#endif
using Solver = Dune::RestartedGMResSolver<Vector>;
return IterativePreconditionedSolverImpl::template solveWithILU0PrecGMRes<Preconditioner, Solver>(*this, A, x, b, this->paramGroup());
......@@ -698,7 +719,11 @@ public:
template<int precondBlockLevel = 1, class Matrix, class Vector>
bool solve(const Matrix& A, Vector& x, const Vector& b)
{
#if DUNE_VERSION_NEWER(DUNE_ISTL,2,6)
using Preconditioner = Dune::SeqILU<Matrix, Vector, Vector, precondBlockLevel>;
#else
using Preconditioner = Dune::SeqILUn<Matrix, Vector, Vector, precondBlockLevel>;
#endif
using Solver = Dune::RestartedGMResSolver<Vector>;
return IterativePreconditionedSolverImpl::template solveWithGMRes<Preconditioner, Solver>(*this, A, x, b, this->paramGroup());
......
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