OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #1852 from atgeirr/flexible-linear-solver

Browse Source 
Flexible linear solver
This commit is contained in:
Atgeirr Flø Rasmussen
2019-06-06 14:21:51 +02:00
committed by GitHub
parent 8ae671eb86 d5c3c93c7d
commit 2c5ef367dd
25 changed files with 2241 additions and 16 deletions
Download Patch File Download Diff File Expand all files Collapse all files

66
tests/matr33.txt Normal file
View File

@@ -0,0 +1,66 @@
%%MatrixMarket matrix coordinate real general
% ISTL_STRUCT blocked 3 3
9 9 63
1 1 1.63148e-07
1 2 0.126774
1 3 -6.09774e-05
2 1 4.64465e-10
2 2 -0.0987407
2 3 -5.86038e-05
3 1 5.66966e-08
3 2 -12.0531
3 3 0.0462478
1 4 -8.84841e-11
1 5 -5.64705e-05
1 6 0
2 4 -3.93176e-10
2 5 0
2 6 0
3 4 -4.79944e-08
3 5 0
3 6 0
4 1 -8.89921e-11
4 2 -0.00786145
4 3 0
5 1 -3.53139e-10
5 2 0.0283977
5 3 -3.79289e-05
6 1 -4.31072e-08
6 2 3.46646
6 3 -0.00944335
4 4 1.20862e-10
4 5 0.0829656
4 6 0
5 4 7.40377e-10
5 5 -0.0929483
5 6 -0.0967963
6 4 4.64141e-07
6 5 -9.10276
6 6 16.3772
4 7 -2.168e-11
4 8 -1.38362e-05
4 9 0
5 7 -7.8331e-10
5 8 0
5 9 0
6 7 -1.20687e-07
6 8 0
6 9 0
7 4 -2.17548e-11
7 5 -0.00183118
7 6 0
8 4 -8.00212e-10
8 5 0.0276779
8 6 0.0320899
9 4 -1.56667e-07
9 5 2.7106
9 6 -7.32919
7 7 4.43266e-09
7 8 0.0744616
7 9 0
8 7 1.59375e-07
8 8 -0.0899451
8 9 -0.0940009
9 7 2.86852e-05
9 8 -5.63337
9 9 13.3423

33
tests/options_flexiblesolver.json Normal file
View File

@@ -0,0 +1,33 @@
{
"tol": "0.5",
"maxiter": "20",
"preconditioner": {
"type": "cpr",
"finesmoother": {
"type": "ILU0",
"relaxation": "1.0"
},
"coarsesolver": {
"tol": "0.5",
"maxiter": "20",
"preconditioner": {
"type": "amg",
"maxlevel": "5",
"coarsenTarget": "1000",
"smoother": "ILU0",
"alpha": "0.2",
"beta": "0.0001",
"verbosity": "0",
"iterations": "1",
"relaxation": "1"
},
"verbosity": "0",
"solver": "bicgstab"
},
"verbosity": "11",
"weights_filename" : "weight_cpr.txt",
"pressure_var_index" : "1"
},
"verbosity": "10",
"solver": "bicgstab"
}

9
tests/options_flexiblesolver_simple.json Normal file
View File

@@ -0,0 +1,9 @@
{
"tol": "1e-12",
"maxiter": "200",
"verbosity": "0",
"solver": "bicgstab",
"preconditioner": {
"type": "nothing"
}
}

12
tests/rhs3.txt Normal file
View File

@@ -0,0 +1,12 @@
%%MatrixMarket matrix array real general
% ISTL_STRUCT blocked 3 1
9 1
-4.48332e-07
3.53746e-07
4.31812e-05
1.48051e-07
-5.76325e-07
-0.000243496
-3.81652e-08
1.21761e-06
-9.97615e-05

123
tests/test_flexiblesolver.cpp Normal file
View File

@@ -0,0 +1,123 @@
/*
Copyright 2019 SINTEF Digital, Mathematics and Cybernetics.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <config.h>
#define BOOST_TEST_MODULE OPM_test_FlexibleSolver
#include <boost/test/unit_test.hpp>
#include <dune/common/version.hh>
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
#include <opm/simulators/linalg/FlexibleSolver.hpp>
#include <boost/property_tree/json_parser.hpp>
#include <boost/property_tree/ptree.hpp>
#include <fstream>
#include <iostream>
template <int bz>
Dune::BlockVector<Dune::FieldVector<double, bz>>
testSolver(const boost::property_tree::ptree& prm, const std::string& matrix_filename, const std::string& rhs_filename)
{
using Matrix = Dune::BCRSMatrix<Dune::FieldMatrix<double, bz, bz>>;
using Vector = Dune::BlockVector<Dune::FieldVector<double, bz>>;
Matrix matrix;
{
std::ifstream mfile(matrix_filename);
if (!mfile) {
throw std::runtime_error("Could not read matrix file");
}
readMatrixMarket(matrix, mfile);
}
Vector rhs;
{
std::ifstream rhsfile(rhs_filename);
if (!rhsfile) {
throw std::runtime_error("Could not read rhs file");
}
readMatrixMarket(rhs, rhsfile);
}
Dune::FlexibleSolver<Matrix, Vector> solver(prm, matrix);
Vector x(rhs.size());
Dune::InverseOperatorResult res;
solver.apply(x, rhs, res);
return x;
}
BOOST_AUTO_TEST_CASE(TestFlexibleSolver)
{
namespace pt = boost::property_tree;
pt::ptree prm;
// Read parameters.
{
std::ifstream file("options_flexiblesolver.json");
pt::read_json(file, prm);
// pt::write_json(std::cout, prm);
}
// Test with 1x1 block solvers.
{
const int bz = 1;
auto sol = testSolver<bz>(prm, "matr33.txt", "rhs3.txt");
Dune::BlockVector<Dune::FieldVector<double, bz>> expected {-1.62493,
-1.76435e-06,
1.86991e-10,
-458.542,
2.28308e-06,
-2.45341e-07,
-1.48005,
-5.02264e-07,
-1.049e-05};
BOOST_REQUIRE_EQUAL(sol.size(), expected.size());
for (size_t i = 0; i < sol.size(); ++i) {
for (int row = 0; row < bz; ++row) {
BOOST_CHECK_CLOSE(sol[i][row], expected[i][row], 1e-3);
}
}
}
// Test with 3x3 block solvers.
{
const int bz = 3;
auto sol = testSolver<bz>(prm, "matr33.txt", "rhs3.txt");
Dune::BlockVector<Dune::FieldVector<double, bz>> expected {{-1.62493, -1.76435e-06, 1.86991e-10},
{-458.542, 2.28308e-06, -2.45341e-07},
{-1.48005, -5.02264e-07, -1.049e-05}};
BOOST_REQUIRE_EQUAL(sol.size(), expected.size());
for (size_t i = 0; i < sol.size(); ++i) {
for (int row = 0; row < bz; ++row) {
BOOST_CHECK_CLOSE(sol[i][row], expected[i][row], 1e-3);
}
}
}
}
#else
// Do nothing if we do not have at least Dune 2.6.
BOOST_AUTO_TEST_CASE(DummyTest)
{
BOOST_REQUIRE(true);
}
#endif

229
tests/test_preconditionerfactory.cpp Normal file
View File

@@ -0,0 +1,229 @@
/*
Copyright 2019 SINTEF Digital, Mathematics and Cybernetics.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <config.h>
#define BOOST_TEST_MODULE OPM_test_PreconditionerFactory
#include <boost/test/unit_test.hpp>
#include <dune/common/version.hh>
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
#include <opm/simulators/linalg/PreconditionerFactory.hpp>
#include <opm/simulators/linalg/FlexibleSolver.hpp>
#include <dune/common/fvector.hh>
#include <dune/istl/bvector.hh>
#include <dune/istl/bcrsmatrix.hh>
#include <dune/istl/matrixmarket.hh>
#include <dune/istl/solvers.hh>
#include <boost/property_tree/json_parser.hpp>
#include <boost/property_tree/ptree.hpp>
#include <fstream>
#include <iostream>
template <class X>
class NothingPreconditioner : public Dune::Preconditioner<X, X>
{
public:
virtual void pre(X&, X&) override
{
}
virtual void apply(X& v, const X& d) override
{
v = d;
}
virtual void post(X&) override
{
}
virtual Dune::SolverCategory::Category category() const override
{
return Dune::SolverCategory::sequential;
}
};
template <int bz>
Dune::BlockVector<Dune::FieldVector<double, bz>>
testPrec(const boost::property_tree::ptree& prm, const std::string& matrix_filename, const std::string& rhs_filename)
{
using Matrix = Dune::BCRSMatrix<Dune::FieldMatrix<double, bz, bz>>;
using Vector = Dune::BlockVector<Dune::FieldVector<double, bz>>;
Matrix matrix;
{
std::ifstream mfile(matrix_filename);
if (!mfile) {
throw std::runtime_error("Could not read matrix file");
}
readMatrixMarket(matrix, mfile);
}
Vector rhs;
{
std::ifstream rhsfile(rhs_filename);
if (!rhsfile) {
throw std::runtime_error("Could not read rhs file");
}
readMatrixMarket(rhs, rhsfile);
}
using Operator = Dune::MatrixAdapter<Matrix, Vector, Vector>;
Operator op(matrix);
using PrecFactory = Dune::PreconditionerFactory<Operator>;
auto prec = PrecFactory::create(op, prm.get_child("preconditioner"));
Dune::BiCGSTABSolver<Vector> solver(op, *prec, prm.get<double>("tol"), prm.get<int>("maxiter"), prm.get<int>("verbosity"));
Vector x(rhs.size());
Dune::InverseOperatorResult res;
solver.apply(x, rhs, res);
return x;
}
namespace pt = boost::property_tree;
void test1(const pt::ptree& prm)
{
const int bz = 1;
auto sol = testPrec<bz>(prm, "matr33.txt", "rhs3.txt");
Dune::BlockVector<Dune::FieldVector<double, bz>> expected {-1.62493,
-1.76435e-06,
1.86991e-10,
-458.542,
2.28308e-06,
-2.45341e-07,
-1.48005,
-5.02264e-07,
-1.049e-05};
BOOST_REQUIRE_EQUAL(sol.size(), expected.size());
for (size_t i = 0; i < sol.size(); ++i) {
for (int row = 0; row < bz; ++row) {
BOOST_CHECK_CLOSE(sol[i][row], expected[i][row], 1e-3);
}
}
}
void test3(const pt::ptree& prm)
{
const int bz = 3;
auto sol = testPrec<bz>(prm, "matr33.txt", "rhs3.txt");
Dune::BlockVector<Dune::FieldVector<double, bz>> expected {{-1.62493, -1.76435e-06, 1.86991e-10},
{-458.542, 2.28308e-06, -2.45341e-07},
{-1.48005, -5.02264e-07, -1.049e-05}};
BOOST_REQUIRE_EQUAL(sol.size(), expected.size());
for (size_t i = 0; i < sol.size(); ++i) {
for (int row = 0; row < bz; ++row) {
BOOST_CHECK_CLOSE(sol[i][row], expected[i][row], 1e-3);
}
}
}
BOOST_AUTO_TEST_CASE(TestDefaultPreconditionerFactory)
{
pt::ptree prm;
// Read parameters.
{
std::ifstream file("options_flexiblesolver.json");
pt::read_json(file, prm);
}
// Test with 1x1 block solvers.
test1(prm);
// Test with 3x3 block solvers.
test3(prm);
}
template <int bz>
using M = Dune::BCRSMatrix<Dune::FieldMatrix<double, bz, bz>>;
template <int bz>
using V = Dune::BlockVector<Dune::FieldVector<double, bz>>;
template <int bz>
using O = Dune::MatrixAdapter<M<bz>, V<bz>, V<bz>>;
template <int bz>
using PF = Dune::PreconditionerFactory<O<bz>>;
BOOST_AUTO_TEST_CASE(TestAddingPreconditioner)
{
namespace pt = boost::property_tree;
pt::ptree prm;
// Read parameters.
{
std::ifstream file("options_flexiblesolver_simple.json"); // Requests "nothing" for preconditioner type.
pt::read_json(file, prm);
}
// Test with 1x1 block solvers.
{
const int bz = 1;
BOOST_CHECK_THROW(testPrec<bz>(prm, "matr33.txt", "rhs3.txt"), std::runtime_error);
}
// Test with 3x3 block solvers.
{
const int bz = 3;
BOOST_CHECK_THROW(testPrec<bz>(prm, "matr33.txt", "rhs3.txt"), std::runtime_error);
}
// Add preconditioner to factory for block size 1.
PF<1>::addCreator("nothing", [](const O<1>&, const pt::ptree&) {
return Dune::wrapPreconditioner<NothingPreconditioner<V<1>>>();
});
// Test with 1x1 block solvers.
test1(prm);
// Test with 3x3 block solvers.
{
const int bz = 3;
BOOST_CHECK_THROW(testPrec<bz>(prm, "matr33.txt", "rhs3.txt"), std::runtime_error);
}
// Add preconditioner to factory for block size 3.
PF<3>::addCreator("nothing", [](const O<3>&, const pt::ptree&) {
return Dune::wrapPreconditioner<NothingPreconditioner<V<3>>>();
});
// Test with 1x1 block solvers.
test1(prm);
// Test with 3x3 block solvers.
test3(prm);
}
#else
// Do nothing if we do not have at least Dune 2.6.
BOOST_AUTO_TEST_CASE(DummyTest)
{
BOOST_REQUIRE(true);
}
#endif