kemccusker / scmcoat

lightweight wrapper around simple climate models (namely, FaIR)

scmcoat : Simple Climate Model wrapper

Wrapper tool to run simple climate models with CIL settings and climate uncertainty.

SCMs included and tested:

• FaIR v1.*

FaIR accepts an array of greenhouse gas emissions and prognostically predicts atmospheric concentrations, change in radiative forcings, and change in global mean temperature, and other variables under certain settings. We frequently run this model under very specific settings and this package is meant to ease that process.

This package is very new. This may change in breaking ways.

Inputs

Emissions inputs to FaIR are a specific format. FaIR v1.* expects 39-species of GHGs if run in MultiGas mode.

Demo emissions (see code example below) are an xarray.DataArray with dimensions year x gas. Note that there are 40 elements in gas: the first element is "year", and the rest are gases.

The `gas` coordinate:
array(['year', 'CO2_Fossil', 'CO2_Land', 'CH4', 'N2O', 'SOx', 'CO', 'NMVOC',
       'NOx', 'BC', 'OC', 'NH3', 'CF4', 'C2F6', 'C6F14', 'HFC23', 'HFC32',
       'HFC4310mee', 'HFC125', 'HFC134a', 'HFC143a', 'HFC227ea', 'HFC245fa',
       'SF6', 'CFC11', 'CFC12', 'CFC113', 'CFC114', 'CFC115', 'CARB_TET',
       'MCF', 'HCFC22', 'HCFC141b', 'HCFC142b', 'Halon1211', 'Halon1202',
       'Halon1301', 'Halon2402', 'CH3Br', 'CH3Cl'], dtype='<U10')
The `year` coordinate:      
array([1765., 1766., 1767., ..., 2498., 2499., 2500.]) # year runs from 1765-2500 or 1751-2500 

TODO: add more detail about year coord

Demo emissions are SSP2-4.5 from the CMIP6 era. Standardized CMIP6 emissions can be obtained like:

import scmcoat as sc

# demo emissions
emissions = sc.get_test_emissions() # returns xr.DataArray of SSP2-4.5

# or get any CMIP6 emissions pathway
emissions = sc.utils.rcmip_emissions("ssp370") # returns xr.DataArray

FairModel.run() will also accept emissions that are pandas.DataFrame with columns that match the gas coord and index that matches year, as described above. E.g. a 40-element gas dimension and a 736- or 751-element year dimension.

Outputs

FaIR outputs are returned as an xarray.Dataset with the following data_vars:

Data variables:
    concentration  (year, gas) float64 278.0 722.0 273.0 ... 13.22 620.6
    forcing        (year, forcing_type) float64 2.623e-05 0.0 0.0 ... 0.0 0.0
    temperature    (year) float64 0.005061 0.009262 0.01363 ... 3.244 3.245
    ocean_heat_content (year) float64 4.358e+21 7.791e+21 ... 5.822e+24
    simulation     () <U7 'default'

Here a simulation equal to "default" indicates it's not running an ensemble of simulations with climate parameters but rather is running with default FaIR settings.

Examples:

Run v1.* FaIR with demo emissions (CMIP6 SSP2-4.5)

import scmcoat as sc

# initialize the model
fm = sc.FairModel()

# Open up the test emissions to use as a demo
emissions = fm.get_test_emissions()

# run FaIR with these emissions under its default settings
response = fm.run(emissions)  # an xr.Dataset

response.temperature.plot()

OR

import scmcoat as sc

# initialize the model
fm = sc.FairModel()

# run FaIR with test emissions under its default settings
response = fm.test_run()  # an xr.Dataset

response.temperature.plot()

Run v1.* FaIR with demo emissions (CMIP6 SSP2-4.5) and climate parameters

When FairModel.run() is called, the argument,simid, specifies how to use the climate parameters. The simid default value is "default".

• If ClimateParams have not been set and simid is "default", then FairModel is run with FaIR's default settings (see above example).
• If ClimateParams have been set but simid is "default", then FairModel will be run with "median" climate parameters. This is our usual "point estimate" setup.
• If ClimateParams have been set, then simid can be set to "median" or a scalar integer that indexes the climate parameters.

import scmcoat as sc

params_filepath = "<path_to_params_file.nc>" # TODO add pointer to the file or add to repo

# set up ClimateParams
cp = sc.core.ClimateParams(params=xr.open_dataset(params_filepath))

# initialize the model with climate parameters
fm = sc.core.FairModel(cp)

# Open up the test emissions to use as a demo
emissions = fm.get_test_emissions()

# run FaIR with these emissions and median climate parameters
response = fm.run(emissions, simid="median")  # an xr.Dataset

response.temperature.plot()

# the model can be run with different climate parameters
response2 = fm.run(emissions, simid=150)

_ = xr.concat([response,response2], dim="simulation").temperature.plot.line(x="year")

FaIR can run any emissions pathway in this format, although may become unstable under radically different emissions. Note this has not been tested with different length time series of emissions and will likely throw errors.

scmcoat's Issues

Error when running with test or built-in emissions and "default" `simid`

If climate params are not set, and test emissions or built-in CMIP6 emissions are used to run, fair throws an exception:

sc.core.FairModel(
    debug=False
).run(emissda=sc.utils.rcmip_emissions("ssp245"))

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[10], line 1
----> 1 sc.core.FairModel(
      2     debug=False).run(emissda=emissionsda)#,simid="median")

File ~/code/scmcoat/src/scmcoat/core.py:278, in FairModel.run(self, emissda, simid, emissions_driven, useMultigas)
    275 else:
    276     raise NotImplementedError("emissions time dimension is not recognized")
--> 278 ret = self._run(emiss=emiss.values, 
    279                 useMultigas=useMultigas)
    281 # Prep the FaIR outputs into xarray objects
    282 
    283 # TODO generalize the years
    284 years = np.arange(reference_year, 2501)  # gets range of years in desired period

File ~/code/scmcoat/src/scmcoat/core.py:147, in FairModel._run(self, emiss, useMultigas)
    145         print(f"Running default FaIR, v{fair.__version__}")
    146     # ignore simid
--> 147     ret = fair.forward.fair_scm(emissions=emiss, useMultigas=useMultigas)
    148     return ret
    149 else:

File /srv/conda/envs/notebook/lib/python3.10/site-packages/fair/forward.py:249, in fair_scm(emissions, emissions_driven, C, other_rf, q, tcrecs, d, F2x, tcr_dbl, a, tau, r0, rc, rt, iirf_max, iirf_h, C_pi, E_pi, restart_in, restart_out, F_tropO3, F_aerosol, F_volcanic, F_solar, F_contrails, F_bcsnow, F_landuse, aviNOx_frac, F_ref_aviNOx, E_ref_aviNOx, F_ref_BC, E_ref_BC, fossilCH4_frac, natural, efficacy, scale, oxCH4_frac, ghg_forcing, scale_F2x, stwv_from_ch4, b_aero, b_tro3, pi_tro3, ghan_params, stevens_params, ref_isSO2, useMultigas, tropO3_forcing, ozone_feedback, lifetimes, aerosol_forcing, scaleAerosolAR5, fixPre1850RCP, useTropO3TFeedback, scaleHistoricalAR5, contrail_forcing, kerosene_supply, landuse_forcing, aCO2land, ariaci_out, bcsnow_forcing, diagnostics, gir_carbon_cycle, temperature_function, lambda_global, ocean_heat_capacity, ocean_heat_exchange, deep_ocean_efficacy)
    247     elif natural.ndim==2:
    248         if natural.shape[1]!=2 or natural.shape[0]!=nt:
--> 249             raise ValueError(
    250               "natural emissions should be a 2-element or nt x 2 " +
    251               "array")
    252 else:
    253     raise ValueError(
    254       "natural emissions should be a scalar, 2-element, or nt x 2 " +
    255       "array")

ValueError: natural emissions should be a 2-element or nt x 2 array

Generalize emissions time dimension

Currently, the FairModel._run() function in src/core.py will work with emissions that are of specific lengths (assuming they are standard CMIP5 or CMIP6 emissions pathways). This should be generalized to work with other emissions pathway lengths if the FaIR model itself can handle it.

Add prep of AR6 concentration to utils

Right now, utils.py only has the functionality to obtain and prepare AR6 scenario emissions, but code exists to also obtain and prep concentrations. We don't have a use for concentrations at the moment, but this would be nice to have.

Update test emissions

Test emissions are currently rcp45 emissions from CMIP5. They were prepared straight from fair v1.6.2, however this era of emissions does not work well with the AR6 climate parameters.

update climate parameter handling to correctly use trimmed down JSON files

The current approach relies on the user having access to a preprocessed climate parameter netcdf file (created from v1.0 JSON file of climate parameters).

PR #14 outlined a way to use the trimmed down v1.1 JSON files but left the implementation of preparing the params for a later PR. In looking into that implementation, I discovered that some v1.1 parameters have a different structure than the v1.0 parameters and so the implementation is not as straightforward as copying the code used to generate the netcdf file.

This Issue is meant to capture adding an implementation that can prepare the climate parameters xr.Dataset from v1.1 JSON files.