Store netCDF output

Store netCDF output#

import matplotlib.pyplot as plt
import netCDF4
from typhon import plots

import konrad


plots.styles.use("typhon")

During the configuration of an RCE simulation we can specify an output netCDF file. Konrad will create the file and store a snapshot of the instantaneous model state at a given output frequency.

plev, phlev = konrad.utils.get_pressure_grids(1000e2, 1, 128)
atmosphere = konrad.atmosphere.Atmosphere(phlev)

rce = konrad.RCE(
    atmosphere,
    surface=konrad.surface.FixedTemperature(temperature=288.),  # Run with a fixed surface temperature.
    timestep='12h',  # Set timestep in model time.
    max_duration='100d',  # Set maximum runtime.
    outfile="my_rce_output.nc",  # Specifiy the output file
    writeevery="5d",  # Set the output frequency
)
rce.run()

Read and visualize model output#

After the run has finished, we can use the netCDF4 package to read the data again. Here, we plot the temporal evolution of the temperature profile from the inital US standard atmosphere towards the equilibrium RCE state.

with netCDF4.Dataset("my_rce_output.nc", "r") as root:
    p = root["plev"][:]
    T = root["atmosphere/T"][:]

sm = plt.cm.ScalarMappable(norm=plt.Normalize(0, T.shape[0]), cmap="magma")

fig, ax = plt.subplots()
for i, temp in enumerate(T):
    plots.profile_p_log(p, temp, linestyle="solid", color=sm.to_rgba(i))
ax.set_xlabel(r"$T$ / K")
ax.set_ylabel(r"$p$ / hPa")

Text(0, 0.5, '$p$ / hPa')
_images/13eaa4827126685afb1873620adb52861d10a9019d9b36e20e86a73638d64d36.png

Restart from a stored model state#

It is possible to initialize the atmospheric and surface state from a previous run.

ncfile = "my_rce_output.nc"
atmosphere = konrad.atmosphere.Atmosphere.from_netcdf(ncfile)
surface = konrad.surface.FixedTemperature.from_netcdf(ncfile)

This allows the user to spin-up the model to a reference state before performing a perturbation. For example, one can double the CO2 concentration in an otherwise equilibrated atmosphere to observe the stratospheric temperature adjustment.

# Plot reference state
fig, ax = plt.subplots()
plots.profile_p_log(atmosphere["plev"], atmosphere["T"][-1])
ax.set_xlabel(r"$T$ / K")
ax.set_ylabel(r"$p$ / hPa")

atmosphere["CO2"][:] *= 2. # Double the CO2 concentration

rce = konrad.RCE(
    atmosphere,
    surface=surface,
    timestep='12h',  # Set timestep in model time.
    max_duration='100d',  # Set maximum runtime.
)
rce.run()

# Plot adjusted state
plots.profile_p_log(rce.atmosphere["plev"], rce.atmosphere["T"][-1])

[<matplotlib.lines.Line2D at 0x7f05adf01dd0>]
_images/4882a812bd754d055b0414437f8974779d602432e4d0db308b2c66c01e1aa076.png
Contents