"""Calculates CO2 response."""
import logging
from typing import cast
import numpy as np
from .construct_conc import calc_inv_sums, interp_bg_conc
from .utils import tgco2_to_tgc
logger = logging.getLogger(__name__)
# CONSTANTS
#
# alpha_j (list): [ppbv/Tg(C)] alpha_j coefficients of the impulse response function G_C
# for CO2 (e.g. from Table 1 of Sausen & Schumann (2000))
# with j being the eigenmodes of the impulse response
# https://doi.org/10.1023/A:1005579306109
ALPHA_ARR = [0.067, 0.1135, 0.152, 0.0970, 0.041]
# m_j (list): [1/yr] inverse of tau_j coefficients of the impulse response function G_C
# for CO2 (e.g. from Table I of Sausen & Schumann (2000))
M_ARR = [0.0, 1.0 / 313.8, 1.0 / 79.8, 1.0 / 18.8, 1.0 / 1.7]
# pre-industrial CO2 and N2O concentration used as reference
# these values are from the SSP scenarios
CO2_0 = 284.32 # [ppm]
N2O_0 = 273.87 # [ppb]
CO2_VALIDITY_MIN = 180.0 # lower bound of Etminan et al. (2016) validity range [ppm]
CO2_VALIDITY_MAX = 2000.0 # upper bound of Etminan et al. (2016) validity range [ppm]
[docs]
def get_co2_emissions(inv_dict: dict) -> dict:
"""Get total CO2 emissions in Tg for each inventory.
Args:
inv_dict (dict): Dictionary of emission inventory xarrays,
keys are inventory years
Returns:
dict: Dictionary with array of CO2 emissions in Tg
"""
# Sum up CO2 emissions in inventories, converted to Tg
_inv_years, emis_co2 = calc_inv_sums("CO2", inv_dict, target_units="Tg")
emis_co2_dict = {"CO2": emis_co2}
return emis_co2_dict
[docs]
def greens_c(time: float) -> float:
"""Green's function / Impulse response for CO2 concentration.
After equation (5) in Sausen & Schumann (2000).
Args:
time (float): Time
Returns:
float: Impulse response for a certain point in time
"""
return np.sum(ALPHA_ARR * np.exp(np.multiply(M_ARR, -time)))
[docs]
def calc_co2_concentration(config: dict, emis_dict: dict) -> dict:
"""Calculates the CO2 concentration values in ppmv for emitted CO2 in Tg.
Gets the method from config and executes the corresponding subroutine.
Args:
config (dict): Configuration dictionary from config
emis_dict (dict): Dictionary with arrays of emissions
for time range as defined in config, keys are species
Returns:
dict: Dictionary with array of CO2 concentration in ppmv
for time range as defined in config, key is species CO2
"""
method = config["responses"]["CO2"]["conc"]["method"]
if method == "Sausen&Schumann":
conc_co2_dict = calc_co2_ss(config, emis_dict)
return conc_co2_dict
else:
raise ValueError("CO2.conc.method in config file is invalid.")
[docs]
def calc_co2_ss(config: dict, emis_dict: dict) -> dict:
"""Calculates the CO2 concentration values in ppmv for emitted CO2 in Tg.
After Sausen & Schumann (2000), formulas (4) and (5).
Args:
config (dict): Configuration dictionary from config
emis_dict (dict): Dictionary with arrays of emissions
for time range as defined in config, keys are species
Returns:
dict: Dictionary with array of CO2 concentration in ppmv
for time range as defined in config, key is species CO2
"""
time_config = config["time"]["range"]
time_range = np.arange(time_config[0], time_config[1], time_config[2], dtype=int)
delta_t = time_config[2]
# Convert Tg CO2 to Tg C
emis_co2_arr = cast(np.ndarray, tgco2_to_tgc(emis_dict["CO2"]))
conc_co2_arr = np.zeros(len(time_range))
i = 0
for year in time_range:
j = 0
conc_co2 = 0
for year_dash in time_range[: (i + 1)]:
if emis_co2_arr[j] != 0.0: # optimize code
# (4) in Sausen & Schumann, 2000
conc_co2 = (
conc_co2 + greens_c(year - year_dash) * emis_co2_arr[j] * delta_t
)
j = j + 1
conc_co2_arr[i] = conc_co2 / 1000.0 # convert ppbv -> ppmv
i = i + 1
return {"CO2": conc_co2_arr}
[docs]
def calc_co2_rf(conc_dict: dict, config: dict) -> dict:
"""Calculates the radiative forcing values for emitted CO2 concentrations.
The CO2 method is taken from the config file.
Args:
conc_dict (dict): Dictionary with array of concentrations (not including
background) between the starting and ending years, keys is species
config (dict): Configuration dictionary from config
Raises:
ValueError: if CO2.rf.method not valid
Returns:
dict: Dictionary with :class:`numpy.ndarray` of CO2 radiative forcing
values between the starting and ending years, key is species CO2
"""
method = config["responses"]["CO2"]["rf"]["method"]
if method == "IPCC_2001_1":
rf_dict = calc_co2_rf_ipcc_2001_1(conc_dict)
return rf_dict
if method == "IPCC_2001_2":
rf_dict = calc_co2_rf_ipcc_2001_2(conc_dict)
return rf_dict
if method == "IPCC_2001_3":
rf_dict = calc_co2_rf_ipcc_2001_3(conc_dict)
return rf_dict
if method == "Etminan_2016":
conc_n2o_bg_dict = interp_bg_conc(config, "N2O")
rf_dict = calc_co2_rf_etminan_2016(conc_dict, conc_n2o_bg_dict)
return rf_dict
# unknown CO2 RF method
raise ValueError("CO2.rf.method in config file is invalid.")
[docs]
def calc_co2_rf_ipcc_2001_1(conc_dict: dict) -> dict:
"""Calculates the radiative forcing values for emitted CO2 concentrations.
After IPCC 2001, Table 6.2, first row.
Args:
conc_dict (dict): Dictionary with array of concentrations (not including
background) between the starting and ending years, keys is species
Returns:
dict: Dictionary with array of CO2 radiative forcing values
between the starting and ending years, key is species CO2
"""
conc_co2_arr = conc_dict["CO2"]
rf_co2_arr = 5.35 * np.log(1 + conc_co2_arr / CO2_0)
return {"CO2": rf_co2_arr}
[docs]
def calc_co2_rf_ipcc_2001_2(conc_dict: dict) -> dict:
"""Calculates the radiative forcing values for emitted CO2 concentrations.
After IPCC 2001, Table 6.2, second row.
Args:
conc_dict (dict): Dictionary with array of concentrations (not including
background) between the starting and ending years, keys is species
Returns:
dict: Dictionary with array of CO2 radiative forcing values
between the starting and ending years, key is species CO2
"""
conc_co2_arr = conc_dict["CO2"]
rf_co2_arr = 4.841 * np.log(1 + conc_co2_arr / CO2_0) + 0.0906 * (
np.sqrt(conc_co2_arr + CO2_0) - np.sqrt(CO2_0)
)
return {"CO2": rf_co2_arr}
[docs]
def calc_co2_rf_ipcc_2001_3(conc_dict: dict) -> dict:
"""Calculates the radiative forcing values for emitted CO2 concentrations.
After IPCC 2001, Table 6.2, third row.
Args:
conc_dict (dict): Dictionary with array of concentrations (not including
background) between the starting and ending years, keys is species
Returns:
dict: Dictionary with array of CO2 radiative forcing values
between the starting and ending years, key is species CO2
"""
def g(conc: float | np.ndarray) -> float | np.ndarray:
return np.log(1.0 + 1.2 * conc + 0.005 * conc**2 + 1.4e-6 * conc**3)
conc_co2_arr = conc_dict["CO2"]
rf_co2_arr = 3.35 * (g(conc_co2_arr + CO2_0) - g(CO2_0))
return {"CO2": rf_co2_arr}
[docs]
def calc_co2_rf_etminan_2016(conc_dict: dict, conc_n2o_bg_dict: dict) -> dict:
"""Calculates the CO2 radiative forcing values for emitted concentrations.
Reference: Etminan, M., Myhre, G., Highwood, E. J., & Shine, K. P. (2016).
Radiative forcing of carbon dioxide, methane, and nitrous oxide: A
significant revision of the methane radiative forcing. Geophysical
Research Letters, 43(24), 12-614.
https://doi.org/10.1002/2016GL071930.
Args:
conc_dict (dict): Dictionary with array of concentrations (not including
background) between the starting and ending years, keys is species
conc_n2o_bg_dict (dict): Dictionary of :class:`numpy.ndarray` of
background N2O concentrations between the starting and ending
years, key is species
Returns:
dict: Dictionary with :class:`numpy.ndarray` of CO2 radiative forcing
values between the starting and ending years, key is species CO2
"""
# concentrations
d_co2_conc = conc_dict["CO2"] # ΔCO2 concentration (compared to background)
co2_conc = d_co2_conc + CO2_0 # CO2 concentration on top of background
n2o_conc = conc_n2o_bg_dict["N2O"]
n2o_conc_mean = 0.5 * (n2o_conc + N2O_0)
# check validity range: 180-2000 ppm for CO2 from Etminan et al. (2016)
if np.any((co2_conc < CO2_VALIDITY_MIN) | (co2_conc > CO2_VALIDITY_MAX)):
logger.warning(
"CO2 concentration is outside of the validity range 180 - 2000 ppm"
"given by Etminan et al. (2016)."
)
# coefficients
a1 = -2.4e-7 # W/m²/ppm
b1 = 7.2e-4 # W/m²/ppm
c1 = -2.1e-4 # W/m²/ppb
# calculate RF
x = a1 * d_co2_conc**2.0 + b1 * np.abs(d_co2_conc) + c1 * n2o_conc_mean + 5.36
y = np.log(co2_conc / CO2_0)
rf_co2_arr = x * y
return {"CO2": rf_co2_arr}
[docs]
def calc_co2_drf_dconc(conc_dict: dict, config: dict) -> dict:
"""Calculates the derivative of CO2 radiative forcing w.r.t. concentration.
This is used for the differential and marginal RF attribution methods.
The CO2 method is taken from the config file.
Args:
conc_dict (dict): Dictionary with array of concentrations
between the starting and ending years, key is species
config (dict): Configuration dictionary from config
Raises:
ValueError: if CO2.rf.method not valid or the derivative undefined
Returns:
dict: Dictionary with :class:`numpy.ndarray` of CO2 radiative forcing
derivative values between the starting and ending years, key is
species CO2
"""
method = config["responses"]["CO2"]["rf"]["method"]
if method == "Etminan_2016":
conc_n2o_bg_dict = interp_bg_conc(config, "N2O")
drf_dconc_dict = calc_co2_drf_dconc_etminan_2016(conc_dict, conc_n2o_bg_dict)
return drf_dconc_dict
raise ValueError(
"CO2.rf.method does not have a valid concentration derivative method"
"and thus cannot be used with the selected attribution method."
)
[docs]
def calc_co2_drf_dconc_etminan_2016(conc_dict: dict, conc_n2o_bg_dict: dict) -> dict:
"""Calculates the derivative of CO2 radiative forcing w.r.t. concentration.
Reference: Etminan, M., Myhre, G., Highwood, E. J., & Shine, K. P. (2016).
Radiative forcing of carbon dioxide, methane, and nitrous oxide: A
significant revision of the methane radiative forcing. Geophysical
Research Letters, 43(24), 12-614.
https://doi.org/10.1002/2016GL071930.
Args:
conc_dict (dict): Dictionary with array of concentrations (not including
background) between the starting and ending years, keys is species
conc_n2o_bg_dict (dict): Dictionary of :class:`numpy.ndarray` of
background N2O concentrations between the starting and ending
years, key is species
Returns:
dict: Dictionary with :class:`numpy.ndarray` of dRF(CO2)/dconc values
between the starting and ending years, key is species CO2
"""
# concentrations
d_co2_conc = conc_dict["CO2"] # ΔCO2 concentration (compared to background)
co2_conc = d_co2_conc + CO2_0 # CO2 concentration on top of background
n2o_conc = conc_n2o_bg_dict["N2O"]
n2o_conc_mean = 0.5 * (n2o_conc + N2O_0)
# coefficients
a1 = -2.4e-7 # W/m²/ppm
b1 = 7.2e-4 # W/m²/ppm
c1 = -2.1e-4 # W/m²/ppb
# calculate derivative of RF w.r.t. concentration using product rule
x = a1 * d_co2_conc**2.0 + b1 * np.abs(d_co2_conc) + c1 * n2o_conc_mean + 5.36
x_prime = 2.0 * a1 * d_co2_conc + b1 * np.sign(d_co2_conc)
y = np.log(co2_conc / CO2_0)
y_prime = 1.0 / co2_conc
drf_dconc_dict = x * y_prime + x_prime * y
return {"CO2": drf_dconc_dict}