Skip to content

Commit

Permalink
adjust remaining aerosol and volcanic RFs and add to their tests
Browse files Browse the repository at this point in the history
  • Loading branch information
ptrscll committed Jul 2, 2024
1 parent b58ecf6 commit 17eca52
Show file tree
Hide file tree
Showing 3 changed files with 72 additions and 14 deletions.
17 changes: 13 additions & 4 deletions src/forcing_component.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -443,27 +443,35 @@ void ForcingComponent::run(const double runToDate) {
double E_OC =
core->sendMessage(M_GETDATA, D_EMISSIONS_OC, message_data(runToDate))
.value(U_TG);
double foc = rho_oc * E_OC;
// The 0.2 value comes from equally distributing the alpha scalar to all 5
// aerosol RF types.
double foc = 0.2 * alpha.value(U_UNITLESS) * rho_oc * E_OC;
forcings[D_RF_OC].set(foc, U_W_M2);

// ---------- Sulphate Aerosols ----------
unitval SO2_emission = core->sendMessage(M_GETDATA, D_EMISSIONS_SO2,
message_data(runToDate));
double fso2 = rho_so2 * SO2_emission.value(U_GG_S);
// The 0.2 value comes from equally distributing the alpha scalar to all 5
// aerosol RF types.
double fso2 = 0.2 * alpha.value(U_UNITLESS) * rho_so2 * SO2_emission.value(U_GG_S);
forcings[D_RF_SO2].set(fso2, U_W_M2);

// ---------- NH3 ----------
double E_NH3 =
core->sendMessage(M_GETDATA, D_EMISSIONS_NH3, message_data(runToDate))
.value(U_TG);
double fnh3 = rho_nh3 * E_NH3;
// The 0.2 value comes from equally distributing the alpha scalar to all 5
// aerosol RF types.
double fnh3 = 0.2 * alpha.value(U_UNITLESS) * rho_nh3 * E_NH3;
forcings[D_RF_NH3].set(fnh3, U_W_M2);

// ---------- RFaci ----------
// ERF from aerosol-cloud interactions (RFaci)
// Based on Equation 7.SM.1.2 from IPCC AR6 where
// The 0.2 value comes from equally distributing the alpha scalar to all 5
// aerosol RF types.
double aci_rf =
-1 * aci_beta *
0.2 * alpha.value(U_UNITLESS) * -1 * aci_beta *
log(1 + (SO2_emission / s_SO2) + ((E_BC + E_OC) / s_BCOC));
forcings[D_RF_ACI].set(aci_rf, U_W_M2);
}
Expand All @@ -478,6 +486,7 @@ void ForcingComponent::run(const double runToDate) {
if (core->checkCapability(D_VOLCANIC_SO2)) {
// The volcanic forcings are read in from an ini file.
forcings[D_RF_VOL] =
volscl.value(U_UNITLESS) *
core->sendMessage(M_GETDATA, D_VOLCANIC_SO2, message_data(runToDate));
}

Expand Down
2 changes: 1 addition & 1 deletion src/temperature_component.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -463,7 +463,7 @@ void TemperatureComponent::run(const double runToDate) {
const double ftot = core->sendMessage(M_GETDATA, D_RF_TOTAL, message_data(runToDate)).value(U_W_M2);
forcing[tstep] =
double(ftot) -
(1.0) * aero_forcing - (1.0) * volcanic_forcing;
(0.0) * aero_forcing - (0.0) * volcanic_forcing;

// Initialize variables for time-stepping through the model
double DQ1 = 0.0;
Expand Down
67 changes: 58 additions & 9 deletions tests/testthat/test_parameters.R
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,11 @@ sampledir <- system.file("output", package = "hector")
testvars <- c(CONCENTRATIONS_CO2(), RF_TOTAL(), GLOBAL_TAS())

dates <- 1750:2100

# Limit the test dates to the future, where the historical
# variability won't impact the temp.
tdates <- 2000:2100

ssp245 <- file.path(inputdir, "hector_ssp245.ini")


Expand Down Expand Up @@ -87,9 +92,6 @@ test_that("Lowering initial CO2 lowers output CO2", {
shutdown(hc)
})

# Limit the test dates to the future, where the historical
# variability won't impact the temp.
tdates <- 2000:2100

test_that("Lowering ECS lowers output Temperature", {

Expand Down Expand Up @@ -176,7 +178,7 @@ test_that("Lowering diffusivity increases temperature", {
test_that("Lowering aerosol forcing scaling factor increases temperature", {

# Relevant vars to save and test.
vars <- c(GLOBAL_TAS())
vars <- c(GLOBAL_TAS(), RF_BC(), RF_OC(), RF_SO2(), RF_NH3(), RF_ACI(), RF_TOTAL())

# Define Hector core.
hc <- newcore(ssp245, suppresslogging = TRUE)
Expand All @@ -193,11 +195,47 @@ test_that("Lowering aerosol forcing scaling factor increases temperature", {
setvar(hc, NA, AERO_SCALE(), new_alpha, getunits(AERO_SCALE()))
reset(hc, hc$reset_date)
run(hc, 2100)
dd2 <- fetchvars(hc, tdates, GLOBAL_TAS())
dd2 <- fetchvars(hc, tdates, vars)

# Check to make sure that the temp and RF have changed.
diff <- dd2$value - dd1$value
expect_gt(min(diff), 0.0)

# Checking that temperatures increase
temp1 <- dd1[dd1$variable == GLOBAL_TAS(),]
temp2 <- dd2[dd2$variable == GLOBAL_TAS(),]
temp_diff <- temp2$value - temp1$value
expect_gt(min(temp_diff), 0.0)

# Checking that black carbon RF decreases
bc1 <- dd1[dd1$variable == RF_BC(),]
bc2 <- dd2[dd2$variable == RF_BC(),]
bc_diff <- bc2$value - bc1$value
expect_lt(max(bc_diff), 0.0)

# Checking organic carbon RF magnitude decreases
oc1 <- dd1[dd1$variable == RF_OC(),]
oc2 <- dd2[dd2$variable == RF_OC(),]
oc_diff <- abs(oc2$value) - abs(oc1$value)
expect_lt(max(oc_diff), 0.0)

# Checking that NH3 RF increases (becomes less negative)
nh3_1 <- dd1[dd1$variable == RF_NH3(),]
nh3_2 <- dd2[dd2$variable == RF_NH3(),]
nh3_diff <- nh3_2$value - nh3_1$value
expect_gt(min(nh3_diff), 0.0)

# Checking that SO2 RF increases (becomes less negative)
so2_1 <- dd1[dd1$variable == RF_SO2(),]
so2_2 <- dd2[dd2$variable == RF_SO2(),]
so2_diff <- so2_2$value - so2_1$value
expect_gt(min(so2_diff), 0.0)

# Checking that ACI RF increases (becomes less negative)
aci1 <- dd1[dd1$variable == RF_ACI(),]
aci2 <- dd2[dd2$variable == RF_ACI(),]
aci_diff <- aci2$value - aci1$value
expect_gt(min(aci_diff), 0.0)



shutdown(hc)
})
Expand All @@ -207,7 +245,7 @@ test_that("Increasing volcanic forcing scaling factor increases the effect of vo
## Because the volcanic forcing scaling factor only has an impact during the
## the volcanic events. Only check the temp during those years.
tdates <- c(1960, 1965)
vars <- GLOBAL_TAS()
vars <- c(GLOBAL_TAS(), RF_VOL())

# Set up and run Hector
hc <- newcore(ssp245, suppresslogging = TRUE)
Expand All @@ -222,7 +260,18 @@ test_that("Increasing volcanic forcing scaling factor increases the effect of vo
run(hc)
new_out <- fetchvars(hc, tdates, vars)

expect_true(all(out$value != new_out$value))
# Getting temperature and RF data
temps <- out[out$variable == GLOBAL_TAS(),]
rfs <- out[out$variable == RF_VOL(),]

new_temps <- new_out[new_out$variable == GLOBAL_TAS(),]
new_rfs <- new_out[new_out$variable == RF_VOL(),]

expect_true(all(temps$value != new_temps$value))

# Volcanic RF should decrease (become more negative)
rf_diff <- new_rfs$value - rfs$value
expect_lt(max(rf_diff), 0.0)
})

test_that("Decreasing vegetation NPP fraction has down stream impacts", {
Expand Down

0 comments on commit 17eca52

Please sign in to comment.