How PetroScript works

Mental model

word problem
  -> the facts the reservoir states: grid, rock, fluid, initial state, wells, schedule
  -> ps.Model(...)              or ps.load("CASE.DATA") for a deck you already have
  -> m.check()                  PetroScript's own checks: no simulator, never raises
  -> ps.compile(m)              a Deck in memory: deck.files, deck.check()
  -> deck.write(folder)         or deck.run(folder=...) to run OPM Flow (or JutulDarcy) on it
  -> ps.results.open(folder)    the run, read by PetroScript's names

PetroScript does not ask you to memorize deck keywords. You state the reservoir by keyword argument, in the model’s own unit system, and the compile writes the cards.

The engineer’s pipeline

Stage	Where it is stated	Typical question
dimensions	ps.Model(name=, units=, start=); capacity=ps.Capacity(...) only to reserve more than the model holds	What is the run, in which units, from which date?
grid	grid=ps.grid.Cartesian(...), CornerPoint, Radial, Spider	What cells exist and where are they?
rock	rock=ps.Rock(...); region maps in regions=ps.Regions(...)	What porosity and permeability does Flow see?
physics	fluid=ps.fluid.*, its tables in ps.pvt; relperm=[ps.relperm.*]; family packages in physics=[...]	Which phases flow, with which PVT and saturation functions?
initial	initial=ps.initial.Equilibrium(...), ExplicitState, Restart	What pressure and contacts start the run?
summary	summary=ps.Summary(...); what else a run writes in output=ps.Output(...)	What should the simulator report?
schedule	wells=[ps.Well(...)], groups=[ps.Group(...)], schedule=ps.Schedule(events=[ps.Event(...)])	What changes over time?

The model counts its own sizes (wells, connections, table rows), so a model states no dimensioning cards. A loaded deck keeps the sizes it states, in m.capacity.

From scratch

A water injector supporting an oil producer, in field units:

import petroscript as ps

flood = ps.Model(
    name="FLOOD",
    units="field",                     # a bare number is in field units: psia, ft, stb/day
    start="2020-01-01",
    grid=ps.grid.Cartesian(shape=(5, 5, 1), dx=300.0, dy=300.0, dz=30.0, top=8000.0),
    rock=ps.Rock(porosity=0.25, permeability=200.0, vertical_ratio=0.1,   # kz = 0.1 kx
                 compressibility=4.0e-6, reference_pressure=3500.0),
    fluid=ps.fluid.OilWater(
        oil=ps.pvt.DeadOil(pressure=[1000.0, 3000.0, 5000.0],
                           oil_fvf=[1.25, 1.22, 1.20],
                           oil_viscosity=[1.0, 1.1, 1.2]),
        water=ps.pvt.Water(reference_pressure=3500.0, fvf=1.01, compressibility=3.0e-6,
                           viscosity=0.5, viscosibility=0.0),
        oil_surface_density=53.0,
        water_surface_density=64.0,
    ),
    relperm=[ps.relperm.OilWater(sw=[0.2, 0.5, 0.8], krw=[0.0, 0.15, 0.5],
                                 krow=[1.0, 0.3, 0.0], pcow=[0.0, 0.0, 0.0])],
    initial=ps.initial.Equilibrium(datum_depth=8015.0, datum_pressure=3500.0,
                                   oil_water_contact=8100.0),
    wells=[ps.Well(name="P1", i=5, j=5), ps.Well(name="I1", i=1, j=1)],
    schedule=ps.Schedule(
        end="2021-01-01",
        report="monthly",
        events=[ps.Event("2020-01-01",
                         ps.Perforate(well="P1", k=(1, 1), diameter=0.5),
                         ps.Perforate(well="I1", k=(1, 1), diameter=0.5),
                         ps.Produce(well="P1", oil_rate=500.0, bhp_limit=1500.0),
                         ps.Inject(well="I1", phase="water", rate=600.0, bhp_limit=5000.0))],
    ),
    summary=ps.Summary(field=["oil_production_rate", "water_cut", "pressure"],
                       wells={"*": ["bottom_hole_pressure"]}),
)

A model that contradicts itself is refused as it is built, with a ps.ModelError whose rows name the code, the address and the fix. Leave out vertical_ratio= above and the model refuses: [rock.permeability_z.unstated], with a fix line naming vertical_ratio= and permeability_z=. A model that builds is then checked, compiled and written:

assert [row for row in flood.check() if row.severity == "error"] == []
deck = ps.compile(flood)              # in memory: nothing is written or run
text = deck.files["FLOOD.DATA"]
assert "WCONINJE" in text and "\nOIL\nWATER\n" in text
deck.write("flood")                   # flood/FLOOD.DATA

deck.run(folder=...) runs OPM Flow on the files and returns the run; QUICKSTART section 4 runs one and reads it.

A second simulator: JutulDarcy

ps.compile(m, "jutuldarcy") writes the same model for JutulDarcy 0.4.0, an open-source simulator written in Julia. Its capability table checks the model first and refuses, at the address that states it, whatever JutulDarcy does not honour exactly: DRSDT, a THP limit, a group control, a well’s efficiency, and every keyword never measured on it. A group or region summary quantity is a warning: JutulDarcy does not report it, and the run says results.not_read for it. The compile for OPM Flow is untouched.

jutul = ps.compile(flood, "jutuldarcy")    # the same model, as JutulDarcy reads it
assert (jutul.target, list(jutul.files)) == ("jutuldarcy", ["FLOOD.DATA"])

PETROSCRIPT_JUTULDARCY names the runner: a program that runs Julia with JutulDarcy and GeoEnergyIO loadable, or a relay into a container that holds them. With it set, jutul.check() asks JutulDarcy to parse and set up the deck, and jutul.run(folder=...) runs it and returns a run read by PetroScript’s names, as a Flow run is: run.field("oil_production_rate"), run.wells("bottom_hole_pressure"), in the deck’s units with Flow’s unit words, one row per report. The reader serves the field’s and the wells’ rates, totals, ratios, pressure and volumes in place from JutulDarcy’s own summary, with PetroScript’s definitions: field pressure weighted by the current hydrocarbon pore volume, and each total integrated from its rate over every solver step (JutulDarcy 0.4.0’s own totals are 86,400 times too large). Cell arrays, the grid and the log are not read yet, and each says so. examples/two_targets.py compiles SPE1 for both simulators, runs both and compares the runs by name.

Phase sets

The fluid’s class is its formulation, never a card you type. ps.fluid.OilWater flows oil and water and ps.fluid.GasWater gas and water; ps.fluid.BlackOil is every other phase set (three phases, oil and gas, water alone, or CO2 or H2 storage with a ps.storage.Package in physics=[...]). A live oil (ps.pvt.LiveOil) states the gas dissolved in it and a wet gas (ps.pvt.WetGas) the oil vaporized in it. The compile writes the phase cards from the fluid:

assert "\nGAS\n" not in text           # OilWater: no gas phase

Edit with verbs

A model is one frozen value. Every verb returns a new model:

later = flood.until("2022-01-01", report="monthly")            # move the run's end
later = later.at("2021-06-01", ps.Adjust(well="P1", oil_rate=400.0))
later = later.at("2021-06-01", ps.Adjust(well="I1", rate=900.0))      # an injector's target
later = later.set("rock/porosity", 0.22)                       # any fact by its address
assert [change.address for change in flood.diff(later)][:1] == ["rock/porosity"]

Each statement carries its own replace, typed, so a misspelt field is refused:

faster = later.control("I1", at="2021-06-01")
assert isinstance(faster, ps.Inject)
faster = faster.replace(rate=950.0)            # not dataclasses.replace: no statement is a dataclass

Every summary quantity has a name and an OPM Flow vector, and ps.help answers either way: ps.help("FGPT") names gas_production_total, and ps.help("gas_production_total") names FGPT, WGPT and the rest. ps.help("WBHP") prints:

WBHP: the well quantity bottom_hole_pressure
  ps.Summary(wells={"*": ["bottom_hole_pressure"]}); run.wells("bottom_hole_pressure")

QUICKSTART section 5 has the edits that come up again and again: one quantity from a date on, a well’s later plan restated, a new well mid-run.

Look it up

ps.help(ps.Inject)        # a class: its fields, the card it writes, where the model holds it
ps.help("WCONINJE")       # a deck keyword: the class that writes it
ps.help("units")          # a topic; ps.help() lists every class, namespace and topic

On nearly every class a misspelt field is refused with the nearest spelling (did you mean 'bhp_limit'?), so read the refusal before guessing again.