diff --git a/opm/models/blackoil/blackoilconvectivemixingmodule.hh b/opm/models/blackoil/blackoilconvectivemixingmodule.hh index e57f51032..ac46d3bfd 100644 --- a/opm/models/blackoil/blackoilconvectivemixingmodule.hh +++ b/opm/models/blackoil/blackoilconvectivemixingmodule.hh @@ -28,15 +28,16 @@ #ifndef EWOMS_CONVECTIVEMIXING_MODULE_HH #define EWOMS_CONVECTIVEMIXING_MODULE_HH -#include -#include -#include -#include - +#include "opm/material/common/MathToolbox.hpp" #include -#include -#include +#include +#include + +#include + +#include +#include namespace Opm { @@ -119,8 +120,8 @@ public: const Scalar, const ConvectiveMixingModuleParam&) {} - }; + template class BlackOilConvectiveMixingModule { @@ -150,7 +151,10 @@ public: }; #if HAVE_ECL_INPUT - static void beginEpisode(const EclipseState& eclState, const Schedule& schedule, const int episodeIdx, ConvectiveMixingModuleParam& info) + static void beginEpisode(const EclipseState& eclState, + const Schedule& schedule, + const int episodeIdx, + ConvectiveMixingModuleParam& info) { // check that Xhi and Psi didn't change std::size_t numRegions = eclState.runspec().tabdims().getNumPVTTables(); @@ -211,8 +215,10 @@ public: const auto salt_ex = Toolbox::value(intQuantsEx.fluidState().saltConcentration()); const auto bLiquidEx = (FluidSystem::phaseIsActive(waterPhaseIdx)) ? - FluidSystem::waterPvt().inverseFormationVolumeFactor(intQuantsEx.pvtRegionIndex(), t_ex, p_ex, 0.0, salt_ex) : - FluidSystem::oilPvt().inverseFormationVolumeFactor(intQuantsEx.pvtRegionIndex(), t_ex, p_ex, 0.0); + FluidSystem::waterPvt().inverseFormationVolumeFactor(intQuantsEx.pvtRegionIndex(), + t_ex, p_ex, Scalar{0.0}, salt_ex) : + FluidSystem::oilPvt().inverseFormationVolumeFactor(intQuantsEx.pvtRegionIndex(), + t_ex, p_ex, Scalar{0.0}); const auto& refDensityLiquidEx = (FluidSystem::phaseIsActive(waterPhaseIdx)) ? FluidSystem::waterPvt().waterReferenceDensity(intQuantsEx.pvtRegionIndex()) : @@ -373,12 +379,9 @@ public: } } } - }; - + } }; } + #endif - - - diff --git a/opm/models/blackoil/blackoilmicpmodules.hh b/opm/models/blackoil/blackoilmicpmodules.hh index ac56397fa..da49371a1 100644 --- a/opm/models/blackoil/blackoilmicpmodules.hh +++ b/opm/models/blackoil/blackoilmicpmodules.hh @@ -127,7 +127,13 @@ public: MICPpara.getMaximumUreaConcentration(), MICPpara.getToleranceBeforeClogging()); // obtain the porosity for the clamp in the blackoilnewtonmethod - params_.phi_ = eclState.fieldProps().get_double("PORO"); + if constexpr (std::is_same_v) { + const auto phi = eclState.fieldProps().get_double("PORO"); + params_.phi_.resize(phi.size()); + std::copy(phi.begin(), phi.end(), params_.phi_.begin()); + } else { + params_.phi_ = eclState.fieldProps().get_double("PORO"); + } } #endif diff --git a/opm/models/blackoil/blackoilnewtonmethod.hh b/opm/models/blackoil/blackoilnewtonmethod.hh index 847bb7c25..b160dc62f 100644 --- a/opm/models/blackoil/blackoilnewtonmethod.hh +++ b/opm/models/blackoil/blackoilnewtonmethod.hh @@ -324,13 +324,13 @@ protected: else if (enableExtbo && pvIdx == Indices::zFractionIdx) { // z fraction updates are also subject to the Appleyard chop const auto& curr = currentValue[Indices::zFractionIdx]; // or currentValue[pvIdx] given the block condition - delta = std::clamp(delta, curr - 1.0, curr); + delta = std::clamp(delta, curr - Scalar{1.0}, curr); } else if (enablePolymerWeight && pvIdx == Indices::polymerMoleWeightIdx) { const double sign = delta >= 0. ? 1. : -1.; // maximum change of polymer molecular weight, the unit is MDa. // applying this limit to stabilize the simulation. The value itself is still experimental. - const double maxMolarWeightChange = 100.0; + const Scalar maxMolarWeightChange = 100.0; delta = sign * std::min(std::abs(delta), maxMolarWeightChange); delta *= satAlpha; } @@ -341,8 +341,8 @@ protected: else if (enableBrine && pvIdx == Indices::saltConcentrationIdx && enableSaltPrecipitation && currentValue.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp) { - const double maxSaltSaturationChange = 0.1; - const double sign = delta >= 0. ? 1. : -1.; + const Scalar maxSaltSaturationChange = 0.1; + const Scalar sign = delta >= 0. ? 1. : -1.; delta = sign * std::min(std::abs(delta), maxSaltSaturationChange); } @@ -352,19 +352,19 @@ protected: // keep the solvent saturation between 0 and 1 if (enableSolvent && pvIdx == Indices::solventSaturationIdx) { if (currentValue.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss ) - nextValue[pvIdx] = std::min(std::max(nextValue[pvIdx], 0.0), 1.0); + nextValue[pvIdx] = std::min(std::max(nextValue[pvIdx], Scalar{0.0}), Scalar{1.0}); } // keep the z fraction between 0 and 1 if (enableExtbo && pvIdx == Indices::zFractionIdx) - nextValue[pvIdx] = std::min(std::max(nextValue[pvIdx], 0.0), 1.0); + nextValue[pvIdx] = std::min(std::max(nextValue[pvIdx], Scalar{0.0}), Scalar{1.0}); // keep the polymer concentration above 0 if (enablePolymer && pvIdx == Indices::polymerConcentrationIdx) - nextValue[pvIdx] = std::max(nextValue[pvIdx], 0.0); + nextValue[pvIdx] = std::max(nextValue[pvIdx], Scalar{0.0}); if (enablePolymerWeight && pvIdx == Indices::polymerMoleWeightIdx) { - nextValue[pvIdx] = std::max(nextValue[pvIdx], 0.0); + nextValue[pvIdx] = std::max(nextValue[pvIdx], Scalar{0.0}); const double polymerConcentration = nextValue[Indices::polymerConcentrationIdx]; if (polymerConcentration < 1.e-10) nextValue[pvIdx] = 0.0; @@ -372,15 +372,15 @@ protected: // keep the foam concentration above 0 if (enableFoam && pvIdx == Indices::foamConcentrationIdx) - nextValue[pvIdx] = std::max(nextValue[pvIdx], 0.0); + nextValue[pvIdx] = std::max(nextValue[pvIdx], Scalar{0.0}); if (enableBrine && pvIdx == Indices::saltConcentrationIdx) { // keep the salt concentration above 0 if (!enableSaltPrecipitation || (enableSaltPrecipitation && currentValue.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Cs)) - nextValue[pvIdx] = std::max(nextValue[pvIdx], 0.0); + nextValue[pvIdx] = std::max(nextValue[pvIdx], Scalar{0.0}); // keep the salt saturation below upperlimit if ((enableSaltPrecipitation && currentValue.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp)) - nextValue[pvIdx] = std::min(nextValue[pvIdx], 1.0-1.e-8); + nextValue[pvIdx] = std::min(nextValue[pvIdx], Scalar{1.0-1.e-8}); } // keep the temperature within given values @@ -398,15 +398,15 @@ protected: // concentration (the urea concentration has been scaled by 10). For // the biofilm and calcite, we set this value equal to the porosity minus the clogging tolerance. if (enableMICP && pvIdx == Indices::microbialConcentrationIdx) - nextValue[pvIdx] = std::clamp(nextValue[pvIdx], 0.0, MICPModule::densityBiofilm()); + nextValue[pvIdx] = std::clamp(nextValue[pvIdx], Scalar{0.0}, MICPModule::densityBiofilm()); if (enableMICP && pvIdx == Indices::oxygenConcentrationIdx) - nextValue[pvIdx] = std::clamp(nextValue[pvIdx], 0.0, MICPModule::maximumOxygenConcentration()); + nextValue[pvIdx] = std::clamp(nextValue[pvIdx], Scalar{0.0}, MICPModule::maximumOxygenConcentration()); if (enableMICP && pvIdx == Indices::ureaConcentrationIdx) - nextValue[pvIdx] = std::clamp(nextValue[pvIdx], 0.0, MICPModule::maximumUreaConcentration()); + nextValue[pvIdx] = std::clamp(nextValue[pvIdx], Scalar{0.0}, MICPModule::maximumUreaConcentration()); if (enableMICP && pvIdx == Indices::biofilmConcentrationIdx) - nextValue[pvIdx] = std::clamp(nextValue[pvIdx], 0.0, MICPModule::phi()[globalDofIdx] - MICPModule::toleranceBeforeClogging()); + nextValue[pvIdx] = std::clamp(nextValue[pvIdx], Scalar{0.0}, MICPModule::phi()[globalDofIdx] - MICPModule::toleranceBeforeClogging()); if (enableMICP && pvIdx == Indices::calciteConcentrationIdx) - nextValue[pvIdx] = std::clamp(nextValue[pvIdx], 0.0, MICPModule::phi()[globalDofIdx] - MICPModule::toleranceBeforeClogging()); + nextValue[pvIdx] = std::clamp(nextValue[pvIdx], Scalar{0.0}, MICPModule::phi()[globalDofIdx] - MICPModule::toleranceBeforeClogging()); } // switch the new primary variables to something which is physically meaningful. diff --git a/opm/models/blackoil/blackoilprimaryvariables.hh b/opm/models/blackoil/blackoilprimaryvariables.hh index e36b0c1cd..e9496bad5 100644 --- a/opm/models/blackoil/blackoilprimaryvariables.hh +++ b/opm/models/blackoil/blackoilprimaryvariables.hh @@ -849,7 +849,7 @@ public: soMax); Scalar rs = (*this)[Indices::compositionSwitchIdx]; - if (rs > std::min(rsMax, rsSat*(1.0 + eps))) { + if (rs > std::min(rsMax, rsSat*(Scalar{1.0} + eps))) { // the gas phase appears, i.e., switch the primary variables to GasMeaning::Sg setPrimaryVarsMeaningGas(GasMeaning::Sg); (*this)[Indices::compositionSwitchIdx] = 0.0; // hydrocarbon gas saturation @@ -876,7 +876,7 @@ public: soMax); Scalar rv = (*this)[Indices::compositionSwitchIdx]; - if (rv > std::min(rvMax, rvSat*(1.0 + eps))) { + if (rv > std::min(rvMax, rvSat*(Scalar{1.0} + eps))) { // switch to phase equilibrium mode because the oil phase appears. here // we also need the capillary pressures to calculate the oil phase // pressure using the gas phase pressure @@ -925,10 +925,10 @@ public: ssol =(*this) [Indices::solventSaturationIdx]; Scalar so = 1.0 - sw - sg - ssol; - sw = std::min(std::max(sw,0.0),1.0); - so = std::min(std::max(so,0.0),1.0); - sg = std::min(std::max(sg,0.0),1.0); - ssol = std::min(std::max(ssol,0.0),1.0); + sw = std::min(std::max(sw, Scalar{0.0}), Scalar{1.0}); + so = std::min(std::max(so, Scalar{0.0}), Scalar{1.0}); + sg = std::min(std::max(sg, Scalar{0.0}), Scalar{1.0}); + ssol = std::min(std::max(ssol, Scalar{0.0}), Scalar{1.0}); Scalar st = sw + so + sg + ssol; sw = sw/st; sg = sg/st; @@ -981,7 +981,7 @@ public: template void serializeOp(Serializer& serializer) { - using FV = Dune::FieldVector()>; + using FV = Dune::FieldVector()>; serializer(static_cast(*this)); serializer(primaryVarsMeaningWater_); serializer(primaryVarsMeaningPressure_); diff --git a/opm/models/discretization/common/fvbasediscretization.hh b/opm/models/discretization/common/fvbasediscretization.hh index 610609a30..2765164cf 100644 --- a/opm/models/discretization/common/fvbasediscretization.hh +++ b/opm/models/discretization/common/fvbasediscretization.hh @@ -1633,7 +1633,7 @@ public: } // do the normalization of the relative error - Scalar alpha = std::max(1e-20, + Scalar alpha = std::max(Scalar{1e-20}, std::max(std::abs(maxRelErr), std::abs(minRelErr))); for (unsigned globalIdx = 0; globalIdx < numGridDof; ++ globalIdx) diff --git a/opm/models/nonlinear/newtonmethod.hh b/opm/models/nonlinear/newtonmethod.hh index 222a05a72..aa99ae354 100644 --- a/opm/models/nonlinear/newtonmethod.hh +++ b/opm/models/nonlinear/newtonmethod.hh @@ -454,13 +454,13 @@ public: if (numIterations_ > targetIterations_()) { Scalar percent = Scalar(numIterations_ - targetIterations_())/targetIterations_(); Scalar nextDt = std::max(problem().minTimeStepSize(), - oldDt/(1.0 + percent)); + oldDt / (Scalar{1.0} + percent)); return nextDt; } Scalar percent = Scalar(targetIterations_() - numIterations_)/targetIterations_(); Scalar nextDt = std::max(problem().minTimeStepSize(), - oldDt*(1.0 + percent/1.2)); + oldDt*(Scalar{1.0} + percent / Scalar{1.2})); return nextDt; }