examples.well_control
Rate control, BHP control, and the well model that connects them.
A well can be told what to do in two ways: give it a rate (rate, ref
Wells.rates; negative to produce) and let its pressure follow, or give it a
pressure (bhp, ref Wells.bhp) and let its rate follow -- as ResSim.wells
accepts either. Both need a well model -- the well index (peaceman_WI,
whence the rw of the record) -- because a well is far smaller than the cell
that holds it, so its cell pressure is not its wellbore pressure.
That distinction is the first thing shown here: the cell pressure is a grid artefact, whereas the bottom-hole pressure obtained from it is not.
The two modes are then contrasted on the same closed reservoir under primary
depletion -- production without injection, which only ct > 0 allows, the
deficit being drawn from storage -- where they behave quite differently:
- at constant rate, the pressure declines linearly, as material balance dictates: $ d\bar{p}/dt = -q / (c_t V_p) $;
- at constant BHP, the rate declines exponentially, $ q ∝ e^{-t/τ} $ with $ τ = c_t V_p / J $, since combining that same material balance with the well model, $ q = J (\bar{p} - p_\mathrm{bh}) $, gives a linear ODE. (The productivity index $ J $ differs from $ WI λ_t $ by the geometry between the well and the average pressure.)
Neither mode alone is how a well is actually run: the industry standard is a
rate target with a BHP limit, i.e. whichever of the two currently binds.
The model does not switch modes natively, but ResSim.well_controls returns
both controls, so an override can switch between them -- lagged by the step
whose pressure it must judge from. That is the third case shown here.
Finally, the two are shown to be one and the same model, seen from either end: prescribing the BHP that the rate-controlled run reported recovers that run exactly (to ~1e-15).
In the figures:
- "diagnostic": the drawdown grows throughout the transient -- the pressure disturbance has not yet reached the closed boundary, so the reservoir behaves as if infinite -- and settles once the flow is boundary-dominated, the pressure profile then frozen in shape and merely subsiding, at about the time $ r^2/η $ that the disturbance needs to reach the boundary (dotted). Refining 32² → 64² moves the producer's cell pressure by a lot (left), while the bottom-hole pressure inferred from it barely moves (right). Only the latter is a property of the well.
- "modes" (left): the rate is flat by construction under rate control, and decays under BHP control -- as a straight line on the log axis, i.e. exponentially, with the analytic slope (dashed). The rate-with-a-limit case traces the former until the limit binds, and joins the latter thereafter (above it, having drained less by then, hence at a higher pressure).
- "modes" (right): the mirror image. The BHP falls linearly under rate control (material balance), and is flat by construction under BHP control -- while the limited well does both in turn, its corner marking the switch.
- "duality": feeding the left run's BHP back in as the control reproduces its rate to machine precision.
1"""Rate control, BHP control, and the well model that connects them. 2 3A well can be told *what to do* in two ways: give it a rate (`rate`, ref 4`Wells.rates`; negative to produce) and let its pressure follow, or give it a 5pressure (`bhp`, ref `Wells.bhp`) and let its rate follow -- as `ResSim.wells` 6accepts either. Both need a **well model** -- the well index (`peaceman_WI`, 7whence the `rw` of the record) -- because a well is far smaller than the cell 8that holds it, so its cell pressure is not its wellbore pressure. 9 10That distinction is the first thing shown here: the cell pressure is a *grid 11artefact*, whereas the bottom-hole pressure obtained from it is not. 12 13The two modes are then contrasted on the same closed reservoir under *primary 14depletion* -- production without injection, which only `ct > 0` allows, the 15deficit being drawn from storage -- where they behave quite differently: 16 17- at constant **rate**, the pressure declines linearly, as material balance 18 dictates: $ d\\bar{p}/dt = -q / (c_t V_p) $; 19- at constant **BHP**, the rate declines *exponentially*, 20 $ q ∝ e^{-t/τ} $ with $ τ = c_t V_p / J $, since combining that same material 21 balance with the well model, $ q = J (\\bar{p} - p_\\mathrm{bh}) $, gives a 22 linear ODE. (The productivity index $ J $ differs from $ WI λ_t $ by the 23 geometry between the well and the average pressure.) 24 25Neither mode alone is how a well is actually run: the industry standard is a 26rate *target* with a BHP *limit*, i.e. whichever of the two currently binds. 27The model does not switch modes natively, but `ResSim.well_controls` returns 28both controls, so an override can switch between them -- lagged by the step 29whose pressure it must judge from. That is the third case shown here. 30 31Finally, the two are shown to be one and the same model, seen from either end: 32prescribing the BHP that the rate-controlled run *reported* recovers that run 33exactly (to ~1e-15). 34 35In the figures: 36 37- "diagnostic": the drawdown grows throughout the *transient* -- the pressure 38 disturbance has not yet reached the closed boundary, so the reservoir behaves 39 as if infinite -- and settles once the flow is *boundary-dominated*, the 40 pressure profile then frozen in shape and merely subsiding, at about the time 41 $ r^2/η $ that the disturbance needs to reach the boundary (dotted). Refining 42 32² → 64² moves the producer's cell pressure by a lot (left), while the 43 bottom-hole pressure inferred from it barely moves (right). Only the latter 44 is a property of the *well*. 45- "modes" (left): the rate is flat by construction under rate control, and 46 decays under BHP control -- as a straight line on the log axis, i.e. 47 exponentially, with the analytic slope (dashed). The rate-with-a-limit case 48 traces the former until the limit binds, and joins the latter thereafter 49 (above it, having drained less by then, hence at a higher pressure). 50- "modes" (right): the mirror image. The BHP falls linearly under rate control 51 (material balance), and is flat by construction under BHP control -- while the 52 limited well does both in turn, its corner marking the switch. 53- "duality": feeding the left run's BHP back in as the control reproduces its 54 rate to machine precision. 55""" 56 57from mpl_tools.place import freshfig 58import numpy as np 59 60from minires import ResSim 61from minires.plotting import show 62 63## Setup 64q = .25 # the rate-controlled rate 65p_bh = .5 # the BHP-controlled pressure 66rw = 1e-3 # well radius 67ct = .1 68dt, nSteps = 2e-3, 150 69tt = dt*np.arange(1, nSteps + 1) 70 71def depleter(N=32, cls=ResSim, **control): 72 """A single producer at the centre of a closed square (cf. `examples.buildup`). 73 74 The `control` is a `rate` and/or a `bhp`; `rw` is what gives it a well model. 75 """ 76 return cls(Lx=1, Ly=1, Nx=N, Ny=N, ct=ct, 77 wells=[dict(xy=[.5, .5], rw=rw, **control)]) 78 79class Limited(ResSim): 80 """Rate control with a BHP limit, by overriding `ResSim.well_controls`. 81 82 The rate target is held for as long as it can be delivered without drawing 83 the well below `p_bh`; thereafter the well switches to BHP control at it. 84 The switch is judged from the *previous* step's pressure, since the new one 85 is not yet known (indeed it depends on the choice) -- so the limit is 86 breached for the one step in which it comes to bind. 87 """ 88 89 def well_controls(self, S, P, k): 90 ctrl = super().well_controls(S, P, k) 91 if P is None: 92 return ctrl # nothing to switch on 93 would = self.bhp(S, P, ctrl["rates"]) # if rate-controlled 94 ctrl["bhp"] = np.where(would < p_bh, p_bh, np.nan) 95 return ctrl 96 97def run(model): 98 SS, PP = model.sim(dt, nSteps, np.zeros(model.Nxy), 99 P0=np.ones(model.Nxy), pbar=False) 100 assert SS.max() == 0, "No water is injected, so none should appear." 101 return PP 102 103## Plot: the diagnostic -- cell pressure is a grid artefact, bottom-hole is not 104fig, (ax1, ax2) = freshfig("Well control -- diagnostic", ncols=2, figsize=(10, 4), 105 sharey=True) 106for N in [32, 64]: 107 model = depleter(N, rate=-q) 108 PP = run(model) 109 pbar = PP.mean(axis=1) 110 ax1.plot(tt, (pbar[1:] - PP[1:, model.xy2ind(*model.wells.xy[0])]), label=f"{N}²") 111 ax2.plot(tt, (pbar[1:] - model.wells.actual_bhp[0]), label=f"{N}²") 112ax1.set(title="Cell drawdown, $\\bar{p} - p_\\mathrm{cell}$", 113 xlabel="Time", ylabel="$\\Delta p$", xscale="log") 114ax2.set(title="Bottom-hole drawdown, $\\bar{p} - p_\\mathrm{bh}$", xlabel="Time", 115 xscale="log") 116eta = 1/ct # Diffusivity (K = por = 1, and λ = 1, there being no water) 117for ax in (ax1, ax2): 118 ax.axvline(.5**2/eta, c="k", ls=":", lw=1, label="$r^2/η$, $r$ = ½ (to the boundary)") 119 ax.legend(title="Grid", fontsize="small") 120fig.tight_layout() 121 122## Simulate: the same reservoir, under either mode of control 123by_rate = depleter(rate=-q) 124PP_rate = run(by_rate) 125by_bhp = depleter(bhp=p_bh) # NB: no `rate` given 126PP_bhp = run(by_bhp) 127limited = depleter(cls=Limited, rate=-q) # ... rate, but limited by p_bh 128run(limited) # (its rate/BHP is the interest) 129 130# Production, i.e. the negated (signed) rates 131prod_rate, prod_bhp, prod_lim = [-m.wells.actual_rates[0] 132 for m in [by_rate, by_bhp, limited]] 133 134# Analytic decline: q = J (pbar - p_bh) with material balance ct Vp dpbar/dt = -q 135Vp = by_bhp.h2 * by_bhp.por.sum() 136J = prod_bhp[-1] / (PP_bhp[-1].mean() - p_bh) 137tau = ct*Vp/J 138 139## Plot: rate and BHP, under either mode 140fig, (ax1, ax2) = freshfig("Well control -- modes", ncols=2, figsize=(10, 4)) 141 142ax1.plot(tt, prod_rate, label="Rate-controlled") 143ax1.plot(tt, prod_bhp, label="BHP-controlled") 144ax1.plot(tt, prod_lim, ":", lw=2, label="Rate, limited") 145ax1.plot(tt, prod_bhp[-1]*np.exp((tt[-1] - tt)/tau), "k--", 146 lw=1, label=f"$\\propto e^{{-t/\\tau}}$, $\\tau = c_t V_p / J$ = {tau:.3f}") 147ax1.set(title="Production rate", xlabel="Time", ylabel="q", yscale="log") 148ax1.legend(fontsize="small") 149 150ax2.plot(tt, by_rate.wells.actual_bhp[0], label="Rate-controlled") 151ax2.plot(tt, by_bhp .wells.actual_bhp[0], label="BHP-controlled") 152ax2.plot(tt, limited.wells.actual_bhp[0], ":", lw=2, label="Rate, limited") 153ax2.plot(tt, 1 - q*tt/(ct*Vp) - (PP_rate[-1].mean() - by_rate.wells.actual_bhp[0, -1]), 154 "k--", lw=1, label="$p_0 - qt/(c_t V_p) - \\Delta p$") 155ax2.set(title="Bottom-hole pressure", xlabel="Time", ylabel="$p_\\mathrm{bh}$") 156ax2.legend(fontsize="small") 157fig.tight_layout() 158 159## The duality: prescribe the BHP that the rate-controlled run reported 160replay = depleter(bhp=by_rate.wells.actual_bhp[0]) 161PP_replay = run(replay) 162err_P = np.abs(PP_replay - PP_rate).max() 163err_q = np.abs(replay.wells.actual_rates + q).max() 164assert err_P < 1e-12 and err_q < 1e-12, "The two controls are not each other's inverse!" 165 166fig, ax = freshfig("Well control -- duality", figsize=(6, 4)) 167ax.plot(tt, prod_rate, lw=4, alpha=.4, label="Rate-controlled: $q$") 168ax.plot(tt, -replay.wells.actual_rates[0], "k--", lw=1, 169 label="BHP-controlled by its own reported $p_\\mathrm{bh}$") 170ax.set(title=f"The same well, controlled from either end (max err {err_q:.0e})", 171 xlabel="Time", ylabel="q", ylim=(0, 2*q)) 172ax.legend() 173fig.tight_layout() 174 175# Regression values, checked by `tests/test_examples.py`. 176# NB: the production rates are negated, preserving the pre-v0.3 references. 177__digest__ = dict(rate_of_bhp_ctrl = prod_bhp, 178 bhp_of_rate_ctrl = by_rate.wells.actual_bhp[0], 179 p_last_bhp_ctrl = PP_bhp[-1], 180 rate_of_limited = prod_lim, 181 bhp_of_limited = limited.wells.actual_bhp[0]) 182 183if __name__ == "__main__": 184 show()