Running a simulation

This section describes some of the details of the various methods and options for setting-up QDYN simulations. For practical examples, see the Tutorials section.

Using the wrapper

The core of QDYN is a Fortran code. While the format of its input and output files is well defined, we find it more convenient to set up the input parameters, perform simulations, and read the output data within a MATLAB (qdyn.m) or Python (pyqdyn.py) environment. The operation of both these wrappers is detailed below.

MATLAB wrapper

You first need to set in MATLAB the full path to the src directory, for instance:

addpath ~/qdyn-read-only/src

Tip for Mac users: if you get an error message related to gfortran libraries (e.g. libgfortran.3.dylib) when running qdyn in MATLAB, execute:

setenv('DYLD_LIBRARY_PATH', '/usr/local/bin')

The second argument should be the path to your gfortran libraries (sometimes set to /opt/local/lib/libgcc).

The general usage syntax is:

[pars, ot, ox] = qdyn(mode, [parsin], ['Property', Value, ...])

The default input values can be listed by executing:

pars = qdyn('set')

Permitted arguments for mode are:

The [parsin] list is a structure to override the default parameters, using Property-Value pairs. See the section on input parameters below.

Python wrapper

As an alternative to the MATLAB wrapper, users can chose to generate their input files through a Python wrapper, called pyqdyn.py. While the variable naming conventions are the same as in the MATLAB wrapper, the data structures and operation modes are somewhat different. The Python wrapper relies heavily on dictionary objects for both the simulation parameters and the generated mesh. A basic operation workflow looks as follows:

1. After importing the pyqdyn.py module, instantiate QDYN as: p = qdyn()

2. Create a dictionary with simulation parameters to override the default values, for example:

   # General simulation parameters
   set_dict = {
       "FAULT_TYPE": 1,
       "ACC": 1e-10,
       "SOLVER": 2,
       "TMAX": 20,
       "MESHDIM": 1,
       "SIGMA": 5e6,
       "L": 100e3,
       "W": 100e5,
       "FRICTION_MODEL": "RSF",
   }
      
   # RSF parameters
   set_dict_RSF = {
       "RNS_LAW": 0,
       "THETA_LAW": 1,
       "A": 0.001,
       "B": 0.0015,
       "DC": 1e-5,
       "MU_SS": 0.6,
       "V_SS": 1e-6,
   }
      
   # Merge dicts
   set_dict["SET_DICT_RSF"] = set_dict_RSF
   

3. Pass the parameter values to the QDYN instance and generate the mesh:

   p.settings(set_dict)
   p.render_mesh()
   

4. The mesh is now populated the values defined in set_dict (along with the default settings, if left unchanged). Individual mesh node parameters can be accessed through p.mesh_dict["PARAM"][element], where "PARAM" is a given variable name (see input parameters), and element is the index of the node (note that Python index count starts at 0, and not at 1 like in MATLAB). The mesh is stored as a set of NumPy arrays, so slicing and vector operations are permitted (for instance: p.mesh_dict["SIGMA"][50:100] = 10e6).

5. Write the input file and run the simulation as:

   p.write_input()
   p.run()
   

6. During and after the simulation, the current/final simulation output can be read using p.read_output(). This creates two new objects contained in the QDYN instance: p.ot contains the time series data, and p.ox contains the fault snapshots. See the section “Simulation output structure” for a list of output parameters. An example of how to access output data:

   time = p.ot_vmax["t"]  # Vector of simulation time steps
   Vmax = p.ot_vmax["v"]  # Vector of maximum slip velocity
      
   # Import matplotlib and plot results
   import matplotlib.pyplot as plt
   plt.plot(time, Vmax)
   plt.show()
   

Simulation input parameters

Parameters defining the geometry of the problem and loading:

Simulation features (only available via the Python wrapper):

QDYN offers various optional simulation features. Set the following parameters to 1 to enable, or 0 to disable (default).

Rate-and-state friction (RSF) parameters:

Note that with the Python wrapper, the rheological parameters below are set in a dictionary SET_DICT_RSF nested within set_dict.

Chen-Niemeijer-Spiers (CNS) model parameters (only available via the Python wrapper):

Note that with the Python wrapper, the rheological parameters below are set in a dictionary SET_DICT_CNS nested within set_dict.

Thermal pressurisation model parameters (only available via Python wrapper):

Note that with the Python wrapper, the rheological parameters below are set in a dictionary SET_DICT_TP nested within set_dict.

Initial conditions:

Discretization and accuracy parameters:

Output control parameters:

Parameters for integration with dynamic code:

NOTE: the QDYN-SPECFEM3D bridge (QSB) is currently incompatible with the input/output structures of QDYN. Proper functioning of this feature is not guaranteed, and use of the QSB is not advised.

Simulation output structure

Time-series output (output_ot_x)

Time-series output at Vmax (output_vmax)

Snapshot output (output_ox, output_ox_dyn)

Examples

There are a number of (self-documented) examples provided in the examples/ directory. Interactive Jupyter Notebooks are found in examples/notebooks/ (see also the Tutorials section).