Merge pull request #6091 from daavid00/cprBug

Bug fix in linear solver
This commit is contained in:
Atgeirr Flø Rasmussen authored and GitHub committed 2025-03-19 15:17:38 +01:00
commit d6cbee1380
3 files changed
+7 -16

No files matched your search

+3 -10
View File
@@ -161,9 +161,8 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
using ElementMapper = GetPropType<TypeTag, Properties::ElementMapper>;
constexpr static std::size_t pressureIndex = GetPropType<TypeTag, Properties::Indices>::pressureSwitchIdx;
enum { enableMICP = getPropValue<TypeTag, Properties::EnableMICP>() };
enum { enablePolymerMolarWeight = getPropValue<TypeTag, Properties::EnablePolymerMW>() };
constexpr static bool isIncompatibleWithCprw = enableMICP || enablePolymerMolarWeight;
constexpr static bool isIncompatibleWithCprw = enablePolymerMolarWeight;
#if HAVE_MPI
using CommunicationType = Dune::OwnerOverlapCopyCommunication<int,int>;
@@ -218,16 +217,10 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
OPM_TIMEBLOCK(IstlSolver);
if (isIncompatibleWithCprw) {
// Some model variants are incompatible with the CPRW linear solver.
// Polymer injectivity is incompatible with the CPRW linear solver.
if (parameters_[0].linsolver_ == "cprw" || parameters_[0].linsolver_ == "hybrid") {
std::string incompatible_model = "Unknown";
if (enableMICP) {
incompatible_model = "MICP";
} else if (enablePolymerMolarWeight) {
incompatible_model = "Polymer injectivity";
}
OPM_THROW(std::runtime_error,
incompatible_model + " model is incompatible with the CPRW linear solver.\n"
"The polymer injectivity model is incompatible with the CPRW linear solver.\n"
"Choose a different option, for example --linear-solver=ilu0");
}
}
@@ -376,13 +376,13 @@ extractCPRPressureMatrix(PressureMatrix& jacobian,
bweights[0][blockSz-1] = 1.0;
diagElem = 1.0; // better scaling could have used the calculation below if weights were calculated
} else {
for (std::size_t i = 0; i < cell_weights.size(); ++i) {
for (std::size_t i = 0; i < blockSz; ++i) {
bweights[0][i] = cell_weights[i];
}
bweights[0][blockSz-1] = 0.0;
diagElem = 0.0;
const auto& locmat = duneD_[0][0];
for (std::size_t i = 0; i < cell_weights.size(); ++i) {
for (std::size_t i = 0; i < blockSz; ++i) {
diagElem += locmat[i][bhp_var_index] * cell_weights[i];
}
+2 -4
View File
@@ -1193,8 +1193,7 @@ add_test_compareECLFiles(CASENAME micp
SIMULATOR flow
ABS_TOL ${abs_tol}
REL_TOL ${rel_tol}
DIR micp
TEST_ARGS --linear-solver=ilu0)
DIR micp)
add_test_compareECLFiles(CASENAME 0_base_model6
FILENAME 0_BASE_MODEL6
@@ -1223,8 +1222,7 @@ add_test_compareECLFiles(CASENAME base_wt_tracer
ABS_TOL ${abs_tol}
REL_TOL ${rel_tol}
DIR tracer
RESTART_STEP 1,3,7)
RESTART_STEP 1,3,7)
add_test_compareECLFiles(CASENAME min_bhp_1
FILENAME MIN_BHP_1