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.252 │       8.94326 │ 1261 (0) │
 Precond apply  │   20.5041 │       17.8865 │ 2522 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2543 │     5.45 % │ 0.1040 │
 Equations     │ 0.2267 │     6.54 % │ 0.1247 │
 Assembly      │ 0.2861 │     8.25 % │ 0.1574 │
 Linear solve  │ 0.3102 │     6.65 % │ 0.1269 │
 Linear setup  │ 1.8837 │    40.40 % │ 0.7704 │
 Precond apply │ 0.1988 │    26.28 % │ 0.5013 │
 Update        │ 0.0853 │     1.83 % │ 0.0349 │
 Convergence   │ 0.0788 │     2.27 % │ 0.0434 │
 Input/Output  │ 0.0479 │     0.35 % │ 0.0067 │
 Other         │ 0.0919 │     1.97 % │ 0.0376 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 4.6631 │   100.00 % │ 1.9072 │
╰───────────────┴────────┴────────────┴────────╯

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  224 ministeps    (wasted) 
├────────────────┼───────────┼───────────────┼────────────┤
 Newton         │   9.54472 │       5.24107 │ 1174 (375) │
 Linearization  │   11.3659 │       6.24107 │ 1398 (400) │
 Linear solver  │   17.9268 │       9.84375 │ 2205 (401) │
 Precond apply  │   35.8537 │       19.6875 │ 4410 (802) │
╰────────────────┴───────────┴───────────────┴────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2557 │     5.54 % │ 0.3001 │
 Equations     │ 0.2510 │     6.48 % │ 0.3509 │
 Assembly      │ 0.2948 │     7.61 % │ 0.4121 │
 Linear solve  │ 0.2516 │     5.45 % │ 0.2953 │
 Linear setup  │ 2.0219 │    43.83 % │ 2.3737 │
 Precond apply │ 0.1923 │    15.66 % │ 0.8481 │
 Update        │ 0.0906 │     1.96 % │ 0.1064 │
 Convergence   │ 0.4194 │    10.83 % │ 0.5864 │
 Input/Output  │ 0.0406 │     0.17 % │ 0.0091 │
 Other         │ 0.1138 │     2.47 % │ 0.1336 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 4.6130 │   100.00 % │ 5.4157 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬─────────────╮
 Iteration type   Avg/step   Avg/ministep        Total 
 123 steps  382 ministeps     (wasted) 
├────────────────┼───────────┼───────────────┼─────────────┤
 Newton         │   21.5854 │       6.95026 │ 2655 (1200) │
 Linearization  │   24.6911 │       7.95026 │ 3037 (1280) │
 Linear solver  │   31.4472 │       10.1257 │ 3868 (1254) │
 Precond apply  │   62.8943 │       20.2513 │ 7736 (2508) │
╰────────────────┴───────────┴───────────────┴─────────────╯
╭───────────────┬────────┬────────────┬─────────╮
 Timing type      Each    Relative    Total 
     ms  Percentage        s 
├───────────────┼────────┼────────────┼─────────┤
 Properties    │ 0.2502 │     6.54 % │  0.6643 │
 Equations     │ 0.2870 │     8.58 % │  0.8716 │
 Assembly      │ 0.2871 │     8.59 % │  0.8720 │
 Linear solve  │ 0.2317 │     6.06 % │  0.6153 │
 Linear setup  │ 1.8832 │    49.24 % │  4.9999 │
 Precond apply │ 0.1955 │    14.90 % │  1.5128 │
 Update        │ 0.0833 │     2.18 % │  0.2212 │
 Convergence   │ 0.0815 │     2.44 % │  0.2476 │
 Input/Output  │ 0.0317 │     0.12 % │  0.0121 │
 Other         │ 0.0516 │     1.35 % │  0.1370 │
├───────────────┼────────┼────────────┼─────────┤
 Total         │ 3.8244 │   100.00 % │ 10.1537 │
╰───────────────┴────────┴────────────┴─────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬───────────╮
 Iteration type   Avg/step   Avg/ministep      Total 
 123 steps  164 ministeps   (wasted) 
├────────────────┼───────────┼───────────────┼───────────┤
 Newton         │   4.88618 │       3.66463 │  601 (30) │
 Linearization  │   6.21951 │       4.66463 │  765 (32) │
 Linear solver  │    12.065 │       9.04878 │ 1484 (30) │
 Precond apply  │   24.1301 │       18.0976 │ 2968 (60) │
╰────────────────┴───────────┴───────────────┴───────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2567 │     5.40 % │ 0.1543 │
 Equations     │ 0.2472 │     6.62 % │ 0.1891 │
 Assembly      │ 0.2887 │     7.73 % │ 0.2209 │
 Linear solve  │ 0.3584 │     7.54 % │ 0.2154 │
 Linear setup  │ 2.2089 │    46.45 % │ 1.3275 │
 Precond apply │ 0.1946 │    20.21 % │ 0.5777 │
 Update        │ 0.0890 │     1.87 % │ 0.0535 │
 Convergence   │ 0.0845 │     2.26 % │ 0.0647 │
 Input/Output  │ 0.0441 │     0.25 % │ 0.0072 │
 Other         │ 0.0800 │     1.68 % │ 0.0481 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 4.7559 │   100.00 % │ 2.8583 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────────╮
 Iteration type   Avg/step   Avg/ministep         Total 
 123 steps  536 ministeps      (wasted) 
├────────────────┼───────────┼───────────────┼──────────────┤
 Newton         │   34.3008 │       7.87127 │  4219 (2565) │
 Linearization  │   38.6585 │       8.87127 │  4755 (2736) │
 Linear solver  │   44.6585 │       10.2481 │  5493 (2621) │
 Precond apply  │   89.3171 │       20.4963 │ 10986 (5242) │
╰────────────────┴───────────┴───────────────┴──────────────╯
╭───────────────┬────────┬────────────┬─────────╮
 Timing type      Each    Relative    Total 
     ms  Percentage        s 
├───────────────┼────────┼────────────┼─────────┤
 Properties    │ 0.2489 │     6.92 % │  1.0503 │
 Equations     │ 0.2413 │     7.56 % │  1.1472 │
 Assembly      │ 0.2939 │     9.20 % │  1.3974 │
 Linear solve  │ 0.2552 │     7.09 % │  1.0768 │
 Linear setup  │ 1.7696 │    49.17 % │  7.4661 │
 Precond apply │ 0.1915 │    13.86 % │  2.1042 │
 Update        │ 0.0815 │     2.26 % │  0.3439 │
 Convergence   │ 0.0798 │     2.50 % │  0.3793 │
 Input/Output  │ 0.0289 │     0.10 % │  0.0155 │
 Other         │ 0.0482 │     1.34 % │  0.2034 │
├───────────────┼────────┼────────────┼─────────┤
 Total         │ 3.5990 │   100.00 % │ 15.1840 │
╰───────────────┴────────┴────────────┴─────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬───────────╮
 Iteration type   Avg/step   Avg/ministep      Total 
 123 steps  157 ministeps   (wasted) 
├────────────────┼───────────┼───────────────┼───────────┤
 Newton         │   4.37398 │       3.42675 │  538 (15) │
 Linearization  │   5.65041 │       4.42675 │  695 (16) │
 Linear solver  │    13.935 │       10.9172 │ 1714 (45) │
 Precond apply  │   27.8699 │       21.8344 │ 3428 (90) │
╰────────────────┴───────────┴───────────────┴───────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2581 │     5.02 % │ 0.1389 │
 Equations     │ 0.2529 │     6.36 % │ 0.1758 │
 Assembly      │ 0.2921 │     7.34 % │ 0.2030 │
 Linear solve  │ 0.3315 │     6.45 % │ 0.1783 │
 Linear setup  │ 2.3025 │    44.79 % │ 1.2387 │
 Precond apply │ 0.1951 │    24.19 % │ 0.6689 │
 Update        │ 0.0903 │     1.76 % │ 0.0486 │
 Convergence   │ 0.0860 │     2.16 % │ 0.0598 │
 Input/Output  │ 0.0463 │     0.26 % │ 0.0073 │
 Other         │ 0.0864 │     1.68 % │ 0.0465 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.1409 │   100.00 % │ 2.7658 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬─────────────╮
 Iteration type   Avg/step   Avg/ministep        Total 
 123 steps  194 ministeps     (wasted) 
├────────────────┼───────────┼───────────────┼─────────────┤
 Newton         │   6.82114 │       4.32474 │   839 (165) │
 Linearization  │   8.39837 │       5.32474 │  1033 (176) │
 Linear solver  │   31.2846 │       19.8351 │  3848 (554) │
 Precond apply  │   62.5691 │       39.6701 │ 7696 (1108) │
╰────────────────┴───────────┴───────────────┴─────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2579 │     4.52 % │ 0.2164 │
 Equations     │ 0.2946 │     6.36 % │ 0.3044 │
 Assembly      │ 0.2913 │     6.28 % │ 0.3009 │
 Linear solve  │ 0.4163 │     7.29 % │ 0.3493 │
 Linear setup  │ 2.3098 │    40.47 % │ 1.9379 │
 Precond apply │ 0.1879 │    30.20 % │ 1.4459 │
 Update        │ 0.0900 │     1.58 % │ 0.0755 │
 Convergence   │ 0.0862 │     1.86 % │ 0.0890 │
 Input/Output  │ 0.0413 │     0.17 % │ 0.0080 │
 Other         │ 0.0724 │     1.27 % │ 0.0608 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.7069 │   100.00 % │ 4.7881 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬─────────────╮
 Iteration type   Avg/step   Avg/ministep        Total 
 123 steps  345 ministeps     (wasted) 
├────────────────┼───────────┼───────────────┼─────────────┤
 Newton         │   19.0976 │        6.8087 │ 2349 (1245) │
 Linearization  │   21.9024 │        7.8087 │ 2694 (1328) │
 Linear solver  │   28.8862 │       10.2986 │ 3553 (1291) │
 Precond apply  │   57.7724 │       20.5971 │ 7106 (2582) │
╰────────────────┴───────────┴───────────────┴─────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2661 │     6.99 % │ 0.6252 │
 Equations     │ 0.2501 │     7.54 % │ 0.6739 │
 Assembly      │ 0.3034 │     9.14 % │ 0.8173 │
 Linear solve  │ 0.2471 │     6.49 % │ 0.5805 │
 Linear setup  │ 1.8438 │    48.44 % │ 4.3311 │
 Precond apply │ 0.1918 │    15.25 % │ 1.3632 │
 Update        │ 0.0831 │     2.18 % │ 0.1952 │
 Convergence   │ 0.0811 │     2.44 % │ 0.2186 │
 Input/Output  │ 0.0330 │     0.13 % │ 0.0114 │
 Other         │ 0.0529 │     1.39 % │ 0.1243 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 3.8061 │   100.00 % │ 8.9406 │
╰───────────────┴────────┴────────────┴────────╯
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.8374 │       3.12583 │  472 (0) │
 Linearization  │   5.06504 │       4.12583 │  623 (0) │
 Linear solver  │    15.122 │       12.3179 │ 1860 (0) │
 Precond apply  │   30.2439 │       24.6358 │ 3720 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2575 │     4.74 % │ 0.1215 │
 Equations     │ 0.2498 │     6.07 % │ 0.1556 │
 Assembly      │ 0.2911 │     7.07 % │ 0.1814 │
 Linear solve  │ 0.3732 │     6.87 % │ 0.1762 │
 Linear setup  │ 2.2941 │    42.23 % │ 1.0828 │
 Precond apply │ 0.1891 │    27.44 % │ 0.7036 │
 Update        │ 0.0889 │     1.64 % │ 0.0420 │
 Convergence   │ 0.0834 │     2.03 % │ 0.0519 │
 Input/Output  │ 0.0462 │     0.27 % │ 0.0070 │
 Other         │ 0.0895 │     1.65 % │ 0.0423 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.4326 │   100.00 % │ 2.5642 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬────────────╮
 Iteration type   Avg/step   Avg/ministep       Total 
 123 steps  180 ministeps    (wasted) 
├────────────────┼───────────┼───────────────┼────────────┤
 Newton         │   5.88618 │       4.02222 │   724 (75) │
 Linearization  │   7.34959 │       5.02222 │   904 (80) │
 Linear solver  │   14.5366 │       9.93333 │ 1788 (111) │
 Precond apply  │   29.0732 │       19.8667 │ 3576 (222) │
╰────────────────┴───────────┴───────────────┴────────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2569 │     5.43 % │ 0.1860 │
 Equations     │ 0.2495 │     6.59 % │ 0.2255 │
 Assembly      │ 0.2917 │     7.70 % │ 0.2637 │
 Linear solve  │ 0.3476 │     7.35 % │ 0.2516 │
 Linear setup  │ 2.2059 │    46.64 % │ 1.5971 │
 Precond apply │ 0.1957 │    20.44 % │ 0.6999 │
 Update        │ 0.0882 │     1.86 % │ 0.0638 │
 Convergence   │ 0.0837 │     2.21 % │ 0.0756 │
 Input/Output  │ 0.0430 │     0.23 % │ 0.0077 │
 Other         │ 0.0739 │     1.56 % │ 0.0535 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 4.7301 │   100.00 % │ 3.4246 │
╰───────────────┴────────┴────────────┴────────╯
Jutul: Simulating 9 years, 44.69 weeks as 123 report steps
╭────────────────┬───────────┬───────────────┬──────────╮
 Iteration type   Avg/step   Avg/ministep     Total 
 123 steps  149 ministeps  (wasted) 
├────────────────┼───────────┼───────────────┼──────────┤
 Newton         │   3.89431 │       3.21477 │  479 (0) │
 Linearization  │   5.10569 │       4.21477 │  628 (0) │
 Linear solver  │   11.0976 │       9.16107 │ 1365 (0) │
 Precond apply  │   22.1951 │       18.3221 │ 2730 (0) │
╰────────────────┴───────────┴───────────────┴──────────╯
╭───────────────┬────────┬────────────┬────────╮
 Timing type      Each    Relative   Total 
     ms  Percentage       s 
├───────────────┼────────┼────────────┼────────┤
 Properties    │ 0.2580 │     5.10 % │ 0.1236 │
 Equations     │ 0.2493 │     6.47 % │ 0.1566 │
 Assembly      │ 0.2903 │     7.53 % │ 0.1823 │
 Linear solve  │ 0.3074 │     6.08 % │ 0.1472 │
 Linear setup  │ 2.3740 │    46.95 % │ 1.1371 │
 Precond apply │ 0.1946 │    21.94 % │ 0.5313 │
 Update        │ 0.0884 │     1.75 % │ 0.0423 │
 Convergence   │ 0.0833 │     2.16 % │ 0.0523 │
 Input/Output  │ 0.0472 │     0.29 % │ 0.0070 │
 Other         │ 0.0875 │     1.73 % │ 0.0419 │
├───────────────┼────────┼────────────┼────────┤
 Total         │ 5.0561 │   100.00 % │ 2.4219 │
╰───────────────┴────────┴────────────┴────────╯

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 69.731553055 seconds to complete.

This page was generated using Literate.jl.