diff --git a/opm/simulators/wells/MultisegmentWellAssemble.cpp b/opm/simulators/wells/MultisegmentWellAssemble.cpp index 3b16ef468..78cd2be51 100644 --- a/opm/simulators/wells/MultisegmentWellAssemble.cpp +++ b/opm/simulators/wells/MultisegmentWellAssemble.cpp @@ -172,6 +172,34 @@ assemblePressureLoss(const int seg, } } +template +void MultisegmentWellAssemble:: +assemblePressureEq(const int seg, + const int seg_upwind, + const int outlet_segment_index, + const EvalWell& pressure_equation, + const EvalWell& outlet_pressure, + Equations& eqns, + bool wfrac, + bool gfrac) const +{ + eqns.resWell_[seg][SPres] = pressure_equation.value(); + eqns.duneD_[seg][seg][SPres][SPres] += pressure_equation.derivative(SPres + Indices::numEq); + eqns.duneD_[seg][seg][SPres][WQTotal] += pressure_equation.derivative(WQTotal + Indices::numEq); + if (wfrac) { + eqns.duneD_[seg][seg_upwind][SPres][WFrac] += pressure_equation.derivative(WFrac + Indices::numEq); + } + if (gfrac) { + eqns.duneD_[seg][seg_upwind][SPres][GFrac] += pressure_equation.derivative(GFrac + Indices::numEq); + } + + // contribution from the outlet segment + eqns.resWell_[seg][SPres] -= outlet_pressure.value(); + for (int pv_idx = 0; pv_idx < numWellEq; ++pv_idx) { + eqns.duneD_[seg][outlet_segment_index][SPres][pv_idx] = -outlet_pressure.derivative(pv_idx + Indices::numEq); + } +} + #define INSTANCE(...) \ template class MultisegmentWellAssemble,__VA_ARGS__,double>; diff --git a/opm/simulators/wells/MultisegmentWellAssemble.hpp b/opm/simulators/wells/MultisegmentWellAssemble.hpp index bc88d9247..4ded47a60 100644 --- a/opm/simulators/wells/MultisegmentWellAssemble.hpp +++ b/opm/simulators/wells/MultisegmentWellAssemble.hpp @@ -88,6 +88,16 @@ public: const EvalWell& accelerationPressureLoss, Equations& eqns) const; + //! \brief Assemble pressure terms. + void assemblePressureEq(const int seg, + const int seg_upwind, + const int outlet_segment_index, + const EvalWell& pressure_equation, + const EvalWell& outlet_pressure, + Equations& eqns, + bool wfrac = has_wfrac_variable, + bool gfrac = has_gfrac_variable) const; + private: const WellInterfaceIndices& well_; //!< Reference to well }; diff --git a/opm/simulators/wells/MultisegmentWellEval.cpp b/opm/simulators/wells/MultisegmentWellEval.cpp index a978216de..c8b4643fe 100644 --- a/opm/simulators/wells/MultisegmentWellEval.cpp +++ b/opm/simulators/wells/MultisegmentWellEval.cpp @@ -1180,25 +1180,14 @@ assembleDefaultPressureEq(const int seg, segments.pressure_drop_friction[seg] = friction_pressure_drop.value(); } - linSys_.resWell_[seg][SPres] = pressure_equation.value(); - const int seg_upwind = upwinding_segments_[seg]; - linSys_.duneD_[seg][seg][SPres][SPres] += pressure_equation.derivative(SPres + Indices::numEq); - linSys_.duneD_[seg][seg][SPres][WQTotal] += pressure_equation.derivative(WQTotal + Indices::numEq); - if (has_wfrac_variable) { - linSys_.duneD_[seg][seg_upwind][SPres][WFrac] += pressure_equation.derivative(WFrac + Indices::numEq); - } - if (has_gfrac_variable) { - linSys_.duneD_[seg][seg_upwind][SPres][GFrac] += pressure_equation.derivative(GFrac + Indices::numEq); - } - // contribution from the outlet segment const int outlet_segment_index = this->segmentNumberToIndex(this->segmentSet()[seg].outletSegment()); const EvalWell outlet_pressure = getSegmentPressure(outlet_segment_index); - linSys_.resWell_[seg][SPres] -= outlet_pressure.value(); - for (int pv_idx = 0; pv_idx < numWellEq; ++pv_idx) { - linSys_.duneD_[seg][outlet_segment_index][SPres][pv_idx] = -outlet_pressure.derivative(pv_idx + Indices::numEq); - } + const int seg_upwind = upwinding_segments_[seg]; + MultisegmentWellAssemble(baseif_). + assemblePressureEq(seg, seg_upwind, outlet_segment_index, + pressure_equation, outlet_pressure, linSys_); if (this->accelerationalPressureLossConsidered()) { handleAccelerationPressureLoss(seg, well_state);