Skip to content

Quarter-five-spot example

The quarter-five-spot is a standard test problem that simulates 1/4 of the five spot well pattern by assuming axial symmetry. The problem contains an injector in one corner and the producer in the opposing corner, with a significant volume of fluids injected into the domain.

julia
using JutulDarcy, Jutul
nx = 50;

Setup

We define a function that, for a given porosity field, computes a solution with an estimated permeability field. For assumptions and derivation of the specific form of the Kozeny-Carman relation used in this example, see Lie, Knut-Andreas. An introduction to reservoir simulation using MATLAB/GNU Octave: User guide for the MATLAB Reservoir Simulation Toolbox (MRST). Cambridge University Press, 2019, Section 2.5.2

julia
function perm_kozeny_carman(Φ)
    return ((Φ^3)*(1e-5)^2)/(0.81*72*(1-Φ)^2);
end

function simulate_qfs(porosity = 0.2)
    Dx = 1000.0
    Dz = 10.0
    Darcy = 9.869232667160130e-13
    Darcy, bar, kg, meter, Kelvin, day, sec = si_units(:darcy, :bar, :kilogram, :meter, :Kelvin, :day, :second)

    mesh = CartesianMesh((nx, nx, 1), (Dx, Dx, Dz))
    K = perm_kozeny_carman.(porosity)
    domain = reservoir_domain(mesh, permeability = K, porosity = porosity)
    Inj = setup_vertical_well(domain, 1, 1, name = :Injector);
    Prod = setup_vertical_well(domain, nx, nx, name = :Producer);
    phases = (LiquidPhase(), VaporPhase())
    rhoLS = 1000.0*kg/meter^3
    rhoGS = 700.0*kg/meter^3
    rhoS = [rhoLS, rhoGS]
    sys = ImmiscibleSystem(phases, reference_densities = rhoS)
    model, parameters = setup_reservoir_model(domain, sys, wells = [Inj, Prod])
    c = [1e-6/bar, 1e-6/bar]
    ρ = ConstantCompressibilityDensities(p_ref = 150*bar, density_ref = rhoS, compressibility = c)
    kr = BrooksCoreyRelativePermeabilities(sys, [2.0, 2.0])
    replace_variables!(model, PhaseMassDensities = ρ, RelativePermeabilities = kr);

    state0 = setup_reservoir_state(model, Pressure = 150*bar, Saturations = [1.0, 0.0])
    dt = repeat([30.0]*day, 12*10)
    dt = vcat([0.1, 1.0, 10.0], dt)
    inj_rate = Dx*Dx*Dz*0.2/sum(dt) # 1 PVI if average porosity is 0.2

    rate_target = TotalRateTarget(inj_rate)
    I_ctrl = InjectorControl(rate_target, [0.0, 1.0], density = rhoGS)
    bhp_target = BottomHolePressureTarget(50*bar)
    P_ctrl = ProducerControl(bhp_target)
    controls = Dict()
    controls[:Injector] = I_ctrl
    controls[:Producer] = P_ctrl
    forces = setup_reservoir_forces(model, control = controls)
    return simulate_reservoir(state0, model, dt, parameters = parameters, forces = forces)
end
simulate_qfs (generic function with 2 methods)

Simulate base case

This will give the solution with uniform porosity of 0.2.

julia
ws, states, report_time = simulate_qfs();
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────╮
 Iteration type   Avg/step   Avg/ministep     Total 
 123 steps  141 ministeps  (wasted) 
├────────────────┼───────────┼───────────────┼──────────┤
 Newton         │    3.3252 │       2.90071 │  409 (0) │
 Linearization  │   4.47154 │       3.90071 │  550 (0) │
 Linear solver  │   10.2358 │       8.92908 │ 1259 (0) │
 Precond apply  │   20.4715 │       17.8582 │ 2518 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2460 │     5.57 % │ 0.1006 │
 Equations     │ 0.1916 │     5.83 % │ 0.1054 │
 Assembly      │ 0.2560 │     7.80 % │ 0.1408 │
 Linear solve  │ 0.2870 │     6.50 % │ 0.1174 │
 Linear setup  │ 1.8234 │    41.29 % │ 0.7458 │
 Precond apply │ 0.1965 │    27.39 % │ 0.4947 │
 Update        │ 0.0761 │     1.72 % │ 0.0311 │
 Convergence   │ 0.0665 │     2.02 % │ 0.0366 │
 Input/Output  │ 0.0333 │     0.26 % │ 0.0047 │
 Other         │ 0.0717 │     1.62 % │ 0.0293 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 4.4163 │   100.00 % │ 1.8063 │
╰───────────────┴────────┴────────────┴────────╯

Plot the solution of the base case

We observe a radial flow pattern initially, before coning occurs near the producer well once the fluid has reached the opposite corner. The uniform permeability and porosity gives axial symmetry at x=y.

julia
using GLMakie
to_2d(x) = reshape(vec(x), nx, nx)
get_sat(state) = to_2d(state[:Saturations][2, :])
nt = length(report_time)
fig = Figure()
h = nothing
ax = Axis(fig[1, 1])
h = contourf!(ax, get_sat(states[nt÷3]))
ax = Axis(fig[1, 2])
h = contourf!(ax, get_sat(states[nt]))
Colorbar(fig[1, end+1], h)
fig

Create 10 realizations

We create a small set of realizations of the same model, with porosity that is uniformly varying between 0.05 and 0.3. This is not especially sophisticated geostatistics - for a more realistic approach, take a look at GeoStats.jl. The main idea is to get significantly different flow patterns as the porosity and permeability changes.

julia
N = 10
saturations = []
wells = []
report_step = nt
for i = 1:N
    poro = 0.05 .+ 0.25*rand(Float64, (nx*nx))
    ws_i, states_i, rt = simulate_qfs(poro)
    push!(wells, ws_i)
    push!(saturations, get_sat(states_i[report_step]))
end
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬─────────────╮
 Iteration type   Avg/step   Avg/ministep        Total 
 123 steps  236 ministeps     (wasted) 
├────────────────┼───────────┼───────────────┼─────────────┤
 Newton         │   10.0569 │       5.24153 │   1237 (30) │
 Linearization  │   11.9756 │       6.24153 │   1473 (32) │
 Linear solver  │   43.7073 │       22.7797 │  5376 (142) │
 Precond apply  │   87.4146 │       45.5593 │ 10752 (284) │
╰────────────────┴───────────┴───────────────┴─────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2521 │     4.47 % │ 0.3119 │
 Equations     │ 0.2199 │     4.64 % │ 0.3238 │
 Assembly      │ 0.2638 │     5.56 % │ 0.3885 │
 Linear solve  │ 0.3679 │     6.52 % │ 0.4551 │
 Linear setup  │ 2.1869 │    38.73 % │ 2.7052 │
 Precond apply │ 0.1867 │    28.75 % │ 2.0078 │
 Update        │ 0.0843 │     1.49 % │ 0.1043 │
 Convergence   │ 0.3790 │     7.99 % │ 0.5583 │
 Input/Output  │ 0.0325 │     0.11 % │ 0.0077 │
 Other         │ 0.0980 │     1.74 % │ 0.1212 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.6458 │   100.00 % │ 6.9838 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬────────────╮
 Iteration type   Avg/step   Avg/ministep       Total 
 123 steps  208 ministeps    (wasted) 
├────────────────┼───────────┼───────────────┼────────────┤
 Newton         │   7.90244 │       4.67308 │   972 (15) │
 Linearization  │    9.5935 │       5.67308 │  1180 (16) │
 Linear solver  │    37.439 │       22.1394 │ 4605 (112) │
 Precond apply  │    74.878 │       44.2788 │ 9210 (224) │
╰────────────────┴───────────┴───────────────┴────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2501 │     4.47 % │ 0.2431 │
 Equations     │ 0.2142 │     4.65 % │ 0.2528 │
 Assembly      │ 0.2628 │     5.71 % │ 0.3101 │
 Linear solve  │ 0.3882 │     6.94 % │ 0.3773 │
 Linear setup  │ 2.2278 │    39.85 % │ 2.1655 │
 Precond apply │ 0.1872 │    31.74 % │ 1.7245 │
 Update        │ 0.0814 │     1.46 % │ 0.0791 │
 Convergence   │ 0.0725 │     1.57 % │ 0.0856 │
 Input/Output  │ 0.6856 │     2.62 % │ 0.1426 │
 Other         │ 0.0551 │     0.99 % │ 0.0535 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.5906 │   100.00 % │ 5.4341 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────╮
 Iteration type   Avg/step   Avg/ministep     Total 
 123 steps  150 ministeps  (wasted) 
├────────────────┼───────────┼───────────────┼──────────┤
 Newton         │    3.8374 │       3.14667 │  472 (0) │
 Linearization  │   5.05691 │       4.14667 │  622 (0) │
 Linear solver  │   12.9512 │         10.62 │ 1593 (0) │
 Precond apply  │   25.9024 │         21.24 │ 3186 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2470 │     5.07 % │ 0.1166 │
 Equations     │ 0.2061 │     5.57 % │ 0.1282 │
 Assembly      │ 0.2589 │     7.00 % │ 0.1610 │
 Linear solve  │ 0.3063 │     6.29 % │ 0.1446 │
 Linear setup  │ 2.1880 │    44.92 % │ 1.0327 │
 Precond apply │ 0.1888 │    26.16 % │ 0.6015 │
 Update        │ 0.0774 │     1.59 % │ 0.0365 │
 Convergence   │ 0.0683 │     1.85 % │ 0.0425 │
 Input/Output  │ 0.0313 │     0.20 % │ 0.0047 │
 Other         │ 0.0654 │     1.34 % │ 0.0309 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 4.8711 │   100.00 % │ 2.2992 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬────────────╮
 Iteration type   Avg/step   Avg/ministep       Total 
 123 steps  225 ministeps    (wasted) 
├────────────────┼───────────┼───────────────┼────────────┤
 Newton         │   9.47967 │       5.18222 │ 1166 (330) │
 Linearization  │   11.3089 │       6.18222 │ 1391 (352) │
 Linear solver  │   17.6016 │       9.62222 │ 2165 (348) │
 Precond apply  │   35.2033 │       19.2444 │ 4330 (696) │
╰────────────────┴───────────┴───────────────┴────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2467 │     6.24 % │ 0.2877 │
 Equations     │ 0.2758 │     8.33 % │ 0.3837 │
 Assembly      │ 0.2597 │     7.84 % │ 0.3613 │
 Linear solve  │ 0.2277 │     5.76 % │ 0.2655 │
 Linear setup  │ 1.9001 │    48.09 % │ 2.2155 │
 Precond apply │ 0.1925 │    18.10 % │ 0.8337 │
 Update        │ 0.0791 │     2.00 % │ 0.0923 │
 Convergence   │ 0.0716 │     2.16 % │ 0.0996 │
 Input/Output  │ 0.0311 │     0.15 % │ 0.0070 │
 Other         │ 0.0520 │     1.32 % │ 0.0607 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 3.9511 │   100.00 % │ 4.6070 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬─────────────╮
 Iteration type   Avg/step   Avg/ministep        Total 
 123 steps  425 ministeps     (wasted) 
├────────────────┼───────────┼───────────────┼─────────────┤
 Newton         │   25.3008 │       7.32235 │ 3112 (1620) │
 Linearization  │   28.7561 │       8.32235 │ 3537 (1728) │
 Linear solver  │   37.0081 │       10.7106 │ 4552 (1673) │
 Precond apply  │   74.0163 │       21.4212 │ 9104 (3346) │
╰────────────────┴───────────┴───────────────┴─────────────╯
╭───────────────┬────────┬────────────┬─────────╮
 Timing type      Each    Relative    Total 
     ms  Percentage        s 
├───────────────┼────────┼────────────┼─────────┤
 Properties    │ 0.2427 │     6.84 % │  0.7554 │
 Equations     │ 0.2134 │     6.84 % │  0.7547 │
 Assembly      │ 0.2693 │     8.63 % │  0.9524 │
 Linear solve  │ 0.2113 │     5.96 % │  0.6576 │
 Linear setup  │ 1.7730 │    49.98 % │  5.5175 │
 Precond apply │ 0.1893 │    15.61 % │  1.7231 │
 Update        │ 0.0752 │     2.12 % │  0.2341 │
 Convergence   │ 0.0713 │     2.28 % │  0.2522 │
 Input/Output  │ 0.0261 │     0.10 % │  0.0111 │
 Other         │ 0.0581 │     1.64 % │  0.1809 │
├───────────────┼────────┼────────────┼─────────┤
 Total         │ 3.5472 │   100.00 % │ 11.0389 │
╰───────────────┴────────┴────────────┴─────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────╮
 Iteration type   Avg/step   Avg/ministep     Total 
 123 steps  150 ministeps  (wasted) 
├────────────────┼───────────┼───────────────┼──────────┤
 Newton         │   3.95935 │       3.24667 │  487 (0) │
 Linearization  │   5.17886 │       4.24667 │  637 (0) │
 Linear solver  │   18.6911 │       15.3267 │ 2299 (0) │
 Precond apply  │   37.3821 │       30.6533 │ 4598 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2465 │     4.41 % │ 0.1201 │
 Equations     │ 0.2816 │     6.59 % │ 0.1794 │
 Assembly      │ 0.2575 │     6.03 % │ 0.1640 │
 Linear solve  │ 0.3909 │     7.00 % │ 0.1904 │
 Linear setup  │ 2.2309 │    39.94 % │ 1.0864 │
 Precond apply │ 0.1870 │    31.62 % │ 0.8600 │
 Update        │ 0.0775 │     1.39 % │ 0.0377 │
 Convergence   │ 0.0693 │     1.62 % │ 0.0441 │
 Input/Output  │ 0.0341 │     0.19 % │ 0.0051 │
 Other         │ 0.0678 │     1.21 % │ 0.0330 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.5859 │   100.00 % │ 2.7203 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────╮
 Iteration type   Avg/step   Avg/ministep     Total 
 123 steps  148 ministeps  (wasted) 
├────────────────┼───────────┼───────────────┼──────────┤
 Newton         │   3.74797 │       3.11486 │  461 (0) │
 Linearization  │   4.95122 │       4.11486 │  609 (0) │
 Linear solver  │   16.4715 │       13.6892 │ 2026 (0) │
 Precond apply  │   32.9431 │       27.3784 │ 4052 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2462 │     4.60 % │ 0.1135 │
 Equations     │ 0.2087 │     5.15 % │ 0.1271 │
 Assembly      │ 0.2573 │     6.35 % │ 0.1567 │
 Linear solve  │ 0.3696 │     6.90 % │ 0.1704 │
 Linear setup  │ 2.2246 │    41.55 % │ 1.0256 │
 Precond apply │ 0.1873 │    30.74 % │ 0.7587 │
 Update        │ 0.0776 │     1.45 % │ 0.0358 │
 Convergence   │ 0.0689 │     1.70 % │ 0.0420 │
 Input/Output  │ 0.0370 │     0.22 % │ 0.0055 │
 Other         │ 0.0716 │     1.34 % │ 0.0330 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.3540 │   100.00 % │ 2.4682 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────╮
 Iteration type   Avg/step   Avg/ministep     Total 
 123 steps  148 ministeps  (wasted) 
├────────────────┼───────────┼───────────────┼──────────┤
 Newton         │   3.77236 │       3.13514 │  464 (0) │
 Linearization  │   4.97561 │       4.13514 │  612 (0) │
 Linear solver  │   16.3333 │       13.5743 │ 2009 (0) │
 Precond apply  │   32.6667 │       27.1486 │ 4018 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2497 │     4.61 % │ 0.1158 │
 Equations     │ 0.2727 │     6.64 % │ 0.1669 │
 Assembly      │ 0.2590 │     6.31 % │ 0.1585 │
 Linear solve  │ 0.3705 │     6.84 % │ 0.1719 │
 Linear setup  │ 2.2213 │    41.03 % │ 1.0307 │
 Precond apply │ 0.1862 │    29.78 % │ 0.7481 │
 Update        │ 0.0798 │     1.47 % │ 0.0370 │
 Convergence   │ 0.0709 │     1.73 % │ 0.0434 │
 Input/Output  │ 0.0389 │     0.23 % │ 0.0058 │
 Other         │ 0.0732 │     1.35 % │ 0.0340 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.4139 │   100.00 % │ 2.5121 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────╮
 Iteration type   Avg/step   Avg/ministep     Total 
 123 steps  151 ministeps  (wasted) 
├────────────────┼───────────┼───────────────┼──────────┤
 Newton         │   3.78049 │       3.07947 │  465 (0) │
 Linearization  │   5.00813 │       4.07947 │  616 (0) │
 Linear solver  │   17.2033 │       14.0132 │ 2116 (0) │
 Precond apply  │   34.4065 │       28.0265 │ 4232 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2470 │     4.57 % │ 0.1149 │
 Equations     │ 0.2060 │     5.05 % │ 0.1269 │
 Assembly      │ 0.2578 │     6.31 % │ 0.1588 │
 Linear solve  │ 0.3775 │     6.98 % │ 0.1755 │
 Linear setup  │ 2.2168 │    40.99 % │ 1.0308 │
 Precond apply │ 0.1874 │    31.53 % │ 0.7929 │
 Update        │ 0.0769 │     1.42 % │ 0.0358 │
 Convergence   │ 0.0687 │     1.68 % │ 0.0423 │
 Input/Output  │ 0.0336 │     0.20 % │ 0.0051 │
 Other         │ 0.0686 │     1.27 % │ 0.0319 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.4085 │   100.00 % │ 2.5150 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬────────────╮
 Iteration type   Avg/step   Avg/ministep       Total 
 123 steps  208 ministeps    (wasted) 
├────────────────┼───────────┼───────────────┼────────────┤
 Newton         │   7.92683 │        4.6875 │   975 (30) │
 Linearization  │   9.61789 │        5.6875 │  1183 (32) │
 Linear solver  │   33.6667 │       19.9087 │ 4141 (128) │
 Precond apply  │   67.3333 │       39.8173 │ 8282 (256) │
╰────────────────┴───────────┴───────────────┴────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2464 │     4.72 % │ 0.2403 │
 Equations     │ 0.2089 │     4.85 % │ 0.2471 │
 Assembly      │ 0.2566 │     5.96 % │ 0.3035 │
 Linear solve  │ 0.3555 │     6.81 % │ 0.3466 │
 Linear setup  │ 2.2162 │    42.43 % │ 2.1608 │
 Precond apply │ 0.1860 │    30.25 % │ 1.5402 │
 Update        │ 0.0766 │     1.47 % │ 0.0747 │
 Convergence   │ 0.0700 │     1.63 % │ 0.0828 │
 Input/Output  │ 0.2159 │     0.88 % │ 0.0449 │
 Other         │ 0.0523 │     1.00 % │ 0.0510 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.2226 │   100.00 % │ 5.0920 │
╰───────────────┴────────┴────────────┴────────╯

Plot the oil rate at the producer over the ensemble

julia
using Statistics
fig = Figure()
ax = Axis(fig[1, 1])
for i = 1:N
    ws = wells[i]
    q = -ws[:Producer][:orat]
    lines!(ax, report_time, q)
end
xlims!(ax, [mean(report_time), report_time[end]])
ylims!(ax, 0, 0.0075)
fig

Plot the average saturation over the ensemble

julia
avg = mean(saturations)
fig = Figure()
h = nothing
ax = Axis(fig[1, 1])
h = contourf!(ax, avg)
fig

Plot the isocontour lines over the ensemble

julia
fig = Figure()
h = nothing
ax = Axis(fig[1, 1])
for s in saturations
    contour!(ax, s, levels = 0:0.1:1)
end
fig

Example on GitHub

If you would like to run this example yourself, it can be downloaded from the JutulDarcy.jl GitHub repository as a script, or as a Jupyter Notebook

This example took 54.341694917 seconds to complete.

This page was generated using Literate.jl.