From a10f6fc15580770522f5e54b4116c9a1b192bcfd Mon Sep 17 00:00:00 2001 From: davidanthoff Date: Fri, 14 Oct 2022 13:13:29 +0000 Subject: [PATCH] Format files using DocumentFormat --- src/MimiMAGICC.jl | 24 +++++------ src/components/ch4_cycle.jl | 82 ++++++++++++++++++------------------- src/components/rf_ch4.jl | 12 +++--- src/components/rf_ch4h2o.jl | 8 ++-- src/components/rf_o3.jl | 42 +++++++++---------- test/runtests.jl | 10 ++--- 6 files changed, 89 insertions(+), 89 deletions(-) diff --git a/src/MimiMAGICC.jl b/src/MimiMAGICC.jl index bfed98c..c97fc1c 100644 --- a/src/MimiMAGICC.jl +++ b/src/MimiMAGICC.jl @@ -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. @@ -64,8 +64,8 @@ function get_model(;rcp_scenario::String="RCP85", start_year::Int=1765, end_year set_param!(m, :ch4_cycle, :AVOC, -0.000315) set_param!(m, :ch4_cycle, :TAUINIT, 7.73) set_param!(m, :ch4_cycle, :SCH4, -0.32) - set_param!(m, :ch4_cycle, :TAUSOIL, 160.) - set_param!(m, :ch4_cycle, :TAUSTRAT, 120.) + set_param!(m, :ch4_cycle, :TAUSOIL, 160.0) + set_param!(m, :ch4_cycle, :TAUSTRAT, 120.0) set_param!(m, :ch4_cycle, :GAM, -1.0) set_param!(m, :ch4_cycle, :CH4_natural, 266.5) set_param!(m, :ch4_cycle, :fffrac, 0.18) @@ -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 diff --git a/src/components/ch4_cycle.jl b/src/components/ch4_cycle.jl index 0447fe4..e4e3b6d 100644 --- a/src/components/ch4_cycle.jl +++ b/src/components/ch4_cycle.jl @@ -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 @@ -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. @@ -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 @@ -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. diff --git a/src/components/rf_ch4.jl b/src/components/rf_ch4.jl index d8571d5..ae1228e 100644 --- a/src/components/rf_ch4.jl +++ b/src/components/rf_ch4.jl @@ -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₄ diff --git a/src/components/rf_ch4h2o.jl b/src/components/rf_ch4h2o.jl index 099fb05..43f7d94 100644 --- a/src/components/rf_ch4h2o.jl +++ b/src/components/rf_ch4h2o.jl @@ -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) diff --git a/src/components/rf_o3.jl b/src/components/rf_o3.jl index 427fb99..b9bf8e8 100644 --- a/src/components/rf_o3.jl +++ b/src/components/rf_o3.jl @@ -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) @@ -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] diff --git a/test/runtests.jl b/test/runtests.jl index c637b79..24055b6 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -6,14 +6,14 @@ using MimiMAGICC @testset "MAGICC" begin -#------------------------------------------------------------------------------ -# 1. Carry out test to check that the model runs. -#------------------------------------------------------------------------------ + #------------------------------------------------------------------------------ + # 1. Carry out test to check that the model runs. + #------------------------------------------------------------------------------ @testset "MAGICC-model" begin - m = MimiMAGICC.get_model() - run(m) + m = MimiMAGICC.get_model() + run(m) end # MAGICC-CH4 model run test. end # All MAGICC-CH4 tests.