Analysis of E3SMv2 Model Output#

Overview#

This workflow example showcases how to use UXarray to analyze the unstructured grid output from the Energy Exascale Earth System Model (E3SM) directly without needing to perform any regridding operations.

Imports#

This notebook requires the following packages to be installed in the notebook environment. cmocean is downloaded for pretty colormaps.

mamba install -c conda-forge uxarray cmocean
import warnings

warnings.filterwarnings("ignore")

from dask.distributed import Client, LocalCluster
import uxarray as ux
import cartopy.crs as ccrs
import geoviews as gv
import geoviews.feature as gf
import holoviews as hv
import cmocean
import glob

import sys

import panel as pn
import dask.distributed

# Set-up for HoloViz plots
hv.extension("matplotlib")
features = (
    gf.coastline(scale="110m", projection=ccrs.PlateCarree())
    * gf.borders(scale="110m", projection=ccrs.PlateCarree())
    * gf.states(scale="110m", projection=ccrs.PlateCarree())
)

Set up local Dask cluster#

For details about setting up Dask cluster, check out the Load Input Data in Parallel with Dask and UXarray page in the User Guide section.

# cluster = LocalCluster()
# client = Client(cluster)
client = dask.distributed.Client(n_workers=4, threads_per_worker=2)
client

Client

Client-2591eb31-909c-11ef-9aa8-0040a687eb14

Connection method: Cluster object Cluster type: distributed.LocalCluster
Dashboard: https://jupyterhub.hpc.ucar.edu/stable/user/rtam/proxy/8787/status

Cluster Info

Data#

Data loaded in this notebook is the output from E3SMv2 that is run for 6 simulated years. The case is an atmosphere-only (AMIP) simulation with present-day control forcing (F2010) at a 1-degree horizontal resolution (ne30pg2), where boundary conditions like sea surface temperatures and sea ice set as default as in the E3SMv2 model and repeat every simulated year.

%%time
# Load multiple files with chunking by time #CHUNK FAIL WITH E3SM IN TIME DIMENSION
data_files = "/glade/campaign/cisl/vast/uxarray/data/e3sm_keeling/ENSO_ctl_1std/unstructured/*.nc"
grid_file = (
    "/glade/campaign/cisl/vast/uxarray/data/e3sm_keeling/E3SM_grid/ne30pg2_grd.nc"
)
uxds_e3sm_multi = ux.open_mfdataset(grid_file, glob.glob(data_files), parallel=True)
CPU times: user 11.9 s, sys: 853 ms, total: 12.7 s
Wall time: 25.9 s
uxds_e3sm_multi
<xarray.UxDataset> Size: 37GB
Dimensions:              (time: 72, n_face: 21600, lev: 72, ilev: 73,
                          cosp_prs: 7, nbnd: 2, cosp_tau: 7, cosp_ht: 40,
                          cosp_sr: 15, cosp_htmisr: 16, cosp_tau_modis: 7,
                          cosp_reffice: 6, cosp_reffliq: 6, cosp_sza: 5,
                          cosp_scol: 10)
Coordinates: (12/13)
  * lev                  (lev) float64 576B 0.1238 0.1828 0.2699 ... 993.8 998.5
  * ilev                 (ilev) float64 584B 0.1 0.1477 0.218 ... 997.0 1e+03
  * cosp_prs             (cosp_prs) float64 56B 9e+04 7.4e+04 ... 2.45e+04 9e+03
  * cosp_tau             (cosp_tau) float64 56B 0.15 0.8 2.45 ... 41.5 100.0
  * cosp_scol            (cosp_scol) int32 40B 1 2 3 4 5 6 7 8 9 10
  * cosp_ht              (cosp_ht) float64 320B 1.896e+04 1.848e+04 ... 240.0
    ...                   ...
  * cosp_sza             (cosp_sza) float64 40B 0.0 20.0 40.0 60.0 80.0
  * cosp_htmisr          (cosp_htmisr) float64 128B 0.0 250.0 ... 1.8e+04
  * cosp_tau_modis       (cosp_tau_modis) float64 56B 0.15 0.8 ... 41.5 100.0
  * cosp_reffice         (cosp_reffice) float64 48B 5e-06 1.5e-05 ... 7.5e-05
  * cosp_reffliq         (cosp_reffliq) float64 48B 4e-06 9e-06 ... 2.5e-05
  * time                 (time) object 576B 0001-02-01 00:00:00 ... 0007-01-0...
Dimensions without coordinates: n_face, nbnd
Data variables: (12/471)
    lat                  (time, n_face) float64 12MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    lon                  (time, n_face) float64 12MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    area                 (time, n_face) float64 12MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    hyam                 (time, lev) float64 41kB dask.array<chunksize=(1, 72), meta=np.ndarray>
    hybm                 (time, lev) float64 41kB dask.array<chunksize=(1, 72), meta=np.ndarray>
    P0                   (time) float64 576B 1e+05 1e+05 1e+05 ... 1e+05 1e+05
    ...                   ...
    soa_c1DDF            (time, n_face) float32 6MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    soa_c1SFWET          (time, n_face) float32 6MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    soa_c2DDF            (time, n_face) float32 6MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    soa_c2SFWET          (time, n_face) float32 6MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    soa_c3DDF            (time, n_face) float32 6MB dask.array<chunksize=(1, 21600), meta=np.ndarray>
    soa_c3SFWET          (time, n_face) float32 6MB dask.array<chunksize=(1, 21600), meta=np.ndarray>

Calculating Cloud Radiative Effect (netCRE) with UXarray#

Cloud radiative effect can be calculated as follows:

Shortwave cloud radiative effect can be approximated as the difference beteween all-sky net shortwave flux (FSNT) at the top of the model and the clear-sky net shortwave flux (FSNTC).

\[ SWCRE = FSNT - FSNTC\]

Longwave cloud radiative effect is similar to that for SWCRE, but the all-sky and clear-sky longwave fluxes are applied instead.

\[ LWCRE = FLUT - FLUT\]

Net cloud radiative effect is thus the difference beteween shortwave and longwave cloud radiative effect.

\[ netCRE = SWCRE - LWCRE\]

Calculate Shortwave Cloud Radiative Effect (SWCRE)#

All of the following dataarrays are lazily loaded, as demonstrated below.

# Define function to display variable size in MiB
def var_size(in_var):
    result = sys.getsizeof(in_var) / 1024 / 1024
    return result


# Loading Variables
FSNT = uxds_e3sm_multi["FSNT"]
FSNTC = uxds_e3sm_multi["FSNTC"]

SWCRE = FSNT - FSNTC
SWCRE.name = "SWCRE"  # UXarray requires dataarray to have a name

print(f"Memory size of SWCRE dask array in MB  : {var_size(SWCRE):.2f} MiB")
Memory size of SWCRE dask array in MB  : 0.00 MiB
%output size=450
hv.Layout(
    SWCRE.isel(time=0).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        fig_size=800,
        title="SWCRE in Year 1 Jan",
        clim=(-100, 0),
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + SWCRE.isel(time=3).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        fig_size=800,
        title="SWCRE in Year 1 April",
        clim=(-100, 0),
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + SWCRE.isel(time=6).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        clim=(-100, 0),
        fig_size=800,
        title="SWCRE in Year 1 July",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + SWCRE.isel(time=9).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        clim=(-100, 0),
        fig_size=800,
        title="SWCRE in Year 1 Nov",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
).cols(2)

Calculate Longwave Cloud Radiative Effect (LWCRE)#

FLUT = uxds_e3sm_multi["FLUT"]
FLUTC = uxds_e3sm_multi["FLUTC"]

LWCRE = FLUT - FLUTC
LWCRE.name = "LWCRE"  # UXarray requires dataarray to have a name
%output size=450
hv.Layout(
    LWCRE.isel(time=0).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        fig_size=800,
        title="LWCRE in Year 1 Jan",
        clim=(-100, 0),
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + LWCRE.isel(time=3).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        fig_size=800,
        clim=(-100, 0),
        title="LWCRE in Year 1 April",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + LWCRE.isel(time=6).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        clim=(-100, 0),
        fig_size=800,
        title="LWCRE in Year 1 July",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + LWCRE.isel(time=9).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        clim=(-100, 0),
        fig_size=800,
        title="LWCRE in Year 1 Nov",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
).cols(2)

Calculate Net CRE#

netCRE = SWCRE - LWCRE
netCRE.name = "netCRE"  # UXarray requires dataarray to have a name
%output size=450
hv.Layout(
    netCRE.isel(time=0).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        fig_size=800,
        title="netCRE in Year 1 Jan",
        clim=(-100, 0),
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + netCRE.isel(time=3).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        fig_size=800,
        clim=(-100, 0),
        title="netCRE in Year 1 April",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + netCRE.isel(time=6).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        clim=(-100, 0),
        fig_size=800,
        title="netCRE in Year 1 July",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
    + netCRE.isel(time=9).plot(
        backend="matplotlib",
        cmap=cmocean.cm.haline,
        pixel_ratio=2.0,
        clim=(-100, 0),
        fig_size=800,
        title="netCRE in Year 1 Nov",
        clabel="(W/m2)",
        fontscale=2.5,
    )
    * features.opts(fig_size=800)
).cols(2)