From 080afc453274941e0e869607cad9beb4441f6fa1 Mon Sep 17 00:00:00 2001 From: Kai Bao Date: Fri, 26 Sep 2025 17:18:32 +0200 Subject: [PATCH] making WellInterfaceIndices template parameter numEq --- .../wells/MultisegmentWellAssemble.hpp | 6 +-- opm/simulators/wells/MultisegmentWellEval.cpp | 26 ++++++------- opm/simulators/wells/MultisegmentWellEval.hpp | 6 +-- .../MultisegmentWellPrimaryVariables.hpp | 6 +-- .../wells/MultisegmentWell_impl.hpp | 12 +++--- .../wells/StandardWellConnections.cpp | 2 +- .../wells/StandardWellConnections.hpp | 8 ++-- opm/simulators/wells/StandardWellEval.cpp | 2 +- opm/simulators/wells/StandardWellEval.hpp | 6 +-- .../wells/StandardWellPrimaryVariables.hpp | 6 +-- opm/simulators/wells/StandardWell_impl.hpp | 2 +- opm/simulators/wells/WellInterface.hpp | 4 +- opm/simulators/wells/WellInterfaceIndices.cpp | 37 ++++++++++++------- opm/simulators/wells/WellInterfaceIndices.hpp | 6 +-- opm/simulators/wells/WellInterface_impl.hpp | 2 +- 15 files changed, 71 insertions(+), 60 deletions(-) diff --git a/opm/simulators/wells/MultisegmentWellAssemble.hpp b/opm/simulators/wells/MultisegmentWellAssemble.hpp index 1270f3d9ba5..a17e87797ca 100644 --- a/opm/simulators/wells/MultisegmentWellAssemble.hpp +++ b/opm/simulators/wells/MultisegmentWellAssemble.hpp @@ -35,7 +35,7 @@ template class Mul template class MultisegmentWellPrimaryVariables; class Schedule; class SummaryState; -template class WellInterfaceIndices; +template class WellInterfaceIndices; template class WellState; //! \brief Class handling assemble of the equation system for MultisegmentWell. @@ -58,7 +58,7 @@ class MultisegmentWellAssemble using EvalWell = DenseAd::Evaluation; //! \brief Constructor initializes reference to well. - explicit MultisegmentWellAssemble(const WellInterfaceIndices& well) + explicit MultisegmentWellAssemble(const WellInterfaceIndices& well) : well_(well) {} @@ -136,7 +136,7 @@ class MultisegmentWellAssemble Equations& eqns) const; private: - const WellInterfaceIndices& well_; //!< Reference to well + const WellInterfaceIndices& well_; //!< Reference to well }; } diff --git a/opm/simulators/wells/MultisegmentWellEval.cpp b/opm/simulators/wells/MultisegmentWellEval.cpp index 58ec821e323..1847e27eb7c 100644 --- a/opm/simulators/wells/MultisegmentWellEval.cpp +++ b/opm/simulators/wells/MultisegmentWellEval.cpp @@ -53,7 +53,7 @@ namespace Opm template MultisegmentWellEval:: -MultisegmentWellEval(WellInterfaceIndices& baseif, const ParallelWellInfo& pw_info) +MultisegmentWellEval(WellInterfaceIndices& baseif, const ParallelWellInfo& pw_info) : MultisegmentWellGeneric(baseif) , pw_info_(pw_info) , baseif_(baseif) @@ -214,12 +214,12 @@ assembleAccelerationPressureLoss(const int seg, const int seg_upwind = segments_.upwinding_segment(seg); // acceleration term is *subtracted* from pressure equation - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assembleAccelerationTerm(seg, seg, seg_upwind, signed_velocity_head, linSys_); if (seg != seg_upwind) {// special treatment for reverse flow // extra derivatives are *added* to Jacobian (hence minus) const EvalWell extra_derivatives = -segments_.accelerationPressureLossContribution(seg, seg_area, /*extra_derivatives*/ true); - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assemblePressureEqExtraDerivatives(seg, seg_upwind, extra_derivatives, linSys_); } @@ -232,13 +232,13 @@ assembleAccelerationPressureLoss(const int seg, segments.pressure_drop_accel[seg] -= signed_velocity_head_inlet.value(); const int inlet_upwind = segments_.upwinding_segment(inlet); - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assembleAccelerationTerm(seg, inlet, inlet_upwind, -signed_velocity_head_inlet, linSys_); if (inlet != inlet_upwind) {// special treatment for reverse flow // extra derivatives are *added* to Jacobian (hence minus minus) const EvalWell extra_derivatives_inlet = segments_.accelerationPressureLossContribution(inlet, inlet_area, /*extra_derivatives*/ true); // in this case inlet_upwind = seg - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assemblePressureEqExtraDerivatives(seg, inlet_upwind, extra_derivatives_inlet, linSys_); } } @@ -269,7 +269,7 @@ assembleDefaultPressureEq(const int seg, if (reverseFlow){ // call function once again to obtain/assemble remaining derivatives extra_derivatives = -segments_.getFrictionPressureLoss(seg, /*extra_reverse_flow_derivatives*/ true); - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assemblePressureEqExtraDerivatives(seg, seg_upwind, extra_derivatives, linSys_); } pressure_equation -= friction_pressure_drop; @@ -281,7 +281,7 @@ assembleDefaultPressureEq(const int seg, const int outlet_segment_index = this->segmentNumberToIndex(this->segmentSet()[seg].outletSegment()); const EvalWell outlet_pressure = primary_variables_.getSegmentPressure(outlet_segment_index); - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assemblePressureEq(seg, seg_upwind, outlet_segment_index, pressure_equation, outlet_pressure, linSys_); @@ -305,7 +305,7 @@ assembleICDPressureEq(const int seg, if (const auto& segment = this->segmentSet()[seg]; (segment.segmentType() == Segment::SegmentType::VALVE) && (segment.valve().status() == Opm::ICDStatus::SHUT) ) { // we use a zero rate equation to handle SHUT valve - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assembleTrivialEq(seg, this->primary_variables_.eval(seg)[WQTotal].value(), linSys_); auto& ws = well_state.well(baseif_.indexOfWell()); @@ -351,7 +351,7 @@ assembleICDPressureEq(const int seg, } } if (reverseFlow){ - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assemblePressureEqExtraDerivatives(seg, seg_upwind, extra_derivatives, linSys_); } @@ -363,7 +363,7 @@ assembleICDPressureEq(const int seg, const int outlet_segment_index = this->segmentNumberToIndex(this->segmentSet()[seg].outletSegment()); const EvalWell outlet_pressure = primary_variables_.getSegmentPressure(outlet_segment_index); - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assemblePressureEq(seg, seg_upwind, outlet_segment_index, pressure_equation, outlet_pressure, linSys_); @@ -388,15 +388,15 @@ assembleAccelerationAndHydroPressureLosses(const int seg, auto& ws = well_state.well(baseif_.indexOfWell()); auto& segments = ws.segments; if (!use_average_density){ - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assembleHydroPressureLoss(seg, seg, hydro_pressure_drop_seg, linSys_); segments.pressure_drop_hydrostatic[seg] = hydro_pressure_drop_seg.value(); } else { const int seg_outlet = this->segmentNumberToIndex(this->segmentSet()[seg].outletSegment()); const auto hydro_pressure_drop_outlet = segments_.getHydroPressureLoss(seg, seg_outlet); - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assembleHydroPressureLoss(seg, seg, 0.5*hydro_pressure_drop_seg, linSys_); - MultisegmentWellAssemble(baseif_). + MultisegmentWellAssemble(baseif_). assembleHydroPressureLoss(seg, seg_outlet, 0.5*hydro_pressure_drop_outlet, linSys_); segments.pressure_drop_hydrostatic[seg] = 0.5*hydro_pressure_drop_seg.value() + 0.5*hydro_pressure_drop_outlet.value(); } diff --git a/opm/simulators/wells/MultisegmentWellEval.hpp b/opm/simulators/wells/MultisegmentWellEval.hpp index 43ab1f24d5e..0f903cf146b 100644 --- a/opm/simulators/wells/MultisegmentWellEval.hpp +++ b/opm/simulators/wells/MultisegmentWellEval.hpp @@ -39,7 +39,7 @@ class ConvergenceReport; class Schedule; class SummaryState; -template class WellInterfaceIndices; +template class WellInterfaceIndices; template class WellState; template @@ -73,7 +73,7 @@ class MultisegmentWellEval : public MultisegmentWellGeneric& pw_info_; protected: - MultisegmentWellEval(WellInterfaceIndices& baseif, const ParallelWellInfo& pw_info); + MultisegmentWellEval(WellInterfaceIndices& baseif, const ParallelWellInfo& pw_info); void initMatrixAndVectors(); @@ -137,7 +137,7 @@ class MultisegmentWellEval : public MultisegmentWellGeneric& baseif_; + const WellInterfaceIndices& baseif_; Equations linSys_; //!< The equation system PrimaryVariables primary_variables_; //!< The primary variables diff --git a/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp b/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp index 590f87c7b83..db6563a28c6 100644 --- a/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp +++ b/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp @@ -37,7 +37,7 @@ namespace Opm class DeferredLogger; template class MultisegmentWellGeneric; -template class WellInterfaceIndices; +template class WellInterfaceIndices; template class WellState; template @@ -75,7 +75,7 @@ class MultisegmentWellPrimaryVariables using Equations = MultisegmentWellEquations; using BVectorWell = typename Equations::BVectorWell; - explicit MultisegmentWellPrimaryVariables(const WellInterfaceIndices& well) + explicit MultisegmentWellPrimaryVariables(const WellInterfaceIndices& well) : well_(well) {} @@ -165,7 +165,7 @@ class MultisegmentWellPrimaryVariables //! \details Contains derivatives and are used in AD calculation std::vector> evaluation_; - const WellInterfaceIndices& well_; //!< Reference to well interface + const WellInterfaceIndices& well_; //!< Reference to well interface static constexpr double bhp_lower_limit = 1. * unit::barsa - 1. * unit::Pascal; static constexpr double seg_pres_lower_limit = 0.; diff --git a/opm/simulators/wells/MultisegmentWell_impl.hpp b/opm/simulators/wells/MultisegmentWell_impl.hpp index 81980baff48..3dc09a5660a 100644 --- a/opm/simulators/wells/MultisegmentWell_impl.hpp +++ b/opm/simulators/wells/MultisegmentWell_impl.hpp @@ -70,7 +70,7 @@ namespace Opm const int index_of_well, const std::vector>& perf_data) : Base(well, pw_info, time_step, param, rate_converter, pvtRegionIdx, num_conservation_quantities, num_phases, index_of_well, perf_data) - , MSWEval(static_cast&>(*this), pw_info) + , MSWEval(static_cast&>(*this), pw_info) , regularize_(false) , segment_fluid_initial_(this->numberOfSegments(), std::vector(this->num_conservation_quantities_, 0.0)) { @@ -1865,7 +1865,7 @@ namespace Opm this->connectionRates_[local_perf_index][comp_idx] = Base::restrictEval(cq_s_effective); - MultisegmentWellAssemble(*this). + MultisegmentWellAssemble(*this). assemblePerforationEq(seg, local_perf_index, comp_idx, cq_s_effective, this->linSys_); } } @@ -1894,7 +1894,7 @@ namespace Opm for (int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) { const EvalWell accumulation_term = regularization_factor * (segment_surface_volume * this->primary_variables_.surfaceVolumeFraction(seg, comp_idx) - segment_fluid_initial_[seg][comp_idx]) / dt; - MultisegmentWellAssemble(*this). + MultisegmentWellAssemble(*this). assembleAccumulationTerm(seg, comp_idx, accumulation_term, this->linSys_); } } @@ -1907,7 +1907,7 @@ namespace Opm seg_upwind, comp_idx) * this->well_efficiency_factor_; - MultisegmentWellAssemble(*this). + MultisegmentWellAssemble(*this). assembleOutflowTerm(seg, seg_upwind, comp_idx, segment_rate, this->linSys_); } } @@ -1922,7 +1922,7 @@ namespace Opm inlet_upwind, comp_idx) * this->well_efficiency_factor_; - MultisegmentWellAssemble(*this). + MultisegmentWellAssemble(*this). assembleInflowTerm(seg, inlet, inlet_upwind, comp_idx, inlet_rate, this->linSys_); } } @@ -1933,7 +1933,7 @@ namespace Opm const auto& summaryState = simulator.vanguard().summaryState(); const Schedule& schedule = simulator.vanguard().schedule(); const bool stopped_or_zero_target = this->stoppedOrZeroRateTarget(simulator, well_state, deferred_logger); - MultisegmentWellAssemble(*this). + MultisegmentWellAssemble(*this). assembleControlEq(well_state, group_state, schedule, diff --git a/opm/simulators/wells/StandardWellConnections.cpp b/opm/simulators/wells/StandardWellConnections.cpp index 0bd6c8e2e42..3aac92ebc0f 100644 --- a/opm/simulators/wells/StandardWellConnections.cpp +++ b/opm/simulators/wells/StandardWellConnections.cpp @@ -162,7 +162,7 @@ namespace Opm template StandardWellConnections:: -StandardWellConnections(const WellInterfaceIndices& well) +StandardWellConnections(const WellInterfaceIndices& well) : well_(well) , perf_densities_(well.numLocalPerfs()) , perf_pressure_diffs_(well.numLocalPerfs()) diff --git a/opm/simulators/wells/StandardWellConnections.hpp b/opm/simulators/wells/StandardWellConnections.hpp index 24a532d60c1..40d45c34bc1 100644 --- a/opm/simulators/wells/StandardWellConnections.hpp +++ b/opm/simulators/wells/StandardWellConnections.hpp @@ -36,7 +36,7 @@ namespace Opm class DeferredLogger; enum class Phase; -template class WellInterfaceIndices; +template class WellInterfaceIndices; template class WellState; template class PerfData; @@ -46,7 +46,7 @@ class StandardWellConnections public: using Scalar = typename FluidSystem::Scalar; using IndexTraits = typename FluidSystem::IndexTraitsType; - explicit StandardWellConnections(const WellInterfaceIndices& well); + explicit StandardWellConnections(const WellInterfaceIndices& well); struct Properties { @@ -106,7 +106,7 @@ class StandardWellConnections Scalar pressure_diff(const unsigned perf) const { return perf_pressure_diffs_[perf]; } - using Eval = typename WellInterfaceIndices::Eval; + using Eval = typename WellInterfaceIndices::Eval; using EvalWell = typename StandardWellPrimaryVariables::EvalWell; Eval connectionRateBrine(Scalar& rate, @@ -168,7 +168,7 @@ class StandardWellConnections copyInPerforationRates(const Properties& props, const PerfData& perf_data) const; - const WellInterfaceIndices& well_; //!< Reference to well interface + const WellInterfaceIndices& well_; //!< Reference to well interface std::vector perf_densities_; //!< densities of the fluid in each perforation std::vector perf_pressure_diffs_; //!< // pressure drop between different perforations diff --git a/opm/simulators/wells/StandardWellEval.cpp b/opm/simulators/wells/StandardWellEval.cpp index 084940ad9f4..e18d0a58b40 100644 --- a/opm/simulators/wells/StandardWellEval.cpp +++ b/opm/simulators/wells/StandardWellEval.cpp @@ -46,7 +46,7 @@ namespace Opm { template StandardWellEval:: -StandardWellEval(const WellInterfaceIndices& baseif) +StandardWellEval(const WellInterfaceIndices& baseif) : baseif_(baseif) , primary_variables_(baseif_) , F0_(numWellConservationEq) diff --git a/opm/simulators/wells/StandardWellEval.hpp b/opm/simulators/wells/StandardWellEval.hpp index f5c2372b8bf..46683e0fc4e 100644 --- a/opm/simulators/wells/StandardWellEval.hpp +++ b/opm/simulators/wells/StandardWellEval.hpp @@ -38,7 +38,7 @@ class DeferredLogger; class Schedule; class SummaryState; template class WellContributions; -template class WellInterfaceIndices; +template class WellInterfaceIndices; template class WellState; template @@ -69,9 +69,9 @@ class StandardWellEval { return linSys_; } protected: - explicit StandardWellEval(const WellInterfaceIndices& baseif); + explicit StandardWellEval(const WellInterfaceIndices& baseif); - const WellInterfaceIndices& baseif_; + const WellInterfaceIndices& baseif_; EvalWell extendEval(const Eval& in) const; diff --git a/opm/simulators/wells/StandardWellPrimaryVariables.hpp b/opm/simulators/wells/StandardWellPrimaryVariables.hpp index 16cf78c0612..6d6505b3a16 100644 --- a/opm/simulators/wells/StandardWellPrimaryVariables.hpp +++ b/opm/simulators/wells/StandardWellPrimaryVariables.hpp @@ -33,7 +33,7 @@ namespace Opm { class DeferredLogger; -template class WellInterfaceIndices; +template class WellInterfaceIndices; template class WellState; //! \brief Class holding primary variables for StandardWell. @@ -89,7 +89,7 @@ class StandardWellPrimaryVariables { using BVectorWell = typename StandardWellEquations::BVectorWell; //! \brief Constructor initializes reference to well interface. - explicit StandardWellPrimaryVariables(const WellInterfaceIndices& well) + explicit StandardWellPrimaryVariables(const WellInterfaceIndices& well) : well_(well) {} @@ -171,7 +171,7 @@ class StandardWellPrimaryVariables { //! \details Contain derivatives and are used in AD calculation std::vector evaluation_; - const WellInterfaceIndices& well_; //!< Reference to well interface + const WellInterfaceIndices& well_; //!< Reference to well interface //! \brief Total number of the well equations and primary variables. //! \details There might be extra equations be used, numWellEq will be updated during the initialization diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index c7d12ce4050..6dee29f360a 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -60,7 +60,7 @@ namespace Opm const int index_of_well, const std::vector>& perf_data) : Base(well, pw_info, time_step, param, rate_converter, pvtRegionIdx, num_conservation_quantities, num_phases, index_of_well, perf_data) - , StdWellEval(static_cast&>(*this)) + , StdWellEval(static_cast&>(*this)) , regularize_(false) { assert(this->num_conservation_quantities_ == numWellConservationEq); diff --git a/opm/simulators/wells/WellInterface.hpp b/opm/simulators/wells/WellInterface.hpp index 89d59d3dbe2..aad13e999df 100644 --- a/opm/simulators/wells/WellInterface.hpp +++ b/opm/simulators/wells/WellInterface.hpp @@ -72,10 +72,10 @@ class WellProductionProperties; template class WellInterface : public WellInterfaceIndices, - GetPropType> + GetPropType::numEq> { using Base = WellInterfaceIndices, - GetPropType>; + GetPropType::numEq>; public: using Grid = GetPropType; using Simulator = GetPropType; diff --git a/opm/simulators/wells/WellInterfaceIndices.cpp b/opm/simulators/wells/WellInterfaceIndices.cpp index 7812b38d92e..0eaf9f28a6f 100644 --- a/opm/simulators/wells/WellInterfaceIndices.cpp +++ b/opm/simulators/wells/WellInterfaceIndices.cpp @@ -24,17 +24,13 @@ #include -#include -#include -#include - #include namespace Opm { -template -WellInterfaceIndices:: +template +WellInterfaceIndices:: WellInterfaceIndices(const Well& well, const ParallelWellInfo& parallel_well_info, const int time_step, @@ -58,9 +54,9 @@ WellInterfaceIndices(const Well& well, { } -template -typename WellInterfaceIndices::Scalar -WellInterfaceIndices:: +template +typename WellInterfaceIndices::Scalar +WellInterfaceIndices:: scalingFactor(const int phaseIdx) const { if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && @@ -72,7 +68,8 @@ scalingFactor(const int phaseIdx) const if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx) == phaseIdx) return 0.01; - if (Indices::enableSolvent && phaseIdx == Indices::contiSolventEqIdx ) + const auto& pu = FluidSystem::phaseUsage(); + if (pu.hasSolvent() && phaseIdx == pu.contiSolventEqIdx() ) return 0.01; // we should not come this far @@ -80,12 +77,26 @@ scalingFactor(const int phaseIdx) const return 1.0; } -#include +template +using FS = BlackOilFluidSystem; + +#define INSTANTIATE_TYPE(T, NUMEQ) \ + template class WellInterfaceIndices, NUMEQ>; -INSTANTIATE_TYPE_INDICES(WellInterfaceIndices, double) +INSTANTIATE_TYPE(double, 1) +INSTANTIATE_TYPE(double, 2) +INSTANTIATE_TYPE(double, 3) +INSTANTIATE_TYPE(double, 4) +INSTANTIATE_TYPE(double, 5) +INSTANTIATE_TYPE(double, 6) #if FLOW_INSTANTIATE_FLOAT -INSTANTIATE_TYPE_INDICES(WellInterfaceIndices, float) +INSTANTIATE_TYPE(float, 1) +INSTANTIATE_TYPE(float, 2) +INSTANTIATE_TYPE(float, 3) +INSTANTIATE_TYPE(float, 4) +INSTANTIATE_TYPE(float, 5) +INSTANTIATE_TYPE(float, 6) #endif } // namespace Opm diff --git a/opm/simulators/wells/WellInterfaceIndices.hpp b/opm/simulators/wells/WellInterfaceIndices.hpp index bde6290a340..0ff13fe7a05 100644 --- a/opm/simulators/wells/WellInterfaceIndices.hpp +++ b/opm/simulators/wells/WellInterfaceIndices.hpp @@ -29,7 +29,7 @@ namespace Opm { -template +template class WellInterfaceIndices : public WellInterfaceFluidSystem { public: @@ -37,7 +37,7 @@ class WellInterfaceIndices : public WellInterfaceFluidSystem using WellInterfaceFluidSystem::Oil; using WellInterfaceFluidSystem::Water; using Scalar = typename FluidSystem::Scalar; - using Eval = DenseAd::Evaluation; + using Eval = DenseAd::Evaluation; using ModelParameters = typename WellInterfaceFluidSystem::ModelParameters; Scalar scalingFactor(const int phaseIdx) const; @@ -47,7 +47,7 @@ class WellInterfaceIndices : public WellInterfaceFluidSystem { Eval out = 0.0; out.setValue(in.value()); - for (int eqIdx = 0; eqIdx < Indices::numEq; ++eqIdx) { + for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) { out.setDerivative(eqIdx, in.derivative(eqIdx)); } return out; diff --git a/opm/simulators/wells/WellInterface_impl.hpp b/opm/simulators/wells/WellInterface_impl.hpp index a34ba1749b8..00ca742754f 100644 --- a/opm/simulators/wells/WellInterface_impl.hpp +++ b/opm/simulators/wells/WellInterface_impl.hpp @@ -65,7 +65,7 @@ namespace Opm const int num_phases, const int index_of_well, const std::vector>& perf_data) - : WellInterfaceIndices(well, + : WellInterfaceIndices(well, pw_info, time_step, param,