Running#

Restarting a crashed simulation#

If the reload_results input argument is set to true, an interrupted response-functions mode simulation will not recompute already computed response function simulation and resume where it was interrupted when restarted.

The observation operator plugin autorestart input argument can also be used to resume the simulation that was interrupted.

Exploit results#

The following code examples can be used to read and exploit the response functions results:

\(\mathbf{H}\) matrix#

Reading the \(\mathbf{H}\) matrix (with the control vector coordinates) from the observation vector decomposition.

Note

The response function simulation need to be run with the input argument dump_obsvect_decompostion set to true to the CIF produces the necessary outputs.

import os

import numpy as np
import pandas as pd
import xarray as xr
from pycif.utils.sparse_array import to_dense_dataset

def read_decomp(path, region_names=None, use_sparse_dtype=False):
    with xr.open_dataset(path) as ds:
        if 'sparse_sim_tl' in ds:
            ds = to_dense_dataset(ds)

        ds = ds.isel(vresol=0).drop_vars('vresol')  # Dropping vertical levels

        # Replacing hresol with region names
        if region_names is not None and 'region_id' in ds:
            ds['hresol'] = (['hresol'], region_names[ds.region_id.values])

        decomp = ds.sim_tl.to_series().unstack(['date', 'hresol'])

    if use_sparse_dtype:
        decomp = decomp.astype(pd.SparseDtype("float", 0.0))

    return decomp

def read_h_matrix(workdir, region_names=None):
    h_dict = {}
    dirname = os.path.join(workdir, "base_functions/decomposition")

    for component in os.listdir(dirname):
        comp_dirname = os.path.join(dirname, component)
          for parameter in os.listdir(comp_dirname):
              print(f"reading H matrix for {component}/{parameter}")
              path = os.path.join(comp_dirname, parameter, "obsvect_decomp.nc")
              h_dict[(component, parameter)] = read_decomp(path, region_names, True)

    df = pd.concat(
        h_dict.values(),
        keys=h_dict.keys(),
        names=['component', 'parameter'],
        axis='columns',
    )

    return df


with xr.open_dataset("region_mask.nc") as ds_mask:
    region_names = ds_mask.regions_name.values

h_matrix = read_h_matrix("path/to/workdir", region_names)

Variational inv#

  • TODO: Output

  • TODO: Rerun with another B, R