Setting up simulations using scripts

Because many scripting languages are able to parse JSON files (see JSON files in scripts), it is possible to write scripts to set up Waiwera JSON input files. This can simplify some aspects of simulation setup (for example, setting up boundary conditions, or large numbers of sources).

A script for setting up a simulation must usually use some form of library for interacting with the simulation mesh, for example to identify the correct cells for particular boundary conditions or sources. The appropriate mesh library depends on the type of mesh being used (see Mesh formats).

One option for creating meshes for use with Waiwera (see Creating meshes) is Layermesh, a Python library for creating and manipulating meshes with a layer/column structure. By importing the Layermesh module together with Python’s built-in JSON module, both the mesh and all other aspects of a Waiwera simulation can be set up from a Python script. (This is similar to using the PyTOUGH library to set up a TOUGH2 simulation.)

Example

The example Python script below sets up a simple steady-state non-isothermal 3-D geothermal simulation and saves it to an input JSON file. The simulation uses the default Water and energy (“we”) equation of state.

A simulation mesh is created with 10 × 1 km cells in the x- and y-directions, and 25 layers of varying thicknesses. Scattered surface elevation data are fitted to the mesh to represent a spatially-varying water table. The mesh is written to a Layermesh HDF5 file (which can be used for post-processing) and an ExodusII mesh file for use with Waiwera.

Uniform initial conditions are prescribed for the start of the steady-state run, and atmospheric boundary conditions are applied at the top surface of the model. Two rock types are created, one for the reservoir and another for cap-rock near the top of the model. An upflow of hot water (enthalpy 1300 kJ/kg) is introduced around the centre of the model.

Time-stepping parameters are specified for the simulation, which stops when the time reaches 10¹⁶ seconds, approximating a steady state. Simulation output for the initial state and intermediate time steps is disabled, so that the solution is only output at the end of the steady-state run.

Finally, the simulation is output to a JSON file (with a specified indent size and sorted keys).

After the simulation is run, the Layermesh library can also be used to plot the results (see Simulation output and scripts).

# Waiwera demo problem

import json
import layermesh.mesh as lm

# set up mesh:
dx = [1000] * 10
dy = dx
dz = [10] * 5 + [20] * 5 + [50] * 5 + [100] * 5 + [200] * 5
mesh = lm.mesh(rectangular = (dx, dy, dz))
mesh.translate((0,0,200))
surface_data = [(0,0,0), (10e3,0,0), (10e3,10e3,0), (0,10e3,0,0),
                (5e3,5e3,200)]
mesh.fit_surface(surface_data)
mesh.write('demo_mesh.h5')
mesh_filename = 'demo_mesh.exo'
mesh.export(mesh_filename)

sim = {}
sim['title'] = 'Waiwera demo problem'
sim['mesh'] = {'filename': mesh_filename}

# initial and boundary conditions:
P0, T0 = 1.e5, 20.
sim['initial'] = {'primary': [P0, T0]}
sim['boundaries'] = [{'primary': [P0, T0], 'region': 1,
                    'faces': {'cells': [c.index for c in mesh.surface_cells],
                              'normal': [0, 0, 1]}}]

# rock properties:
sim['mesh']['zones'] = {'cap': {'z': [0, 200]}, 'reservoir': {'-': 'cap'}}
sim['rock'] = {'types': [
    {'name': 'rock', 'permeability': 20.e-15, 'zones': 'reservoir'},
    {'name': 'caprock', 'permeability': 0.5e-15, 'zones': 'cap'}
]}

# upflow:
upflow_cols = mesh.find([(4000,4000), (6000,6000)])
upflow_cells = [col.cell[-1].index for col in upflow_cols]
sim['source'] = [{'rate': 100, 'enthalpy': 1300e3,
                  'cells': upflow_cells}]

# time stepping:
sim['time'] = {
    'step': {
        'size': 1e6,
        'adapt': {'on': True},
        'maximum': {'number': 100}
    },
    'stop': 1e16,
}

sim['output'] = {'initial': False, 'frequency': 0}

json.dump(sim, open('demo.json', 'w'), indent = 2, sort_keys = True)