Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[AUTO] Format files using DocumentFormat #17

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open

[AUTO] Format files using DocumentFormat #17

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 10 additions & 10 deletions src/MimiMAGICC.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,26 +9,26 @@ include("components/rf_ch4.jl")
include("components/rf_ch4h2o.jl")
include("components/rf_o3.jl")

function get_model(;rcp_scenario::String="RCP85", start_year::Int=1765, end_year::Int=2300)
function get_model(; rcp_scenario::String="RCP85", start_year::Int=1765, end_year::Int=2300)

# ---------------------------------------------
# Load and clean up necessary data.
# ---------------------------------------------

# Load RCP emissions and concentration scenario values (RCP options = "RCP26" or "RCP85").
rcp_emissions = DataFrame(load(joinpath(@__DIR__, "..", "data", rcp_scenario*"_EMISSIONS.csv"), skiplines_begin=36))
rcp_concentrations = DataFrame(load(joinpath(@__DIR__, "..", "data", rcp_scenario*"_CONCENTRATIONS.csv"), skiplines_begin=37))
rcp_emissions = DataFrame(load(joinpath(@__DIR__, "..", "data", rcp_scenario * "_EMISSIONS.csv"), skiplines_begin=36))
rcp_concentrations = DataFrame(load(joinpath(@__DIR__, "..", "data", rcp_scenario * "_CONCENTRATIONS.csv"), skiplines_begin=37))

# Load temperature scenario (mean response from SNEASY across 100,000 posterior parameter samples).
rcp_temperature = DataFrame(load(joinpath(@__DIR__, "..", "data", rcp_scenario*"_temperature_sneasy_1765_2300.csv")))
rcp_temperature = DataFrame(load(joinpath(@__DIR__, "..", "data", rcp_scenario * "_temperature_sneasy_1765_2300.csv")))

# Load fossil CO₂ emissions scenario to calculate historic tropospheric O₃ forcing (this scenario is specific to how MAGICC calculates this forcing).
foss_hist_for_O₃ = DataFrame(load(joinpath(@__DIR__, "..", "data", "FOSS_HIST_magicc.csv")))
foss_hist_for_O₃ = DataFrame(load(joinpath(@__DIR__, "..", "data", "FOSS_HIST_magicc.csv")))

# Find indices for start and end years to crop scenarios to correct time frame.
start_index, end_index = findall((in)([start_year, end_year]), rcp_emissions.YEARS)
rcp_emissions = rcp_emissions[start_index:end_index, :]
rcp_temperature = rcp_temperature[start_index:end_index, :]
rcp_emissions = rcp_emissions[start_index:end_index, :]
rcp_temperature = rcp_temperature[start_index:end_index, :]
foss_hist_for_O₃ = foss_hist_for_O₃[start_index:end_index, :]

# Set initial CH₄ and N₂O concentrations to RCP 1765 values.
Expand Down Expand Up @@ -97,11 +97,11 @@ function get_model(;rcp_scenario::String="RCP85", start_year::Int=1765, end_year
# ---------------------------------------------
# Create connections between Mimi components.
# ---------------------------------------------
connect_param!(m, :rf_ch4, :CH₄, :ch4_cycle, :CH₄)
connect_param!(m, :rf_ch4, :CH₄, :ch4_cycle, :CH₄)
connect_param!(m, :rf_ch4h2o, :CH₄, :ch4_cycle, :CH₄)
connect_param!(m, :rf_o3, :CH₄, :ch4_cycle, :CH₄)
connect_param!(m, :rf_o3, :CH₄, :ch4_cycle, :CH₄)

return(m)
return (m)
end

end #module
82 changes: 41 additions & 41 deletions src/components/ch4_cycle.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,34 +4,34 @@

@defcomp ch4_cycle begin

BBCH4 = Parameter() # Tg to ppb unit conversion between emissions and concentrations Note:(1 Teragram = 1 Megatonne).
ANOX = Parameter() # Hydroxyl lifetime coefficient for nitrogen oxides.
ACO = Parameter() # Hydroxyl lifetime coefficient for carbon monoxide.
AVOC = Parameter() # Hydroxyl lifetime coefficient for non-methane volatile organic compounds.
SCH4 = Parameter() # Hydroxyl lifetime coefficient for methane.
TAUINIT = Parameter() # Initial methane tropospheric lifetime due to hydroxyl abundance (years).
TAUSOIL = Parameter() # Methane soil sink lifetime (years).
TAUSTRAT = Parameter() # Methane stratospheric lifetime (years).
GAM = Parameter() # Hydroxyl lifetime parameter.
CH4_natural = Parameter() # Natural methane emissions (Mt yr⁻¹).
fffrac = Parameter() # Fossil fuel fraction of total methane emissions.
CH₄_0 = Parameter() # Atmospheric methane pre-industrial concentration (ppb).
temperature = Parameter(index=[time]) # Global surface temperature anomaly (°C).
CH4_emissions = Parameter(index=[time]) # Global anthropogenic methane emissions (Mt yr⁻¹).
NOx_emissions = Parameter(index=[time]) # Global nitrogen oxides emissions (Mt yr⁻¹).
CO_emissions = Parameter(index=[time]) # Global carbon monoxide emissions (Mt yr⁻¹).
NMVOC_emissions = Parameter(index=[time]) # Global non-methane volatile organic compound emissions (Mt yr⁻¹).

TAU_OH = Variable(index=[time]) # Methane tropospheric lifetime due to hydroxyl abundance (years).
CH₄ = Variable(index=[time]) # Atmospheric methane concetration for current period (ppb).
emeth = Variable(index=[time]) # Methane that has been oxidized to carbon dioxide.
BBCH4 = Parameter() # Tg to ppb unit conversion between emissions and concentrations Note:(1 Teragram = 1 Megatonne).
ANOX = Parameter() # Hydroxyl lifetime coefficient for nitrogen oxides.
ACO = Parameter() # Hydroxyl lifetime coefficient for carbon monoxide.
AVOC = Parameter() # Hydroxyl lifetime coefficient for non-methane volatile organic compounds.
SCH4 = Parameter() # Hydroxyl lifetime coefficient for methane.
TAUINIT = Parameter() # Initial methane tropospheric lifetime due to hydroxyl abundance (years).
TAUSOIL = Parameter() # Methane soil sink lifetime (years).
TAUSTRAT = Parameter() # Methane stratospheric lifetime (years).
GAM = Parameter() # Hydroxyl lifetime parameter.
CH4_natural = Parameter() # Natural methane emissions (Mt yr⁻¹).
fffrac = Parameter() # Fossil fuel fraction of total methane emissions.
CH₄_0 = Parameter() # Atmospheric methane pre-industrial concentration (ppb).
temperature = Parameter(index=[time]) # Global surface temperature anomaly (°C).
CH4_emissions = Parameter(index=[time]) # Global anthropogenic methane emissions (Mt yr⁻¹).
NOx_emissions = Parameter(index=[time]) # Global nitrogen oxides emissions (Mt yr⁻¹).
CO_emissions = Parameter(index=[time]) # Global carbon monoxide emissions (Mt yr⁻¹).
NMVOC_emissions = Parameter(index=[time]) # Global non-methane volatile organic compound emissions (Mt yr⁻¹).

TAU_OH = Variable(index=[time]) # Methane tropospheric lifetime due to hydroxyl abundance (years).
CH₄ = Variable(index=[time]) # Atmospheric methane concetration for current period (ppb).
emeth = Variable(index=[time]) # Methane that has been oxidized to carbon dioxide.


function run_timestep(p, v, d, t)

# Set initial conditions. NOTE: The CH₄ cycle equations require temperature values from the previous two timesteps (t-1 and t-2), so the cycle switches
# on starting in timestep 3. Further, timestep 1 CH₄ concentrations do not affect the results below (but timestep 2 CH₄ concentrations do). We
# therefore treat period 2 CH₄ concentrations as an initial (and uncertain) condition.
# on starting in timestep 3. Further, timestep 1 CH₄ concentrations do not affect the results below (but timestep 2 CH₄ concentrations do). We
# therefore treat period 2 CH₄ concentrations as an initial (and uncertain) condition.
if t <= TimestepIndex(2)
v.CH₄[TimestepIndex(2)] = p.CH₄_0
v.emeth[t] = 0.0
Expand All @@ -48,14 +48,14 @@

# Calculate change in emissions for NOx, CO, and NMVOC (relative to period 3 when CH₄ cycle model engages due to t-1 and t-2 terms above).
DENOX = p.NOx_emissions[t] - p.NOx_emissions[TimestepIndex(3)]
DECO = p.CO_emissions[t] - p.CO_emissions[TimestepIndex(3)]
DECO = p.CO_emissions[t] - p.CO_emissions[TimestepIndex(3)]
DEVOC = p.NMVOC_emissions[t] - p.NMVOC_emissions[TimestepIndex(3)]

# Add natural CH₄ emissions to anthropogenic CH₄ emissions.
Total_CH4emissions = p.CH4_emissions[t] + p.CH4_natural

# Calcualte TAUOTHER term (based on stratospheric and soil sink lifetimes).
TAUOTHER = 1.0 / (1.0/p.TAUSOIL + 1.0/p.TAUSTRAT)
TAUOTHER = 1.0 / (1.0 / p.TAUSOIL + 1.0 / p.TAUSTRAT)

##########################################################
# Iterate to Calculate OH lifetime and CH₄ concentration.
Expand All @@ -69,7 +69,7 @@
B00 = p.CH₄_0 * p.BBCH4

# Account for other gases on OH abundance and tropospheric methane lifetime.
AAA = exp(p.GAM * (p.ANOX*DENOX + p.ACO*DECO + p.AVOC*DEVOC))
AAA = exp(p.GAM * (p.ANOX * DENOX + p.ACO * DECO + p.AVOC * DEVOC))

# Multiply two parameters together.
X = p.GAM * p.SCH4
Expand All @@ -83,39 +83,39 @@
# FIRST ITERATION
#----------------------------------------
BBAR = B
TAUBAR = U * (BBAR/B00) ^ X
TAUBAR = p.TAUINIT / (p.TAUINIT/TAUBAR + 0.0316 * D_Temp)
DB1 = Total_CH4emissions - (BBAR/TAUBAR) - (BBAR/TAUOTHER)
TAUBAR = U * (BBAR / B00)^X
TAUBAR = p.TAUINIT / (p.TAUINIT / TAUBAR + 0.0316 * D_Temp)
DB1 = Total_CH4emissions - (BBAR / TAUBAR) - (BBAR / TAUOTHER)
B1 = B + DB1

#----------------------------------------
# SECOND ITERATION
#----------------------------------------
BBAR = (B + B1) / 2.0
TAUBAR = U * (BBAR/B00) ^ X
TAUBAR = TAUBAR * (1.0 - 0.5 * X * DB1/B)
TAUBAR = p.TAUINIT / (p.TAUINIT/TAUBAR + 0.0316 * D_Temp)
DB2 = Total_CH4emissions - (BBAR/TAUBAR) - (BBAR/TAUOTHER)
TAUBAR = U * (BBAR / B00)^X
TAUBAR = TAUBAR * (1.0 - 0.5 * X * DB1 / B)
TAUBAR = p.TAUINIT / (p.TAUINIT / TAUBAR + 0.0316 * D_Temp)
DB2 = Total_CH4emissions - (BBAR / TAUBAR) - (BBAR / TAUOTHER)
B2 = B + DB2

#----------------------------------------
# THIRD ITERATION
#----------------------------------------
BBAR = (B + B2) / 2.0
TAUBAR = U * (BBAR/B00) ^ X
TAUBAR = TAUBAR * (1.0 - 0.5 * X * DB2/B)
TAUBAR = p.TAUINIT / (p.TAUINIT/TAUBAR + 0.0316 * D_Temp)
DB3 = Total_CH4emissions - (BBAR/TAUBAR) - (BBAR/TAUOTHER)
TAUBAR = U * (BBAR / B00)^X
TAUBAR = TAUBAR * (1.0 - 0.5 * X * DB2 / B)
TAUBAR = p.TAUINIT / (p.TAUINIT / TAUBAR + 0.0316 * D_Temp)
DB3 = Total_CH4emissions - (BBAR / TAUBAR) - (BBAR / TAUOTHER)
B3 = B + DB3

#----------------------------------------
# FOURTH ITERATION
#----------------------------------------
BBAR = (B + B3) / 2.0
TAUBAR = U * (BBAR/B00) ^ X
TAUBAR = TAUBAR * (1.0 - 0.5 * X * DB3/B)
TAUBAR = p.TAUINIT / (p.TAUINIT/TAUBAR + 0.0316 * D_Temp)
DB4 = Total_CH4emissions - (BBAR/TAUBAR) - (BBAR/TAUOTHER)
TAUBAR = U * (BBAR / B00)^X
TAUBAR = TAUBAR * (1.0 - 0.5 * X * DB3 / B)
TAUBAR = p.TAUINIT / (p.TAUINIT / TAUBAR + 0.0316 * D_Temp)
DB4 = Total_CH4emissions - (BBAR / TAUBAR) - (BBAR / TAUOTHER)
B4 = B + DB4

# Calculate CH₄ tropospheric lifetime.
Expand Down
12 changes: 6 additions & 6 deletions src/components/rf_ch4.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,24 +4,24 @@

# Create a function 'interact' to account for the overlap in methane and nitrous oxide in their radiative absoprtion bands.
function interact(M, N)
d = 1.0 + 0.636 * ((M * N)/(10^6))^0.75 + 0.007 * (M/(10^3)) * ((M * N)/(10^6))^1.52
d = 1.0 + 0.636 * ((M * N) / (10^6))^0.75 + 0.007 * (M / (10^3)) * ((M * N) / (10^6))^1.52
return 0.47 * log(d)
end

@defcomp rf_ch4 begin
CH₄_0 = Parameter() # Pre-industrial atmospheric methane concentration (ppb).
N₂O_0 = Parameter() # Pre-industrial atmospheric nitrous oxide concentration (ppb).
CH₄_0 = Parameter() # Pre-industrial atmospheric methane concentration (ppb).
N₂O_0 = Parameter() # Pre-industrial atmospheric nitrous oxide concentration (ppb).
scale_CH₄ = Parameter() # Scaling factor for uncertainty in methane radiative forcing.
CH₄ = Parameter(index=[time]) # Atmospheric methane concentration (ppb).
CH₄ = Parameter(index=[time]) # Atmospheric methane concentration (ppb).

QMeth = Variable(index=[time]) # Direct forcing from atmospheric methane concentrations (Wm⁻²).
QMeth = Variable(index=[time]) # Direct forcing from atmospheric methane concentrations (Wm⁻²).


function run_timestep(p, v, d, t)

if is_first(t)
# Set initial forcings to 0.
v.QMeth[t] = 0.0
v.QMeth[t] = 0.0
else
# Direct radiative forcing from atmospheric CH₄ concentrations.
v.QMeth[t] = (0.036 * (sqrt(p.CH₄[t]) - sqrt(p.CH₄_0)) - (interact(p.CH₄[t], p.N₂O_0) - interact(p.CH₄_0, p.N₂O_0))) * p.scale_CH₄
Expand Down
8 changes: 4 additions & 4 deletions src/components/rf_ch4h2o.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,11 +4,11 @@

@defcomp rf_ch4h2o begin

CH₄_0 = Parameter() # Pre-industrial atmospheric methane concentration (ppb).
STRATH2O = Parameter() # Share of direct methane radiative forcing that represents stratospheric water vapor forcing from methane oxidation.
CH₄ = Parameter(index=[time]) # Atmospheric methane concentration (ppb).
CH₄_0 = Parameter() # Pre-industrial atmospheric methane concentration (ppb).
STRATH2O = Parameter() # Share of direct methane radiative forcing that represents stratospheric water vapor forcing from methane oxidation.
CH₄ = Parameter(index=[time]) # Atmospheric methane concentration (ppb).

QCH4H2O = Variable(index=[time]) # Radiative forcing for stratoshperic water vapor from methane oxidation (Wm⁻²). rf_O₃ = Variable(index=[time]) # Total radiative forcing from tropospheric ozone (Wm⁻²).
QCH4H2O = Variable(index=[time]) # Radiative forcing for stratoshperic water vapor from methane oxidation (Wm⁻²). rf_O₃ = Variable(index=[time]) # Total radiative forcing from tropospheric ozone (Wm⁻²).


function run_timestep(p, v, d, t)
Expand Down
42 changes: 21 additions & 21 deletions src/components/rf_o3.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,23 +4,23 @@

@defcomp rf_o3 begin

CH₄_0 = Parameter() # Pre-industrial atmospheric methane concentration (ppb).
OZ00CH4 = Parameter() # Year 2000 tropospheric ozone forcing due to methane (Wm⁻²).
TROZSENS = Parameter() # Tropospheric ozone radiative efficiency factor.
OZNOX = Parameter() # Sensitivity coefficient of tropospheric ozone to nitrogen oxides emissions.
OZCO = Parameter() # Sensitivity coefficient of tropospheric ozone to carbon monoxide emissions.
OZVOC = Parameter() # Sensitivity coefficient of tropospheric ozone to Non-methane volatile organic compound emissions.
TROZSENS = Parameter() # Tropospheric ozone radiative efficiency factor.
OZCH4 = Parameter() # Sensitivity coefficient of tropospheric ozone to methane concentration.
CH₄ = Parameter(index=[time]) # Atmospheric methane concentration (ppb).
NOx_emissions = Parameter(index=[time]) # Nitrogen oxides emissions (Mt yr⁻¹).
CO_emissions = Parameter(index=[time]) # Carbon monoxide emissions (Mt yr⁻¹).
NMVOC_emissions = Parameter(index=[time]) # Non-methane volatile organic compound emissions (Mt yr⁻¹).
FOSSHIST = Parameter(index=[time]) # Historical fossil carbon dioxide emissions values used to calculate pre-2000 non-methane ozone forcing.

QOZ = Variable(index=[time]) # Tropospheric ozone radiative forcing not attributed to methane (Wm⁻²).
QCH4OZ = Variable(index=[time]) # Indirect forcing effect from methane due to methane-induced ozone production (Wm⁻²).
rf_O₃ = Variable(index=[time]) # Total radiative forcing from tropospheric ozone (Wm⁻²).
CH₄_0 = Parameter() # Pre-industrial atmospheric methane concentration (ppb).
OZ00CH4 = Parameter() # Year 2000 tropospheric ozone forcing due to methane (Wm⁻²).
TROZSENS = Parameter() # Tropospheric ozone radiative efficiency factor.
OZNOX = Parameter() # Sensitivity coefficient of tropospheric ozone to nitrogen oxides emissions.
OZCO = Parameter() # Sensitivity coefficient of tropospheric ozone to carbon monoxide emissions.
OZVOC = Parameter() # Sensitivity coefficient of tropospheric ozone to Non-methane volatile organic compound emissions.
TROZSENS = Parameter() # Tropospheric ozone radiative efficiency factor.
OZCH4 = Parameter() # Sensitivity coefficient of tropospheric ozone to methane concentration.
CH₄ = Parameter(index=[time]) # Atmospheric methane concentration (ppb).
NOx_emissions = Parameter(index=[time]) # Nitrogen oxides emissions (Mt yr⁻¹).
CO_emissions = Parameter(index=[time]) # Carbon monoxide emissions (Mt yr⁻¹).
NMVOC_emissions = Parameter(index=[time]) # Non-methane volatile organic compound emissions (Mt yr⁻¹).
FOSSHIST = Parameter(index=[time]) # Historical fossil carbon dioxide emissions values used to calculate pre-2000 non-methane ozone forcing.

QOZ = Variable(index=[time]) # Tropospheric ozone radiative forcing not attributed to methane (Wm⁻²).
QCH4OZ = Variable(index=[time]) # Indirect forcing effect from methane due to methane-induced ozone production (Wm⁻²).
rf_O₃ = Variable(index=[time]) # Total radiative forcing from tropospheric ozone (Wm⁻²).


function run_timestep(p, v, d, t)
Expand All @@ -41,16 +41,16 @@

else
# Calculate change in emissions for NOX, CO, and NMVOC relative to 2000.
DENOX = p.NOx_emissions[t] - p.NOx_emissions[TimestepValue(2000)]
DECO = p.CO_emissions[t] - p.CO_emissions[TimestepValue(2000)]
DENOX = p.NOx_emissions[t] - p.NOx_emissions[TimestepValue(2000)]
DECO = p.CO_emissions[t] - p.CO_emissions[TimestepValue(2000)]
DEVOC = p.NMVOC_emissions[t] - p.NMVOC_emissions[TimestepValue(2000)]

# Calculate tropospheric O₃ radiative forcing not attributed to CH₄.
v.QOZ[t] = (0.33 - p.OZ00CH4) + (v.QOZ[TimestepValue(2000)-1] - v.QOZ[TimestepValue(2000)-2]) + p.TROZSENS*(p.OZNOX*DENOX + p.OZCO*DECO + p.OZVOC*DEVOC)
v.QOZ[t] = (0.33 - p.OZ00CH4) + (v.QOZ[TimestepValue(2000)-1] - v.QOZ[TimestepValue(2000)-2]) + p.TROZSENS * (p.OZNOX * DENOX + p.OZCO * DECO + p.OZVOC * DEVOC)
end

# Calculate indirect CH₄ forcing due to CH₄-induced O₃ production (based on pure emissions version of CH₄ cycle).
v.QCH4OZ[t] = p.TROZSENS * p.OZCH4 * log(p.CH₄[t]/p.CH₄_0)
v.QCH4OZ[t] = p.TROZSENS * p.OZCH4 * log(p.CH₄[t] / p.CH₄_0)

# Calculate total tropospheric O₃ radiative forcing (from direct O₃ and indirect CH₄ effects).
v.rf_O₃[t] = v.QOZ[t] + v.QCH4OZ[t]
Expand Down
Loading