Implement injection networks - #7345
Closed
atgeirr wants to merge 34 commits into
Closed
Conversation
A network branch lookup clamped the rate and upstream pressure to nothing: VFPHelpers::findInterpData extrapolates linearly past the axis ends, and the flow-line tables are zero-filled where the line cannot deliver the rate, so an injection network with rates beyond the table axis produced node pressures of -260 bar, and a THP below the axis start gave 0. Both were then handed to the wells as THP limits. The branch calculator now clamps the lookup point to the table axes (rates scaled uniformly so prod WFR/GFR are kept) and treats a result <= 1 atm as 'no solution'. Such nodes, and their descendants, are reported by NetworkPressureComputation::invalidNodes(); updatePressures() keeps their previous pressure, counts the network as unbalanced (max update bound), warns once per report step, and does not push a dynamic THP to wells under them. Also make updateActiveStateImpl accumulate instead of clearing active_ on an inactive domain (unreachable today, but wrong). Tests: gas_injection_rate_beyond_flow_axis, gas_injection_zero_cell_region and gas_injection_thp_below_axis in test_networkpressure.cpp; on 0996672 they give -213 bar, 0 bar and 28.9 bar respectively. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The explicit node-pressure iteration P <- P + w (P_computed - P) cannot converge an injection network: with BHP pinned by the reservoir the injectors' THP-mode rate is very steep in THP (~6e4 sm3/d/bar on GNETINJE_GAS-01) and the flow-line pressure falls ~0.5 bar per 1000 sm3/d, a loop gain of ~30, so with w=0.1 the iteration falls into a period-2 limit cycle (leaf inflow 0 <-> 1.5e6, pressure 500 <-> 0). Wells stopped at the wrong moment were then shut for good. NodePressureUpdater (new header) replaces the update per node: keep a sign-change bracket on r(P) = P_computed(P) - P and take Illinois regula-falsi steps inside it; before a bracket exists take a secant step clipped to [P, P_computed], capped at max(25% P, 10 bar), and, when the leaf's wells are on group/rate control (flat response down to their own THP), go straight to just below that kink. Stale bracket ends (from wells that were momentarily stopped, or from the other leaves moving) are dropped when contradicted, when the bracket collapses, or when their Illinois weight is halved away. Convergence is judged on the remaining pressure uncertainty (bracket width), not on the residual, which is amplified by the well response. Enabled by default; the old damped update is the fallback and is kept behind --network-pressure-update-secant=false. Well side: the injector's dynamic THP is clamped to its VFPINJ THP axis instead of being silently left stale; injectors get the previous step's node pressure at report-step start like producers (starting at the WCONINJE THP made them inject at their rate limit and starve the other leaves); the 'Option B' potentials refresh is removed (it wrote well_potentials = surface_rates, so a well that had been stopped restarted every THP solve at zero rate and was stopped again). Result on opm-tests/network/GNETINJE_GAS-01 and _WAT-01: 0 unconverged network messages, no wells shut, and rates/THP/control modes match the E100 reference (opm-tests/eclref) to within 1% / 0.2 bar at every report step. Tests: NodePressureUpdaterTests in test_networkpressure.cpp (stiff response converges where the damped update does not, stale end dropped, plateau floor, step cap, invalid evaluation moves down). Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The physical-shut path is meant for wells that cannot operate under their own THP/BHP limits. An injector whose THP is a network node pressure may be unable to inject at this iteration's pressure and yet inject fine once the network is balanced (the pressure depends on the other wells' rates); shutting it makes the network converge to a wrong state with the remaining wells carrying the whole group target. Keep such wells stopped for the step; the local solve reopens them when the pressure allows. Injectors only for now. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
GAS and WATER injection networks were written into the same leaf-rate map, so a group that is a leaf of both would see the second network's rates, and the injection phase threaded into NetworkPressureComputation was ignored. Key the map by (Phase, group) and use the phase in the calculator; the mock in test_networkpressure.cpp now returns only the network's own phase and its phase-less lookup returns nothing, so the gas/water cases fail without this. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
assignNodeAndBranchValues() only exported the production network. Fill the gas/water injection maps of data::GroupAndNetworkValues (opm-common c4b9da0a4) from the per-domain node pressures and branch data as well. On GNETINJE_GAS-01, GPRG now matches the E100 reference: 340.0 / 209.4 / 209.4 / 204.3 / 204.3 bar (PLAT-A / M5S / G1 / M5N / F1) at day 31. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
… default --network-pressure-update-secant becomes a string: injection (default), all, none. With the bracketing update also on production networks NETWORK-01 gives a slightly different (better converged) solution than the reference (BPR 0.04 bar at one step); keep production networks on the damped update until the references are regenerated. Production-network results are now identical to the pre-change scheme. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
…upported Take GNETINJE off the unsupported-keyword list (it made the decks fail at the default parsing strictness) and register the two opm-tests injection network decks as regression tests (--solver-max-time-step-in-days=1). The references are generated with this branch; both cases agree with the E100 runs in opm-tests/eclref within 2% on node pressures and 3% on well rates at every report step (tools/plot_injection_network.py --check). Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
An injection group target (GCONINJE VREP/REIN) is a function of the produced voidage, so within a network sub-iteration the injectors must see this iteration's production, not the previous one's. The well loop ran in deck order, so that was a matter of luck. Producers are now re-solved first, then the injectors. The group data is refreshed in between only when a deck has both a production and an injection network -- the only case where the producers re-solved here can feed an injection target in the same sub-iteration -- so nothing else pays for the extra update. details::activeNetworks() already returns the production network first, so the pressure computation was in the right order already. No deck in opm-tests has both network kinds, so this is a correctness guard rather than a fix for an observed failure: NETWORK-01, NETWORK-01_STANDARD, NETWORK-01-WTEST, GNETINJE_GAS-01 and GNETINJE_WAT-01 are all bit-identical (compareECL -t SMRY at 1e-9) before and after. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
…t help yet)
An experiment against the per-node bracketing update: Anderson uses the last few
(P, G(P)) pairs of the whole pressure vector of a network, so it sees the coupling
between nodes that a per-node update cannot, without needing any derivative.
Kept deliberately isolated: NetworkAndersonAcceleration.hpp is self-contained (no
OPM dependencies, own tiny least-squares solve), there is one call site in
updatePressures(), and it is off unless
--network-pressure-update-acceleration=anderson is given
(--network-anderson-depth sets the history length, default 4). With the default
the results are bit-identical to the previous commit on all five network decks.
Measured on the opm-tests decks, against the E100 reference:
GNETINJE_WAT-01 bracketing: 3 unconverged steps anderson(4): 0
GNETINJE_GAS-01 bracketing: 0 unconverged steps anderson(4): 18
GNETINJE_GAS-01 WGIR:G-3H vs E100, over all reported steps
bracketing: median 0.24 %, p90 3.4 %, 86 % within 2 %
anderson(4): median 0.02 %, p90 100 %, 63 % within 2 %
So it is better where the response is smooth and much worse where it is not: the
gas case fails the report-step comparison with Eclipse that the bracketing update
passes. That is the expected failure of unsafeguarded extrapolation across the
kinks where wells switch control -- Anderson replaces the bracket, the step caps
and the plateau rule that the per-node updater relies on. It should not be
enabled as it stands; a variant that keeps those safeguards (limit the Anderson
step to the same bound, reject it when the residual grows) is the obvious next
thing to try.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
… yet) Balance the injection networks against each injector's linearised rate response instead of re-solving the well equations, so a residual evaluation costs a VFP lookup. --network-well-proxy=ipr, default none; nothing changes when it is off. Two things had to be fixed before this could even be measured. The explicit IPR (ipr_a_/ipr_b_, WellInterface::updateIPR) is identically zero for injectors: its crossflow guard drops every connection with pressure_diff = p_r - h_perf > 0, which for an injector is the ordinary case. Nothing noticed because the only consumers are the producer operability checks. The implicit IPR is the usable one -- it differentiates the converged well equation through the well Jacobian -- but updateIPRImplicit hard-coded the producer control swap, while the rhs picks out the control equation and so requires that equation to be the bhp one. Made symmetric here; the producer path is unchanged. With that, the response is exact at the point of linearisation (it reproduces each injector's solved rate to four figures), and the proxy is still not good enough: on GNETINJE_GAS-01 and _WAT-01 it takes unconverged network steps from 0/3 to 0/1 but the deviation from the E100 reference goes from 44/1 violations to 620/646, because from about two thirds through both runs it drives every injector to zero rate. dq/dbhp is large enough here that the intersection bhp moves 0.04 bar for a 292 bar change in applied thp, so any node-pressure error lands on the steep part, the linear rate goes negative and the well shuts. A Newton on this Jacobian will need globalisation and the well's own rate/group limits, not a full step. Kept because the linearisation itself is sound and is the building block for that Newton. See injection_network_findings.md. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
tests/test_networksolve.cpp extracts the GNETINJE_GAS-01 network with the wells replaced by their inflow performance plus the control logic, so solution methods can be compared without a reservoir or a well solve. The three VFPINJ tables are the deck's, verbatim; the wells are calibrated to the Eclipse 100 operating point at day 31, and the bench then reproduces it to 0.1 bar. Table 9999 means "no table", so G1 carries M5S's pressure and F1 carries M5N's and the whole problem has two unknowns. That is small enough to look at a Jacobian by hand. From the wells' WCONINJE THP, to 0.01 bar: damped (omega 0.1) FAILED (limit cycle -- the original branch's method) bracketing (shipped) 42 anderson (depth 4) 21 newton (full step) FAILED (overshoots off the plateau to 450/506 bar) newton + line search 14 and over dq/dbhp = 1e4 .. 1e6 sm3/d/bar, bracketing 29/42/30/38 against newton+ls 10/14/21/22. So a Newton on this problem is worth doing and must be globalised: the full step diverges on its own, and the same line search that fixes it also keeps the advantage as the wells stiffen. Note for anyone extending this: VFPInjProperties::addTable stores a reference_wrapper, so the tables have to outlive the properties object and must not be moved. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The Newton direction was already there; what was missing was any measurement of what to do with it. Four strategies now share one Newton driver -- full step, step cap, backtracking line search (plain or Armijo), trust region -- so the only thing that varies between them is how the step is accepted. Iteration count from one good start says little about a globalisation, so the bench sweeps a 23x23 grid of starting pressures across the tables' THP axis and counts how many reach the right answer: bracketing (shipped) 529/529 mean 31 iterations newton, full step 3/529 mean 5 newton, capped step 54/529 mean 8 newton, line search 529/529 mean 12 newton, trust region 529/529 mean 13 Two things worth having measured. Capping the step -- which is what --network-max-pressure-update-in-bars already does -- is not a globalisation: it recovers from a tenth of the space. And a real one costs nothing in robustness, so the choice between the line search and the trust region is a choice about code, and the line search is the smaller piece of code. The trust region needed its collapse case handled to get there: when the Jacobian is taken across a control switch the radius halves away to nothing, and returning the unchanged point stalls the solver. Taking the smallest step and reopening the region fixes it (83% -> 100%). Also covers the group-target branch of the well control logic, which the earlier tests left dead. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
… one The bench solved the network in the node pressures alone, with every rate recovered from them by an inner solve. This adds the same network without the eliminations -- node pressures, branch rates and each well's (rate, bhp) as unknowns, with the eliminations restated as node balances, branch drops, inflow performance and one control equation per well -- and lets the same Newton solve either. 16 unknowns against 2. Over the grid of 529 starting pressures, at the measured stiffness: eliminated, plain newton 3/529 mean 5 iterations eliminated, line search 529/529 mean 12 full, plain newton 529/529 mean 9 full, line search 529/529 mean 9 The eliminated residual needs globalising because the control clamps put kinks in it. The full system holds its controls fixed while the step is taken, so there are no kinks and an unglobalised Newton is already fully robust -- and faster. That is the case for carrying the rates. Three things had to be right for that to hold, and each is a finding in its own right: Controls must be chosen by most-restrictive-wins, matching the clamp the eliminated form applies. A fixed priority chain makes the active set chatter and the Newton never terminates. Globalisation must not veto an active-set change. The residual jumps when a control switches and that is not a failure to make progress; letting the line search reject it stalls the switch instead of resolving it (363/529 -> 529/529). The table limits want to be bounds on the unknowns, not a clamp on the lookup. Outside the tables' box the zero-filled cells extrapolate away and the full system has a spurious root there -- at the softest wells a plain Newton reaches it from a third of the grid, ending with every well at its rate limit and node pressures of -683 and -10975 bar. Clamping the lookup to the axes, which is what the simulator's pressure computation does, removes that root and replaces it with a flat residual the Newton cannot descend (383/529 -> 23/529). Neither is right; the limits belong in the active set beside the well controls. table_bounds_want_to_be_constraints pins both halves down. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
… build Three things, of which the first two were the open questions. The table limits belong on the unknowns. Clamping the lookups to the axes, which is what the simulator's pressure computation does, flattens the residual and leaves a Newton nothing to descend (22/529 starts). Leaving the tables to extrapolate keeps the derivatives but admits a root outside them (412/529). Holding the branch flows inside the box does best (421/529), and only by projecting the offending components -- a scalar fraction-to-boundary lets one binding rate throttle the pressure updates as well, which is worse than doing nothing (299/529). The group target is now an equation with a multiplier, so wells the group does not bind stay on their own controls. It converges in 4 iterations and lands on the target exactly, with the network at higher pressure than it runs free, which is what curtailment should do. One trap: the activation test has to include equality, because at the solution a well's rate is its share exactly and a strict test flips the control every iteration. The case is now data rather than code. A NetworkCase is nodes with a parent and a table, wells hanging off nodes, a terminal pressure and an optional group target; gnetinjeGas() is one instance of it, assembled the way a deck reader would -- topology from GRUPTREE/GNETINJE, tables from the VFPINJ includes, limits from WCONINJE, and the calibration point from the reference summary. Both problems and all the solvers work off that, with no node names or unknown indices written into them, so a second case is a builder function rather than a rewrite. The refactor reproduces every previously measured number. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
A case needs a reference operating point per well -- the rate it takes and the pressure at its node -- because that is what its IPR is a linearisation about. Both ways of supplying one are now there. Reference::set() writes it down. That is what the built-in GNETINJE_GAS-01 case uses, so the bench still needs no files and runs everywhere, and it is also how you would pose a case no reference run covers. Reference::fromSummary() reads it out of a reference run at a given time, paired with fromDeck(), which takes the topology from GRUPTREE/GNETINJE, the tables from the VFPINJ keywords and the limits from WCONINJE. It reads the keywords off the parsed Deck rather than building a Schedule, so it stays usable from a test. built_from_deck_matches_the_builtin_case checks the two agree, node for node and well for well, and that the deck-built case solves to the same place. It skips when opm-tests is absent. Two things that cost time and are worth knowing: WCONINJE's RATE is a UDA item whose unit depends on the injected phase, so the parser leaves it dimensionless and getSI() hands back the raw deck number. Taken at face value the rate limit came out 86400x too large, which silently removes the rate control -- the case still looked right and simply stopped converging. The operating point has to be read at a report step. Off the ministep vectors the same nominal time gives a slightly different state, which is a different case. Correcting the built-in reference against the summary at day 31 moved it 0.6% on the G wells, and the bench now lands on (209.405, 204.240) against E100's (209.409, 204.244) -- previously (209.30, 204.19). Rates are carried as an (aqua, liquid, vapour) triple rather than one number, which is what a production network will need: there the branch rate splits and VFPPROD wants WFR and GFR fractions and an ALQ, and those fractions are extra unknowns per branch with their own mixing equations at the nodes. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Until now the bench solved one hand-picked operating point. This runs the full formulation at every report step of both GNETINJE decks, calibrated at each step to what the Eclipse reference says the wells were doing, with the field target applied where the reference has them on group control. GNETINJE_GAS-01 4/4 report steps, worst 0.05 bar, 4-16 iterations GNETINJE_WAT-01 6/6 report steps, worst 0.39 bar, 4-8 iterations Both group-controlled steps land on the reference exactly. Those are the steps where the method we ship is weakest -- day 91 of the gas case is where it spends two reported steps at the validity floor. Getting there needed two corrections, both of which matter beyond the bench. A well must be calibrated against its own bhp, not the pressure at its node. The two agree only while it is on THP control. At day 91 the wells sit 162 bar below their node (WTHP 202 against GPRG 365), and calibrating there against the node pressure builds a well that does not exist -- which is why that step would not converge at all before. The stiffness knob only transfers between cases as a fraction of the well's own rate. 6e4 sm3/d/bar is reasonable for a gas injector taking 5e5 sm3/d and meaningless for a water injector taking 700; carrying it across is what made three of the water steps diverge. setRelativeStiffness(0.12) reproduces the gas figure and works for both. fromDeck also had to learn that GNETINJE spells the phase WAT where WCONINJE spells it WATER. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
opm/simulators/wells/NetworkSystem.hpp now holds the injection network solved simultaneously in its pressures and its rates: the branch drops, node balances, inflow performance and per-well control equation, an optional group multiplier, the active-set control selection, the branch-flow bounds and the Newton. The bench is now a caller of it rather than a separate implementation. All it still owns is where a well's inflow performance comes from -- a linearisation about a reference operating point, where the simulator will use the well Jacobian -- plus the eliminated formulation and the fixed-point methods it is compared against. Every measured number is unchanged, which is the check that the extraction is faithful. setClampToAxes is on the system only so the comparison that rules it out can still be made, and says so. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
…newton
The simulator now fills in the same NetworkSolve::System the bench does. The
only thing that differs is where a well's rate response comes from: the bench
linearises about a reference operating point, the simulator reads the implicit
IPR out of the well state. Everything downstream -- the equations, the
active-set control selection, the branch-flow bounds, the Newton -- is shared.
Default is unchanged (fixedpoint), and both GNETINJE references still compare
clean. With --network-solver=newton, both decks at 1-day steps:
unconverged network steps deviations from E100
gas fixedpoint 0 44
gas newton 0 10
water fixedpoint 3 4
water newton 0 4
so the gas case lands four times closer to Eclipse and the water case stops
taking unconverged steps, in slightly fewer Newton iterations either way.
The node pressures come back already at the fixed point, so the relaxation that
follows sees no imbalance and stops; the branch data from the ordinary
evaluation is kept for the output. Anything the solve cannot handle -- more than
one root, a well with no usable inflow performance, a network that does not
converge -- returns nullopt and leaves the fixed point in charge.
The group multiplier the system can carry is not used here. The simulator's own
group machinery has already decided each well's share by the time the network is
solved, so a well on GRUP enters as one held at that rate.
updateIPRImplicit is refreshed for injectors before the solve; the well solve
only maintains it for producers.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
--network-analytic-jacobian assembles the Jacobian from the derivatives the VFP
interpolation already computes and throws away, instead of differencing the
residual n+1 times. Everything but the two table lookups is constant. Checked
against the difference quotient in four states -- the ordinary THP one, a
converged one, one with the group holding the wells, and one driven onto the
bhp and rate limits -- and agreeing to about 1e-7 relative. On both decks it
gives the same answer as differencing, which is the point: it is a cost, not an
accuracy, change.
--network-group-control hands the network the group's injection total and lets
it place the split, so a well that runs into its own limit is taken up by the
others rather than the total quietly dropping. On these two decks it changes
nothing, because the group machinery has already made the split and no well is
near a limit, so the multiplier has nothing to redistribute. It is off by
default and unproven; a case where a group-controlled well hits its bhp limit is
what would show whether it earns its place.
Both decks, 1-day steps, against the Eclipse reference:
unconverged deviations
gas fixedpoint 0 44
gas newton 0 10
gas newton + analytic 0 10
gas newton + analytic + grp 0 10
water fixedpoint 3 4
water newton 0 4
water newton + analytic 0 4
water newton + analytic + grp 0 4
Also fixes two things in the previous commit's well setup. A group-controlled
well is held at the rate the group gave it, which had been lost; without it the
gas case went from 10 deviations to 252. And a well with neither a rate nor a
target now takes the whole network back to the fixed point instead of entering
with an invented limit -- std::max(rate, 1) is 86400 sm3/d in SI, which is no
limit at all for a water injector taking 700, and cost the water case 4
deviations against 240.
The solve is now started from the wells' current rates where the caller knows
them, so the first control selection is made on the state the well is actually
in.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
A well's share of a group total should be proportional to what it can inject at the pressure the network gives it. The previous version used its current rate, which is the split being decided -- so the allocation reproduced whatever it already was, which is why --network-group-control changed nothing. System::thpPotential() meets a well's inflow performance with its tubing curve at a given node pressure. The guides are refreshed from it each iteration and have to settle before the solve is called converged, the same way the controls do. It fixes the node pressures at the gas case's day 61 -- GPRG:M5N/F1 were 2.6 % high and are now inside tolerance, taking the gas deviations from 10 to 8 -- and leaves the well rates bit-identical. That last part is the useful half of the result. Under group control the network writes nothing to a well but a dynamic THP limit, and for a group-held well that limit does not bind, so the network cannot move the split however good its guide rates are. The rates come from the group allocation, which runs before the network and does not see it. Closing that needs the allocation to become part of the same solve, not a better guide handed to it. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The simultaneous solve read the rank-local well container and communicated nothing, so under MPI each rank would assemble a network out of its own wells and get its own answer. It was latent because the solver is opt-in and had only been run serially. The relaxed computation stays right by working from group and node quantities, which are the same everywhere by the time it runs. This does the same thing one level down. The well list and everything static about a well -- its node, its vfp table, its bhp and rate limits -- now come from the schedule, which is replicated. Only what a well is currently doing is rank-local: its inflow performance, its rate and its control mode are contributed once by the rank that owns it and summed, after which every rank holds the same numbers and reaches the same decisions, including the decision to give up and leave it to the fixed point. Contributed by the owner rather than by every holder because a distributed well appears on several ranks carrying the same values; summing those would count it twice. On GNETINJE_GAS-01 the node pressures now agree across 1, 2 and 4 ranks to about 1e-5 relative, and are identical to six figures on the group-controlled steps. The deviation count against the Eclipse reference is 8 on all three, and 4 on the water case serial and on two ranks. What is left between serial and parallel is a connection total on a shut well that the default solver differs on in exactly the same way, so it is not this. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Every way out of newtonNodePressures() returned nullopt without a word, so there was no way to tell whether the simultaneous solve had run or the relaxed update had quietly done the work. Each exit now names its cause -- more than one root, no terminal pressure, the phase not active, no injectors on the network, a well with no usable inflow performance, a well with neither a rate nor a target, or simply not converging -- and a successful solve reports its iteration count. Debug level, because this runs once per network sub-iteration; it lands in the .DBG file next to the network trace. It paid for itself immediately. On the two decks the solve is used for most of the run but not all of it: gas 282 solved against 21 handed back (6.9 %), water 344 against 37 (9.7 %), nearly all of them non-convergence in 50 iterations. Neither run's answer changes -- the relaxed update covers those steps, which is what the fallback is for -- but a tenth of the water case was taking a path nobody could see. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
A non-converged solve now carries its residual, whether a control was still switching, whether the group shares were still moving, and one letter per well per iteration for the last eight iterations. The caller prints all of it. That was enough to identify what the fallbacks actually are. Every one of them is an active set that will not settle, and the trace shows a clean period-2 cycle: water TTTT GGGG TTTT GGGG TTTT GGGG TTTT GGGG gas TTRR TTGG TTRR TTGG and RRTR GGTG RRTR GGTG so the wells are flipping between group and thp control, not converging slowly. The share a well is tested against moves with the guides and with the multiplier while its rate is chasing that share, and the two never meet. Two candidate fixes were measured and neither is right yet, so both are written down rather than shipped. Making the bhp and rate activation tests inclusive, as the group one already is, looks correct -- a well at a limit sits exactly on it -- and is not: start() clamps the opening rate to the limit, so an inclusive test latches every well onto rate control at the first iteration and holds it there, which sends the bench to the out-of-table root at -683 bar. The tests stay strict, with a comment saying why. Deciding group membership from the group rather than the well -- it binds when the wells could between them take more than the target -- does cure the cycling, taking the water case from 37 fallbacks to 2 and the gas case from 21 to 7. It also moves the gas answer away from the reference, 8 deviations to 26. Not a trade worth making blind, so it is described in the code where the next person will find it. Nothing here changes an answer: gas 8 deviations and water 4, as before. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The simulator refreshes each well's share of a group total from what it can inject at the network pressure. The bench never turned that on, so the one code path that fails in the simulator was the one path with no coverage. It has coverage now, in the simulator's own group configuration -- the target is the total the wells are already injecting, the guides are their current rates, so the multiplier starts at one and every well starts exactly on its share, right on the activation boundary. It converges, in three iterations, on the same answer as holding the guides fixed. So the GRUP/THP cycling is not the guide logic by itself. The half the bench cannot supply is the inflow performance: the simulator's comes from the well Jacobian at a start-up or post-control-change state, and every one of its failures falls at exactly those. Reproducing one needs the failing system dumped and replayed here, which is the next thing to build rather than the next thing to guess. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
--network-dump-failures=<prefix> writes every network system that fails to converge: nodes, wells with the inflow performance the well Jacobian gave them, limits, group target and starting pressures. Point OPM_NETWORK_DUMP at the directory and tests/test_networksolve.cpp replays each one against the same tables, with no simulator. The failures reproduce exactly, cycle for cycle, in about a second. It found its first bug immediately. The group residual summed every well's rate while the target was accumulated only over the wells the group was actually holding, so the group was being asked to account for rates it does not control -- a constraint nothing can satisfy, which the active set then thrashes against instead of settling. 26 of the 146 captured systems converge once the residual counts only what the group holds. Guide rates are now taken once per solve rather than refreshed each iteration. That is what they are in the simulator -- explicit, set per timestep -- and refreshing them inside the Newton makes each well's share a moving target while its rate is chasing it. It takes the water case with group control from 37 solves handed back to 2. The solve also opens just inside a well's limits rather than exactly on them, so the first control selection is not made on a rate that clamping put there. No answer changes: gas 8 deviations from the reference with group control and 10 without, water 4 either way, both references still comparing clean by default. One thing measured and not taken: making the bhp and rate activation tests inclusive, as the group one is, takes the gas case from 70 solves handed back to 4 and the water case from 39 to 2 -- and sends the bench to the out-of-table root from two thirds of its starting points. The simulator and the bench disagree about it, which is worth understanding before either is believed. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Member
Author
|
Closing in favour of #7348 which is a revised version. |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
This is not ready for review or user testing yet, but is made a draft PR to ease developer testing and improvement.