diff --git a/src/MimiDICE2010.jl b/src/MimiDICE2010.jl index 311f2b5..3a09927 100644 --- a/src/MimiDICE2010.jl +++ b/src/MimiDICE2010.jl @@ -1,7 +1,7 @@ module MimiDICE2010 using Mimi -using XLSX:readxlsx +using XLSX: readxlsx include("parameters.jl") include("marginaldamage.jl") @@ -48,7 +48,7 @@ function get_model(params=nothing) # Socioeconomics connect_param!(m, :grosseconomy, :I, :neteconomy, :I) connect_param!(m, :emissions, :YGROSS, :grosseconomy, :YGROSS) - + # Climate connect_param!(m, :co2cycle, :E, :emissions, :E) connect_param!(m, :radiativeforcing, :MAT, :co2cycle, :MAT) @@ -76,4 +76,4 @@ end construct_dice = get_model # still export the old version of the function name -end # module \ No newline at end of file +end # module diff --git a/src/components/climatedynamics_component.jl b/src/components/climatedynamics_component.jl index 1195fbd..0d826d9 100644 --- a/src/components/climatedynamics_component.jl +++ b/src/components/climatedynamics_component.jl @@ -1,12 +1,12 @@ @defcomp climatedynamics begin - TATM = Variable(index=[time]) # Increase in temperature of atmosphere (degrees C from 1900) - TOCEAN = Variable(index=[time]) # Increase in temperature of lower oceans (degrees C from 1900) + TATM = Variable(index=[time]) # Increase in temperature of atmosphere (degrees C from 1900) + TOCEAN = Variable(index=[time]) # Increase in temperature of lower oceans (degrees C from 1900) - FORC = Parameter(index=[time]) # Increase in radiative forcing (watts per m2 from 1900) - fco22x = Parameter() # Forcings of equilibrium CO2 doubling (Wm-2) - t2xco2 = Parameter() # Equilibrium temp impact (oC per doubling CO2) - tatm0 = Parameter() # Initial atmospheric temp change (C from 1900) in 2005 - tatm1 = Parameter() # Initial atmospheric temp change (C from 1900) in 2015 + FORC = Parameter(index=[time]) # Increase in radiative forcing (watts per m2 from 1900) + fco22x = Parameter() # Forcings of equilibrium CO2 doubling (Wm-2) + t2xco2 = Parameter() # Equilibrium temp impact (oC per doubling CO2) + tatm0 = Parameter() # Initial atmospheric temp change (C from 1900) in 2005 + tatm1 = Parameter() # Initial atmospheric temp change (C from 1900) in 2015 tocean0 = Parameter() # Initial lower stratum temp change (C from 1900) # Transient TSC Correction ("Speed of Adjustment Parameter") @@ -25,10 +25,10 @@ if t.t == 2 v.TATM[t] = p.tatm1 else - v.TATM[t] = v.TATM[t - 1] + p.c1 * ((p.FORC[t] - (p.fco22x / p.t2xco2) * v.TATM[t - 1]) - (p.c3 * (v.TATM[t - 1] - v.TOCEAN[t - 1]))) + v.TATM[t] = v.TATM[t-1] + p.c1 * ((p.FORC[t] - (p.fco22x / p.t2xco2) * v.TATM[t-1]) - (p.c3 * (v.TATM[t-1] - v.TOCEAN[t-1]))) end - v.TOCEAN[t] = v.TOCEAN[t - 1] + p.c4 * (v.TATM[t - 1] - v.TOCEAN[t - 1]) + v.TOCEAN[t] = v.TOCEAN[t-1] + p.c4 * (v.TATM[t-1] - v.TOCEAN[t-1]) end end end diff --git a/src/components/co2cycle_component.jl b/src/components/co2cycle_component.jl index ce43d00..7425eea 100644 --- a/src/components/co2cycle_component.jl +++ b/src/components/co2cycle_component.jl @@ -1,23 +1,23 @@ @defcomp co2cycle begin - MAT = Variable(index=[time]) # Carbon concentration increase in atmosphere (GtC from 1750) - MAT_final = Variable() # MAT calculation one timestep further than the model's index - ML = Variable(index=[time]) # Carbon concentration increase in lower oceans (GtC from 1750) - MU = Variable(index=[time]) # Carbon concentration increase in shallow oceans (GtC from 1750) + MAT = Variable(index=[time]) # Carbon concentration increase in atmosphere (GtC from 1750) + MAT_final = Variable() # MAT calculation one timestep further than the model's index + ML = Variable(index=[time]) # Carbon concentration increase in lower oceans (GtC from 1750) + MU = Variable(index=[time]) # Carbon concentration increase in shallow oceans (GtC from 1750) - E = Parameter(index=[time]) # Total CO2 emissions (GtC per year) - mat0 = Parameter() # Initial Concentration in atmosphere 2000 (GtC) - mat1 = Parameter() # Concentration 2010 (GtC) - ml0 = Parameter() # Initial Concentration in lower strata (GtC) - mu0 = Parameter() # Initial Concentration in upper strata (GtC) + E = Parameter(index=[time]) # Total CO2 emissions (GtC per year) + mat0 = Parameter() # Initial Concentration in atmosphere 2000 (GtC) + mat1 = Parameter() # Concentration 2010 (GtC) + ml0 = Parameter() # Initial Concentration in lower strata (GtC) + mu0 = Parameter() # Initial Concentration in upper strata (GtC) # Parameters for long-run consistency of carbon cycle - b11 = Parameter() # Carbon cycle transition matrix atmosphere to atmosphere - b12 = Parameter() # Carbon cycle transition matrix atmosphere to shallow ocean - b21 = Parameter() # Carbon cycle transition matrix biosphere/shallow oceans to atmosphere - b22 = Parameter() # Carbon cycle transition matrix shallow ocean to shallow oceans - b23 = Parameter() # Carbon cycle transition matrix shallow to deep ocean - b32 = Parameter() # Carbon cycle transition matrix deep ocean to shallow ocean - b33 = Parameter() # Carbon cycle transition matrix deep ocean to deep oceans + b11 = Parameter() # Carbon cycle transition matrix atmosphere to atmosphere + b12 = Parameter() # Carbon cycle transition matrix atmosphere to shallow ocean + b21 = Parameter() # Carbon cycle transition matrix biosphere/shallow oceans to atmosphere + b22 = Parameter() # Carbon cycle transition matrix shallow ocean to shallow oceans + b23 = Parameter() # Carbon cycle transition matrix shallow to deep ocean + b32 = Parameter() # Carbon cycle transition matrix deep ocean to shallow ocean + b33 = Parameter() # Carbon cycle transition matrix deep ocean to deep oceans function run_timestep(p, v, d, t) # Define functions for MU, ML, and MAT @@ -29,15 +29,15 @@ v.MAT[TimestepIndex(1)] = p.mat0 v.MAT[TimestepIndex(2)] = p.mat1 - else + else - v.MU[t] = v.MAT[t - 1] * p.b12 + v.MU[t - 1] * p.b22 + v.ML[t - 1] * p.b32 + v.MU[t] = v.MAT[t-1] * p.b12 + v.MU[t-1] * p.b22 + v.ML[t-1] * p.b32 - v.ML[t] = v.ML[t - 1] * p.b33 + v.MU[t - 1] * p.b23 + v.ML[t] = v.ML[t-1] * p.b33 + v.MU[t-1] * p.b23 + + if !is_last(t) + v.MAT[t+1] = v.MAT[t] * p.b11 + v.MU[t] * p.b21 + p.E[t] * 10 - if ! is_last(t) - v.MAT[t + 1] = v.MAT[t] * p.b11 + v.MU[t] * p.b21 + p.E[t] * 10 - else # last timestep v.MAT_final = v.MAT[t] * p.b11 + v.MU[t] * p.b21 + p.E[t] * 10 end diff --git a/src/components/damages_component.jl b/src/components/damages_component.jl index a34059a..1e370fb 100644 --- a/src/components/damages_component.jl +++ b/src/components/damages_component.jl @@ -1,20 +1,20 @@ @defcomp damages begin DAMFRAC = Variable(index=[time]) # Damages (fraction of gross output) - TATM = Parameter(index=[time]) # Increase temperature of atmosphere (degrees C from 1900) - YGROSS = Parameter(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) - a1 = Parameter() # Damage intercept - a2 = Parameter() # Damage quadratic term - a3 = Parameter() # Damage exponent - - TotSLR = Parameter(index=[time]) # Path of total SLR - b1 = Parameter() # Coefficient on SLR - b2 = Parameter() # Coefficient on quadratic SLR term - b3 = Parameter() # SLR exponent + TATM = Parameter(index=[time]) # Increase temperature of atmosphere (degrees C from 1900) + YGROSS = Parameter(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) + a1 = Parameter() # Damage intercept + a2 = Parameter() # Damage quadratic term + a3 = Parameter() # Damage exponent + + TotSLR = Parameter(index=[time]) # Path of total SLR + b1 = Parameter() # Coefficient on SLR + b2 = Parameter() # Coefficient on quadratic SLR term + b3 = Parameter() # SLR exponent function run_timestep(p, v, d, t) # Define function for DAMFRAC - v.DAMFRAC[t] = p.a1 * p.TATM[t] + p.a2 * p.TATM[t]^p.a3 + p.b1 * p.TotSLR[t] + p.b2 * p.TotSLR[t]^p.b3 + v.DAMFRAC[t] = p.a1 * p.TATM[t] + p.a2 * p.TATM[t]^p.a3 + p.b1 * p.TotSLR[t] + p.b2 * p.TotSLR[t]^p.b3 end end diff --git a/src/components/emissions_component.jl b/src/components/emissions_component.jl index 844e70c..a404277 100644 --- a/src/components/emissions_component.jl +++ b/src/components/emissions_component.jl @@ -1,25 +1,25 @@ @defcomp emissions begin - CCA = Variable(index=[time]) # Cumulative indiustrial emissions - E = Variable(index=[time]) # Total CO2 emissions (GtC per year) - EIND = Variable(index=[time]) # Industrial emissions (GtC per year) + CCA = Variable(index=[time]) # Cumulative indiustrial emissions + E = Variable(index=[time]) # Total CO2 emissions (GtC per year) + EIND = Variable(index=[time]) # Industrial emissions (GtC per year) - etree = Parameter(index=[time]) # Emissions from deforestation - MIU = Parameter(index=[time]) # Emission control rate GHGs - sigma = Parameter(index=[time]) # CO2-equivalent-emissions output ratio - YGROSS = Parameter(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) + etree = Parameter(index=[time]) # Emissions from deforestation + MIU = Parameter(index=[time]) # Emission control rate GHGs + sigma = Parameter(index=[time]) # CO2-equivalent-emissions output ratio + YGROSS = Parameter(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) function run_timestep(p, v, d, t) # Define function for EIND v.EIND[t] = p.sigma[t] * p.YGROSS[t] * (1 - p.MIU[t]) - + # Define function for E v.E[t] = v.EIND[t] + p.etree[t] - + # Define function for CCA if is_first(t) - v.CCA[t] = v.EIND[t] * 10 + v.CCA[t] = v.EIND[t] * 10 else - v.CCA[t] = v.CCA[t - 1] + v.EIND[t] * 10 + v.CCA[t] = v.CCA[t-1] + v.EIND[t] * 10 end end end diff --git a/src/components/grosseconomy_component.jl b/src/components/grosseconomy_component.jl index 310ed9d..8fefc9d 100644 --- a/src/components/grosseconomy_component.jl +++ b/src/components/grosseconomy_component.jl @@ -1,20 +1,20 @@ @defcomp grosseconomy begin - K = Variable(index=[time]) # Capital stock (trillions 2005 US dollars) - YGROSS = Variable(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) + K = Variable(index=[time]) # Capital stock (trillions 2005 US dollars) + YGROSS = Variable(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) - al = Parameter(index=[time]) # Level of total factor productivity - I = Parameter(index=[time]) # Investment (trillions 2005 USD per year) - l = Parameter(index=[time]) # Level of population and labor - dk = Parameter() # Depreciation rate on capital (per year) - gama = Parameter() # Capital elasticity in production function - k0 = Parameter() # Initial capital value (trill 2005 USD) + al = Parameter(index=[time]) # Level of total factor productivity + I = Parameter(index=[time]) # Investment (trillions 2005 USD per year) + l = Parameter(index=[time]) # Level of population and labor + dk = Parameter() # Depreciation rate on capital (per year) + gama = Parameter() # Capital elasticity in production function + k0 = Parameter() # Initial capital value (trill 2005 USD) function run_timestep(p, v, d, t) # Define function for K if is_first(t) v.K[t] = p.k0 else - v.K[t] = v.K[t - 1] * (1 - p.dk)^10 + 10 * p.I[t - 1] + v.K[t] = v.K[t-1] * (1 - p.dk)^10 + 10 * p.I[t-1] end # Define function for YGROSS diff --git a/src/components/neteconomy_component.jl b/src/components/neteconomy_component.jl index 218a3b5..7c9629c 100644 --- a/src/components/neteconomy_component.jl +++ b/src/components/neteconomy_component.jl @@ -1,22 +1,22 @@ @defcomp neteconomy begin - ABATECOST = Variable(index=[time]) # Cost of emissions reductions (trillions 2005 USD per year) - C = Variable(index=[time]) # Consumption (trillions 2005 US dollars per year) - CPC = Variable(index=[time]) # Per capita consumption (thousands 2005 USD per year) - CPRICE = Variable(index=[time]) # Carbon price (2005$ per ton of CO2) - I = Variable(index=[time]) # Investment (trillions 2005 USD per year) - Y = Variable(index=[time]) # Gross world product net of abatement and damages (trillions 2005 USD per year) - YNET = Variable(index=[time]) # Output net of damages equation (trillions 2005 USD per year) - - cost1 = Parameter(index=[time]) # Abatement cost function coefficient - DAMFRAC = Parameter(index=[time]) # Damages (fraction of gross output) - l = Parameter(index=[time]) # Level of population and labor - MIU = Parameter(index=[time]) # Emission control rate GHGs - partfract = Parameter(index=[time]) # Fraction of emissions in control regime - pbacktime = Parameter(index=[time]) # Backstop price - S = Parameter(index=[time]) # Gross savings rate as fraction of gross world product - YGROSS = Parameter(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) - expcost2 = Parameter() # Exponent of control cost function + ABATECOST = Variable(index=[time]) # Cost of emissions reductions (trillions 2005 USD per year) + C = Variable(index=[time]) # Consumption (trillions 2005 US dollars per year) + CPC = Variable(index=[time]) # Per capita consumption (thousands 2005 USD per year) + CPRICE = Variable(index=[time]) # Carbon price (2005$ per ton of CO2) + I = Variable(index=[time]) # Investment (trillions 2005 USD per year) + Y = Variable(index=[time]) # Gross world product net of abatement and damages (trillions 2005 USD per year) + YNET = Variable(index=[time]) # Output net of damages equation (trillions 2005 USD per year) + + cost1 = Parameter(index=[time]) # Abatement cost function coefficient + DAMFRAC = Parameter(index=[time]) # Damages (fraction of gross output) + l = Parameter(index=[time]) # Level of population and labor + MIU = Parameter(index=[time]) # Emission control rate GHGs + partfract = Parameter(index=[time]) # Fraction of emissions in control regime + pbacktime = Parameter(index=[time]) # Backstop price + S = Parameter(index=[time]) # Gross savings rate as fraction of gross world product + YGROSS = Parameter(index=[time]) # Gross world product GROSS of abatement and damages (trillions 2005 USD per year) + expcost2 = Parameter() # Exponent of control cost function function run_timestep(p, v, d, t) # Define function for YNET @@ -37,7 +37,7 @@ elseif t < TimestepIndex(60) v.C[t] = v.Y[t] - v.I[t] elseif t == TimestepIndex(60) - v.C[t] = v.C[t - 1] + v.C[t] = v.C[t-1] end # Define function for CPC @@ -45,7 +45,7 @@ # Define function for CPRICE if t == TimestepIndex(26) - v.CPRICE[t] = v.CPRICE[t - 1] + v.CPRICE[t] = v.CPRICE[t-1] else v.CPRICE[t] = p.pbacktime[t] * 1000 * p.MIU[t]^(p.expcost2 - 1) end diff --git a/src/components/radiativeforcing_component.jl b/src/components/radiativeforcing_component.jl index 2dbaf25..4d9088c 100644 --- a/src/components/radiativeforcing_component.jl +++ b/src/components/radiativeforcing_component.jl @@ -1,15 +1,15 @@ @defcomp radiativeforcing begin - FORC = Variable(index=[time]) # Increase in radiative forcing (watts per m2 from 1900) + FORC = Variable(index=[time]) # Increase in radiative forcing (watts per m2 from 1900) - forcoth = Parameter(index=[time]) # Exogenous forcing for other greenhouse gases - MAT = Parameter(index=[time]) # Carbon concentration increase in atmosphere (GtC from 1750) + forcoth = Parameter(index=[time]) # Exogenous forcing for other greenhouse gases + MAT = Parameter(index=[time]) # Carbon concentration increase in atmosphere (GtC from 1750) MAT_final = Parameter() # MAT calculation one timestep further than the model's index - fco22x = Parameter() # Forcings of equilibrium CO2 doubling (Wm-2) + fco22x = Parameter() # Forcings of equilibrium CO2 doubling (Wm-2) function run_timestep(p, v, d, t) # Define function for FORC - if ! is_last(t) - v.FORC[t] = p.fco22x * (log((((p.MAT[t] + p.MAT[t + 1]) / 2) + 0.000001) / 596.4) / log(2)) + p.forcoth[t] + if !is_last(t) + v.FORC[t] = p.fco22x * (log((((p.MAT[t] + p.MAT[t+1]) / 2) + 0.000001) / 596.4) / log(2)) + p.forcoth[t] else # final timestep # need to use MAT_final, calculated one step further v.FORC[t] = p.fco22x * (log((((p.MAT[t] + p.MAT_final) / 2) + 0.000001) / 596.4) / log(2)) + p.forcoth[t] diff --git a/src/components/sealevelrise_component.jl b/src/components/sealevelrise_component.jl index 78b7d1b..5d2a5d7 100644 --- a/src/components/sealevelrise_component.jl +++ b/src/components/sealevelrise_component.jl @@ -1,48 +1,48 @@ @defcomp sealevelrise begin - ThermSLR = Variable(index=[time]) # Path of SLR from thermal expansion - GSICSLR = Variable(index=[time]) # Path of SLR from G&SIC - GISSLR = Variable(index=[time]) # Path of SLR from GIS - AISSLR = Variable(index=[time]) # Path of SLR from AIS - TotSLR = Variable(index=[time]) # Path of total SLR + ThermSLR = Variable(index=[time]) # Path of SLR from thermal expansion + GSICSLR = Variable(index=[time]) # Path of SLR from G&SIC + GISSLR = Variable(index=[time]) # Path of SLR from GIS + AISSLR = Variable(index=[time]) # Path of SLR from AIS + TotSLR = Variable(index=[time]) # Path of total SLR - TATM = Parameter(index=[time]) # Path of atmospheric temperature anomalies (from climate dynamics?) + TATM = Parameter(index=[time]) # Path of atmospheric temperature anomalies (from climate dynamics?) - therm0 = Parameter() # Initial SLR from thermal expansion (meters above 2000) - gsic0 = Parameter() # Initial SLR from G&SIC - gis0 = Parameter() # Initial SLR from GIS - ais0 = Parameter() # Initial SLR from AIS + therm0 = Parameter() # Initial SLR from thermal expansion (meters above 2000) + gsic0 = Parameter() # Initial SLR from G&SIC + gis0 = Parameter() # Initial SLR from GIS + ais0 = Parameter() # Initial SLR from AIS - therm_asym = Parameter() # Asymptotic rise from thermal expansion - gsic_asym = Parameter() # Asymptotic rise from G&SIC - gis_asym = Parameter() # Asymptotic rise from GIS - ais_asym = Parameter() # Asymptotic rise from AIS + therm_asym = Parameter() # Asymptotic rise from thermal expansion + gsic_asym = Parameter() # Asymptotic rise from G&SIC + gis_asym = Parameter() # Asymptotic rise from GIS + ais_asym = Parameter() # Asymptotic rise from AIS - thermrate = Parameter() # Rate of thermal expansion - gsicrate = Parameter() # Rate of G&SIC - gisrate = Parameter() # Rate of GIS - aisrate = Parameter() # Rate of AIS + thermrate = Parameter() # Rate of thermal expansion + gsicrate = Parameter() # Rate of G&SIC + gisrate = Parameter() # Rate of GIS + aisrate = Parameter() # Rate of AIS slrthreshold = Parameter() # Temperature threshold for AIS function run_timestep(p, v, d, t) if is_first(t) v.ThermSLR[t] = p.therm0 - v.GSICSLR[t] = p.gsic0 - v.GISSLR[t] = p.gis0 - v.AISSLR[t] = p.ais0 + v.GSICSLR[t] = p.gsic0 + v.GISSLR[t] = p.gis0 + v.AISSLR[t] = p.ais0 else - v.ThermSLR[t] = v.ThermSLR[t - 1] + p.thermrate * p.TATM[t] - v.GSICSLR[t] = v.GSICSLR[t - 1] + p.gsicrate * (p.gsic_asym - v.GSICSLR[t - 1]) * p.TATM[t] - v.GISSLR[t] = v.GISSLR[t - 1] + p.gisrate * (p.gis_asym - v.GISSLR[t - 1]) * p.TATM[t] - v.AISSLR[t] = 0 + v.ThermSLR[t] = v.ThermSLR[t-1] + p.thermrate * p.TATM[t] + v.GSICSLR[t] = v.GSICSLR[t-1] + p.gsicrate * (p.gsic_asym - v.GSICSLR[t-1]) * p.TATM[t] + v.GISSLR[t] = v.GISSLR[t-1] + p.gisrate * (p.gis_asym - v.GISSLR[t-1]) * p.TATM[t] + v.AISSLR[t] = 0 if p.TATM[t] > p.slrthreshold - v.AISSLR[t] = v.AISSLR[t - 1] + p.aisrate * (p.ais_asym - v.AISSLR[t - 1]) * p.TATM[t] - end + v.AISSLR[t] = v.AISSLR[t-1] + p.aisrate * (p.ais_asym - v.AISSLR[t-1]) * p.TATM[t] + end end - + v.TotSLR[t] = v.ThermSLR[t] + v.GSICSLR[t] + v.GISSLR[t] + v.AISSLR[t] end end diff --git a/src/components/welfare_component.jl b/src/components/welfare_component.jl index a2cf01c..1af9931 100644 --- a/src/components/welfare_component.jl +++ b/src/components/welfare_component.jl @@ -1,14 +1,14 @@ @defcomp welfare begin - CEMUTOTPER = Variable(index=[time]) # Period utility - PERIODU = Variable(index=[time]) # One period utility function - UTILITY = Variable() # Welfare Function + CEMUTOTPER = Variable(index=[time]) # Period utility + PERIODU = Variable(index=[time]) # One period utility function + UTILITY = Variable() # Welfare Function - CPC = Parameter(index=[time]) # Per capita consumption (thousands 2005 USD per year) - l = Parameter(index=[time]) # Level of population and labor - rr = Parameter(index=[time]) # Average utility social discount rate - elasmu = Parameter() # Elasticity of marginal utility of consumption - scale1 = Parameter() # Multiplicative scaling coefficient - scale2 = Parameter() # Additive scaling coefficient + CPC = Parameter(index=[time]) # Per capita consumption (thousands 2005 USD per year) + l = Parameter(index=[time]) # Level of population and labor + rr = Parameter(index=[time]) # Average utility social discount rate + elasmu = Parameter() # Elasticity of marginal utility of consumption + scale1 = Parameter() # Multiplicative scaling coefficient + scale2 = Parameter() # Additive scaling coefficient function run_timestep(p, v, d, t) # Define function for PERIODU diff --git a/src/marginaldamage.jl b/src/marginaldamage.jl index dc4abfa..209e206 100644 --- a/src/marginaldamage.jl +++ b/src/marginaldamage.jl @@ -32,8 +32,8 @@ function compute_scc_mm(m::Model=get_model(); year::Union{Int,Nothing}=nothing, mm = get_marginal_model(m; year=year, pulse_size=pulse_size) scc = _compute_scc(mm; year=year, last_year=last_year, prtp=prtp, eta=eta) - - return (scc = scc, mm = mm) + + return (scc=scc, mm=mm) end # helper function for computing SCC from a MarginalModel, not to be exported or advertised to users @@ -72,7 +72,7 @@ end """ Adds a marginal emission component to `m` which adds `pulse_size` of additional C emissions over ten years starting in the specified `year`. """ -function add_marginal_emissions!(m::Model, year::Int, pulse_size::Float64) +function add_marginal_emissions!(m::Model, year::Int, pulse_size::Float64) add_comp!(m, Mimi.adder, :marginalemission, before=:co2cycle) time = Mimi.dimension(m, :time) @@ -87,7 +87,7 @@ end # Old available function -function getmarginal_dice_models(;emissionyear=2015) +function getmarginal_dice_models(; emissionyear=2015) DICE = get_model() run(DICE) diff --git a/src/parameters.jl b/src/parameters.jl index 7d9eb94..7068154 100644 --- a/src/parameters.jl +++ b/src/parameters.jl @@ -36,27 +36,27 @@ function dice2010_excel_parameters(filename=datafile; nsteps=nothing) f = readxlsx(filename) # Preferences - p[:elasmu] = read_params(f, "B19") # Elasticity of marginal utility of consumption - p[:rr] = read_params(f, "B18:BI18", nsteps) # Social time preference factor + p[:elasmu] = read_params(f, "B19") # Elasticity of marginal utility of consumption + p[:rr] = read_params(f, "B18:BI18", nsteps) # Social time preference factor # Population and technology - p[:gama] = read_params(f, "B9") # Capital elasticity in production function - p[:l] = read_params(f, "B27:BI27", nsteps) # Population (millions) - p[:dk] = read_params(f, "B10") # Depreciation rate on capital (per year) - p[:k0] = read_params(f, "B13") # Initial capital value (trill 2005 USD) - p[:al] = read_params(f, "B21:BI21", nsteps) # Level of total factor productivity + p[:gama] = read_params(f, "B9") # Capital elasticity in production function + p[:l] = read_params(f, "B27:BI27", nsteps) # Population (millions) + p[:dk] = read_params(f, "B10") # Depreciation rate on capital (per year) + p[:k0] = read_params(f, "B13") # Initial capital value (trill 2005 USD) + p[:al] = read_params(f, "B21:BI21", nsteps) # Level of total factor productivity # Emissions parameters - p[:sigma] = read_params(f, "B46:BI46", nsteps) - p[:etree] = read_params(f, "B52:BI52", nsteps) + p[:sigma] = read_params(f, "B46:BI46", nsteps) + p[:etree] = read_params(f, "B52:BI52", nsteps) p[:MIU] = read_params(f, "B133:BI133", nsteps) # emissions control rate p[:S] = read_params(f, "B132:BI132", nsteps) / 100 # savings rate # Carbon cycle - p[:mat0] = read_params(f, "B57") # Initial Concentration in atmosphere 2000 (GtC) - p[:mat1] = read_params(f, "B58") # Initial Concentration in atmosphere 2010 (GtC) - p[:mu0] = read_params(f, "B59") # Initial Concentration in biosphere/shallow oceans (GtC) - p[:ml0] = read_params(f, "B60") # Initial Concentration in deep oceans (GtC) + p[:mat0] = read_params(f, "B57") # Initial Concentration in atmosphere 2000 (GtC) + p[:mat1] = read_params(f, "B58") # Initial Concentration in atmosphere 2010 (GtC) + p[:mu0] = read_params(f, "B59") # Initial Concentration in biosphere/shallow oceans (GtC) + p[:ml0] = read_params(f, "B60") # Initial Concentration in deep oceans (GtC) # Flow parameters p[:b12] = read_params(f, "B64") / 100 # Carbon cycle transition matrix atmosphere to shallow ocean @@ -70,16 +70,16 @@ function dice2010_excel_parameters(filename=datafile; nsteps=nothing) p[:b33] = read_params(f, "B68") / 100 # Carbon cycle transition matrix deep ocean to deep oceans # Climate model parameters - p[:t2xco2] = read_params(f, "B76") # Equilibrium temp impact (oC per doubling CO2) - p[:tatm0] = read_params(f, "B73") # Initial lower stratum temp change (C from 1900) - p[:tatm1] = read_params(f, "C73") # Initial atmospheric temp change 2015 (C from 1900) + p[:t2xco2] = read_params(f, "B76") # Equilibrium temp impact (oC per doubling CO2) + p[:tatm0] = read_params(f, "B73") # Initial lower stratum temp change (C from 1900) + p[:tatm1] = read_params(f, "C73") # Initial atmospheric temp change 2015 (C from 1900) p[:tocean0] = read_params(f, "B74") # Initial lower stratum temp change (C from 1900) # Transient TSC Correction ("Speed of Adjustment Parameter") - p[:c1] = read_params(f, "B75") # Climate equation coefficient for upper level - p[:c3] = read_params(f, "B78") # Transfer coefficient upper to lower stratum - p[:c4] = read_params(f, "B79") # Transfer coefficient for lower level - p[:fco22x] = read_params(f, "B77") # Forcings of equilibrium CO2 doubling (Wm-2) + p[:c1] = read_params(f, "B75") # Climate equation coefficient for upper level + p[:c3] = read_params(f, "B78") # Transfer coefficient upper to lower stratum + p[:c4] = read_params(f, "B79") # Transfer coefficient for lower level + p[:fco22x] = read_params(f, "B77") # Forcings of equilibrium CO2 doubling (Wm-2) # Climate damage parameters p[:a1] = read_params(f, "B33") # Damage coefficient @@ -87,10 +87,10 @@ function dice2010_excel_parameters(filename=datafile; nsteps=nothing) p[:a3] = read_params(f, "B35") # Damage exponent # Abatement cost - p[:expcost2] = read_params(f, "B44") # Exponent of control cost function - p[:pbacktime] = read_params(f, "B42:BI42", nsteps) # backstop price + p[:expcost2] = read_params(f, "B44") # Exponent of control cost function + p[:pbacktime] = read_params(f, "B42:BI42", nsteps) # backstop price # Adjusted cost for backstop (or: "Abatement cost function coefficient") - p[:cost1] = read_params(f, "B37:BI37", nsteps) + p[:cost1] = read_params(f, "B37:BI37", nsteps) # Scaling and inessential parameters p[:scale1] = read_params(f, "B88") # Multiplicative scaling coefficient @@ -105,24 +105,24 @@ function dice2010_excel_parameters(filename=datafile; nsteps=nothing) # SLR Parameters - p[:b1] = read_params(f, "B51", sheet="Parameters") - p[:b2] = read_params(f, "B52", sheet="Parameters") - p[:b3] = read_params(f, "B53", sheet="Parameters") + p[:b1] = read_params(f, "B51", sheet="Parameters") + p[:b2] = read_params(f, "B52", sheet="Parameters") + p[:b3] = read_params(f, "B53", sheet="Parameters") - p[:therm0] = read_params(f, "B178") # meters above 2000 - p[:gsic0] = read_params(f, "B179") - p[:gis0] = read_params(f, "B180") - p[:ais0] = read_params(f, "B181") + p[:therm0] = read_params(f, "B178") # meters above 2000 + p[:gsic0] = read_params(f, "B179") + p[:gis0] = read_params(f, "B180") + p[:ais0] = read_params(f, "B181") - p[:therm_asym] = read_params(f, "B173") - p[:gsic_asym] = read_params(f, "B174") - p[:gis_asym] = read_params(f, "B175") - p[:ais_asym] = read_params(f, "B176") + p[:therm_asym] = read_params(f, "B173") + p[:gsic_asym] = read_params(f, "B174") + p[:gis_asym] = read_params(f, "B175") + p[:ais_asym] = read_params(f, "B176") - p[:thermrate] = read_params(f, "C173") - p[:gsicrate] = read_params(f, "C174") - p[:gisrate] = read_params(f, "C175") - p[:aisrate] = read_params(f, "C176") + p[:thermrate] = read_params(f, "C173") + p[:gsicrate] = read_params(f, "C174") + p[:gisrate] = read_params(f, "C175") + p[:aisrate] = read_params(f, "C176") p[:slrthreshold] = read_params(f, "D176") diff --git a/test/runtests.jl b/test/runtests.jl index f19dfb4..c541900 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,7 +1,7 @@ using Test using Mimi using MimiDICE2010 -using XLSX:readxlsx +using XLSX: readxlsx using DataFrames using CSVFiles @@ -42,10 +42,10 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # Climate dynamics tests - TATM = m[:climatedynamics, :TATM] - TOCEAN = m[:climatedynamics, :TOCEAN] + TATM = m[:climatedynamics, :TATM] + TOCEAN = m[:climatedynamics, :TOCEAN] - True_TATM = read_params(f, "B121:BI121", T) + True_TATM = read_params(f, "B121:BI121", T) True_TOCEAN = read_params(f, "B123:BI123", T) @test maximum(abs, TATM .- True_TATM) ≈ 0. atol = Precision @@ -53,15 +53,15 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # CO2 Cycle tests - MAT = m[:co2cycle, :MAT] - MAT_final = m[:co2cycle, :MAT_final] - ML = m[:co2cycle, :ML] - MU = m[:co2cycle, :MU] + MAT = m[:co2cycle, :MAT] + MAT_final = m[:co2cycle, :MAT_final] + ML = m[:co2cycle, :ML] + MU = m[:co2cycle, :MU] - True_MAT = read_params(f, "B112:BI112", T) - True_MAT_final = read_params(f, "BJ112") - True_ML = read_params(f, "B115:BI115", T) - True_MU = read_params(f, "B114:BI114", T) + True_MAT = read_params(f, "B112:BI112", T) + True_MAT_final = read_params(f, "BJ112") + True_ML = read_params(f, "B115:BI115", T) + True_MU = read_params(f, "B114:BI114", T) @test maximum(abs, MAT .- True_MAT) ≈ 0. atol = Precision @test abs(MAT_final - True_MAT_final) ≈ 0. atol = Precision @@ -72,18 +72,18 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # Damages test DAMFRAC = m[:damages, :DAMFRAC] - True_DAMFRAC = read_params(f, "B93:BI93", T) + True_DAMFRAC = read_params(f, "B93:BI93", T) @test maximum(abs, DAMFRAC .- True_DAMFRAC) ≈ 0. atol = Precision # Emissions tests - CCA = m[:emissions, :CCA] - E = m[:emissions, :E] - EIND = m[:emissions, :EIND] + CCA = m[:emissions, :CCA] + E = m[:emissions, :E] + EIND = m[:emissions, :EIND] - True_CCA = read_params(f, "B117:BI117", T) - True_E = read_params(f, "B109:BI109", T) - True_EIND = read_params(f, "B110:BI110", T) + True_CCA = read_params(f, "B117:BI117", T) + True_E = read_params(f, "B109:BI109", T) + True_EIND = read_params(f, "B110:BI110", T) @test maximum(abs, CCA .- True_CCA) ≈ 0. atol = Precision @test maximum(abs, E .- True_E) ≈ 0. atol = Precision @@ -91,10 +91,10 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # Gross Economy tests - K = m[:grosseconomy, :K] - YGROSS = m[:grosseconomy, :YGROSS] + K = m[:grosseconomy, :K] + YGROSS = m[:grosseconomy, :YGROSS] - True_K = read_params(f, "B102:BI102", T) + True_K = read_params(f, "B102:BI102", T) True_YGROSS = read_params(f, "B92:BI92", T) @test maximum(abs, K .- True_K) ≈ 0. atol = 3.0e-11 # Relax the precision just for this variable @@ -102,21 +102,21 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # Net Economy tests - ABATECOST = m[:neteconomy, :ABATECOST] - C = m[:neteconomy, :C] - CPC = m[:neteconomy, :CPC] - CPRICE = m[:neteconomy, :CPRICE] - I = m[:neteconomy, :I] - Y = m[:neteconomy, :Y] - YNET = m[:neteconomy, :YNET] - - True_ABATECOST = read_params(f, "B97:BI97", T) - True_C = read_params(f, "B125:BI125", T) - True_CPC = read_params(f, "B126:BI126", T) - True_CPRICE = read_params(f, "B134:BI134", T) - True_I = read_params(f, "B101:BI101", T) - True_Y = read_params(f, "B98:BI98", T) - True_YNET = read_params(f, "B95:BI95", T) + ABATECOST = m[:neteconomy, :ABATECOST] + C = m[:neteconomy, :C] + CPC = m[:neteconomy, :CPC] + CPRICE = m[:neteconomy, :CPRICE] + I = m[:neteconomy, :I] + Y = m[:neteconomy, :Y] + YNET = m[:neteconomy, :YNET] + + True_ABATECOST = read_params(f, "B97:BI97", T) + True_C = read_params(f, "B125:BI125", T) + True_CPC = read_params(f, "B126:BI126", T) + True_CPRICE = read_params(f, "B134:BI134", T) + True_I = read_params(f, "B101:BI101", T) + True_Y = read_params(f, "B98:BI98", T) + True_YNET = read_params(f, "B95:BI95", T) @test maximum(abs, ABATECOST .- True_ABATECOST) ≈ 0. atol = Precision @test maximum(abs, C .- True_C) ≈ 0. atol = Precision @@ -129,22 +129,22 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # Radiative Forcing test FORC = m[:radiativeforcing, :FORC] - True_FORC = read_params(f, "B122:BI122", T) + True_FORC = read_params(f, "B122:BI122", T) @test maximum(abs, FORC .- True_FORC) ≈ 0. atol = Precision # Sea Level Rise tests - ThermSLR = m[:sealevelrise, :ThermSLR] - GSICSLR = m[:sealevelrise, :GSICSLR] - GISSLR = m[:sealevelrise, :GISSLR] - AISSLR = m[:sealevelrise, :AISSLR] - TotSLR = m[:sealevelrise, :TotSLR] + ThermSLR = m[:sealevelrise, :ThermSLR] + GSICSLR = m[:sealevelrise, :GSICSLR] + GISSLR = m[:sealevelrise, :GISSLR] + AISSLR = m[:sealevelrise, :AISSLR] + TotSLR = m[:sealevelrise, :TotSLR] - True_ThermSLR = read_params(f, "B178:BI178", T) - True_GSICSLR = read_params(f, "B179:BI179", T) - True_GISSLR = read_params(f, "B180:BI180", T) - True_AISSLR = read_params(f, "B181:BI181", T) - True_TotSLR = read_params(f, "B182:BI182", T) + True_ThermSLR = read_params(f, "B178:BI178", T) + True_GSICSLR = read_params(f, "B179:BI179", T) + True_GISSLR = read_params(f, "B180:BI180", T) + True_AISSLR = read_params(f, "B181:BI181", T) + True_TotSLR = read_params(f, "B182:BI182", T) @test maximum(abs, ThermSLR .- True_ThermSLR) ≈ 0. atol = Precision @test maximum(abs, GSICSLR .- True_GSICSLR) ≈ 0. atol = Precision @@ -154,13 +154,13 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # Welfare tests - CEMUTOTPER = m[:welfare, :CEMUTOTPER] - PERIODU = m[:welfare, :PERIODU] - UTILITY = m[:welfare, :UTILITY] + CEMUTOTPER = m[:welfare, :CEMUTOTPER] + PERIODU = m[:welfare, :PERIODU] + UTILITY = m[:welfare, :UTILITY] - True_CEMUTOTPER = read_params(f, "B129:BI129", T) - True_PERIODU = read_params(f, "B128:BI128", T) - True_UTILITY = read_params(f, "B130") + True_CEMUTOTPER = read_params(f, "B129:BI129", T) + True_PERIODU = read_params(f, "B128:BI128", T) + True_UTILITY = read_params(f, "B130") @test maximum(abs, CEMUTOTPER .- True_CEMUTOTPER) ≈ 0. atol = Precision @test maximum(abs, PERIODU .- True_PERIODU) ≈ 0. atol = Precision @@ -172,7 +172,7 @@ using MimiDICE2010: read_params, dice2010_excel_parameters # 3. Run tests to make sure integration version (Mimi v0.5.0) # values match Mimi 0.4.0 values # ------------------------------------------------------------------------------ - + @testset "dice2010-integration" begin Precision = 1.0e-11 @@ -190,7 +190,7 @@ using MimiDICE2010: read_params, dice2010_excel_parameters df = load(filepath) |> DataFrame if typeof(results) <: Number - validation_results = df[1,1] + validation_results = df[1, 1] else validation_results = Matrix(df) @@ -266,9 +266,9 @@ using MimiDICE2010: read_params, dice2010_excel_parameters :last_year => [2295, 2595], :pulse_size => [1e3, 1e7, 1e10] ]) - + results = DataFrame(year=[], eta=[], prtp=[], last_year=[], pulse_size=[], SC=[]) - + for year in specs[:year] for eta in specs[:eta] for prtp in specs[:prtp] @@ -281,7 +281,7 @@ using MimiDICE2010: read_params, dice2010_excel_parameters end end end - + validation_results = load(joinpath(@__DIR__, "..", "data", "SC validation data", "deterministic_sc_values_v1-0-1.csv")) |> DataFrame # diffs = sort(results[!, :SC] - validation_results[!, :SC], rev = true) # println(diffs) diff --git a/test/test_climatedynamics.jl b/test/test_climatedynamics.jl index c8298fd..2dea0c1 100644 --- a/test/test_climatedynamics.jl +++ b/test/test_climatedynamics.jl @@ -12,10 +12,10 @@ include("../src/components/climatedynamics_component.jl") add_comp!(m, climatedynamics, :climatedynamics) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :climatedynamics, :FORC, read_params(f, "B122:BI122", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :climatedynamics, :fco22x, p[:fco22x]) set_param!(m, :climatedynamics, :t2xco2, p[:t2xco2]) @@ -26,19 +26,19 @@ include("../src/components/climatedynamics_component.jl") set_param!(m, :climatedynamics, :c3, p[:c3]) set_param!(m, :climatedynamics, :c4, p[:c4]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables - TATM = m[:climatedynamics, :TATM] - TOCEAN = m[:climatedynamics, :TOCEAN] + # Extract the generated variables + TATM = m[:climatedynamics, :TATM] + TOCEAN = m[:climatedynamics, :TOCEAN] -# Extract the true values - True_TATM = read_params(f, "B121:BI121", T) + # Extract the true values + True_TATM = read_params(f, "B121:BI121", T) True_TOCEAN = read_params(f, "B123:BI123", T) -# Test that the values are the same + # Test that the values are the same @test maximum(abs, TATM .- True_TATM) ≈ 0. atol = Precision @test maximum(abs, TOCEAN .- True_TOCEAN) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_co2cycle.jl b/test/test_co2cycle.jl index 13cca6b..6e32377 100644 --- a/test/test_co2cycle.jl +++ b/test/test_co2cycle.jl @@ -12,10 +12,10 @@ include("../src/components/co2cycle_component.jl") add_comp!(m, co2cycle, :co2cycle) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :co2cycle, :E, read_params(f, "B109:BI109", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :co2cycle, :mat0, p[:mat0]) set_param!(m, :co2cycle, :mat1, p[:mat1]) @@ -29,25 +29,25 @@ include("../src/components/co2cycle_component.jl") set_param!(m, :co2cycle, :b32, p[:b32]) set_param!(m, :co2cycle, :b33, p[:b33]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables - MAT = m[:co2cycle, :MAT] - MAT_final = m[:co2cycle, :MAT_final] - ML = m[:co2cycle, :ML] - MU = m[:co2cycle, :MU] + # Extract the generated variables + MAT = m[:co2cycle, :MAT] + MAT_final = m[:co2cycle, :MAT_final] + ML = m[:co2cycle, :ML] + MU = m[:co2cycle, :MU] -# Extract the true values - True_MAT = read_params(f, "B112:BI112", T) - True_MAT_final = read_params(f, "BJ112") - True_ML = read_params(f, "B115:BI115", T) - True_MU = read_params(f, "B114:BI114", T) + # Extract the true values + True_MAT = read_params(f, "B112:BI112", T) + True_MAT_final = read_params(f, "BJ112") + True_ML = read_params(f, "B115:BI115", T) + True_MU = read_params(f, "B114:BI114", T) -# Test that the values are the same + # Test that the values are the same @test maximum(abs, MAT .- True_MAT) ≈ 0. atol = Precision @test abs(MAT_final - True_MAT_final) ≈ 0. atol = Precision @test maximum(abs, ML .- True_ML) ≈ 0. atol = Precision @test maximum(abs, MU .- True_MU) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_damages.jl b/test/test_damages.jl index acdcb8b..39afa9d 100644 --- a/test/test_damages.jl +++ b/test/test_damages.jl @@ -12,12 +12,12 @@ include("../src/components/damages_component.jl") add_comp!(m, damages, :damages) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :damages, :TATM, read_params(f, "B121:BI121", T)) set_param!(m, :damages, :YGROSS, read_params(f, "B92:BI92", T)) set_param!(m, :damages, :TotSLR, read_params(f, "B182:BI182", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :damages, :a1, p[:a1]) set_param!(m, :damages, :a2, p[:a2]) @@ -26,16 +26,16 @@ include("../src/components/damages_component.jl") set_param!(m, :damages, :b2, p[:b2]) set_param!(m, :damages, :b3, p[:b3]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables + # Extract the generated variables DAMFRAC = m[:damages, :DAMFRAC] -# Extract the true values - True_DAMFRAC = read_params(f, "B93:BI93", T) + # Extract the true values + True_DAMFRAC = read_params(f, "B93:BI93", T) -# Test that the values are the same + # Test that the values are the same @test maximum(abs, DAMFRAC .- True_DAMFRAC) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_emissions.jl b/test/test_emissions.jl index b5f13f2..2e7b922 100644 --- a/test/test_emissions.jl +++ b/test/test_emissions.jl @@ -12,31 +12,31 @@ include("../src/components/emissions_component.jl") add_comp!(m, emissions, :emissions) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :emissions, :YGROSS, read_params(f, "B92:BI92", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :emissions, :sigma, p[:sigma]) set_param!(m, :emissions, :MIU, p[:MIU]) set_param!(m, :emissions, :etree, p[:etree]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables - CCA = m[:emissions, :CCA] - E = m[:emissions, :E] - EIND = m[:emissions, :EIND] + # Extract the generated variables + CCA = m[:emissions, :CCA] + E = m[:emissions, :E] + EIND = m[:emissions, :EIND] -# Extract the true values - True_CCA = read_params(f, "B117:BI117", T) - True_E = read_params(f, "B109:BI109", T) - True_EIND = read_params(f, "B110:BI110", T) + # Extract the true values + True_CCA = read_params(f, "B117:BI117", T) + True_E = read_params(f, "B109:BI109", T) + True_EIND = read_params(f, "B110:BI110", T) -# Test that the values are the same + # Test that the values are the same @test maximum(abs, CCA .- True_CCA) ≈ 0. atol = Precision @test maximum(abs, E .- True_E) ≈ 0. atol = Precision @test maximum(abs, EIND .- True_EIND) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_grosseconomy.jl b/test/test_grosseconomy.jl index c98e3c7..f38cd5f 100644 --- a/test/test_grosseconomy.jl +++ b/test/test_grosseconomy.jl @@ -12,10 +12,10 @@ include("../src/components/grosseconomy_component.jl") add_comp!(m, grosseconomy, :grosseconomy) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :grosseconomy, :I, read_params(f, "B101:BI101", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :grosseconomy, :al, p[:al]) set_param!(m, :grosseconomy, :l, p[:l]) @@ -23,19 +23,19 @@ include("../src/components/grosseconomy_component.jl") set_param!(m, :grosseconomy, :dk, p[:dk]) set_param!(m, :grosseconomy, :k0, p[:k0]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables - K = m[:grosseconomy, :K] - YGROSS = m[:grosseconomy, :YGROSS] + # Extract the generated variables + K = m[:grosseconomy, :K] + YGROSS = m[:grosseconomy, :YGROSS] -# Extract the true values - True_K = read_params(f, "B102:BI102", T) + # Extract the true values + True_K = read_params(f, "B102:BI102", T) True_YGROSS = read_params(f, "B92:BI92", T) -# Test that the values are the same + # Test that the values are the same @test maximum(abs, K .- True_K) ≈ 0. atol = Precision @test maximum(abs, YGROSS .- True_YGROSS) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_neteconomy.jl b/test/test_neteconomy.jl index 73b6e38..5ac6168 100644 --- a/test/test_neteconomy.jl +++ b/test/test_neteconomy.jl @@ -12,11 +12,11 @@ include("../src/components/neteconomy_component.jl") add_comp!(m, neteconomy, :neteconomy) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :neteconomy, :YGROSS, read_params(f, "B92:BI92", T)) set_param!(m, :neteconomy, :DAMFRAC, read_params(f, "B93:BI93", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :neteconomy, :cost1, p[:cost1]) set_param!(m, :neteconomy, :MIU, p[:MIU]) @@ -26,28 +26,28 @@ include("../src/components/neteconomy_component.jl") set_param!(m, :neteconomy, :S, p[:S]) set_param!(m, :neteconomy, :l, p[:l]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables - ABATECOST = m[:neteconomy, :ABATECOST] - C = m[:neteconomy, :C] - CPC = m[:neteconomy, :CPC] - CPRICE = m[:neteconomy, :CPRICE] - I = m[:neteconomy, :I] - Y = m[:neteconomy, :Y] - YNET = m[:neteconomy, :YNET] - -# Extract the true values - True_ABATECOST = read_params(f, "B97:BI97", T) - True_C = read_params(f, "B125:BI125", T) - True_CPC = read_params(f, "B126:BI126", T) - True_CPRICE = read_params(f, "B134:BI134", T) - True_I = read_params(f, "B101:BI101", T) - True_Y = read_params(f, "B98:BI98", T) - True_YNET = read_params(f, "B95:BI95", T) - -# Test that the values are the same + # Extract the generated variables + ABATECOST = m[:neteconomy, :ABATECOST] + C = m[:neteconomy, :C] + CPC = m[:neteconomy, :CPC] + CPRICE = m[:neteconomy, :CPRICE] + I = m[:neteconomy, :I] + Y = m[:neteconomy, :Y] + YNET = m[:neteconomy, :YNET] + + # Extract the true values + True_ABATECOST = read_params(f, "B97:BI97", T) + True_C = read_params(f, "B125:BI125", T) + True_CPC = read_params(f, "B126:BI126", T) + True_CPRICE = read_params(f, "B134:BI134", T) + True_I = read_params(f, "B101:BI101", T) + True_Y = read_params(f, "B98:BI98", T) + True_YNET = read_params(f, "B95:BI95", T) + + # Test that the values are the same @test maximum(abs, ABATECOST .- True_ABATECOST) ≈ 0. atol = Precision @test maximum(abs, C .- True_C) ≈ 0. atol = Precision @test maximum(abs, CPC .- True_CPC) ≈ 0. atol = Precision @@ -56,4 +56,4 @@ include("../src/components/neteconomy_component.jl") @test maximum(abs, Y .- True_Y) ≈ 0. atol = Precision @test maximum(abs, YNET .- True_YNET) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_radiativeforcing.jl b/test/test_radiativeforcing.jl index 5ac888c..0f02482 100644 --- a/test/test_radiativeforcing.jl +++ b/test/test_radiativeforcing.jl @@ -12,25 +12,25 @@ include("../src/components/radiativeforcing_component.jl") add_comp!(m, radiativeforcing, :radiativeforcing) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :radiativeforcing, :MAT, read_params(f, "B112:BI112", T)) set_param!(m, :radiativeforcing, :MAT_final, read_params(f, "BJ112")) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :radiativeforcing, :forcoth, p[:forcoth]) set_param!(m, :radiativeforcing, :fco22x, p[:fco22x]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables + # Extract the generated variables FORC = m[:radiativeforcing, :FORC] -# Extract the true values - True_FORC = read_params(f, "B122:BI122", T) + # Extract the true values + True_FORC = read_params(f, "B122:BI122", T) -# Test that the values are the same + # Test that the values are the same @test maximum(abs, FORC .- True_FORC) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_sealevelrise.jl b/test/test_sealevelrise.jl index 15b4342..3603a84 100644 --- a/test/test_sealevelrise.jl +++ b/test/test_sealevelrise.jl @@ -12,10 +12,10 @@ include("../src/components/sealevelrise_component.jl") add_comp!(m, sealevelrise, :sealevelrise) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :sealevelrise, :TATM, read_params(f, "B121:BI121", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :sealevelrise, :therm0, p[:therm0]) set_param!(m, :sealevelrise, :gsic0, p[:gsic0]) @@ -31,28 +31,28 @@ include("../src/components/sealevelrise_component.jl") set_param!(m, :sealevelrise, :aisrate, p[:aisrate]) set_param!(m, :sealevelrise, :slrthreshold, p[:slrthreshold]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables - ThermSLR = m[:sealevelrise, :ThermSLR] - GSICSLR = m[:sealevelrise, :GSICSLR] - GISSLR = m[:sealevelrise, :GISSLR] - AISSLR = m[:sealevelrise, :AISSLR] - TotSLR = m[:sealevelrise, :TotSLR] - -# Extract the true values - True_ThermSLR = read_params(f, "B178:BI178", T) - True_GSICSLR = read_params(f, "B179:BI179", T) - True_GISSLR = read_params(f, "B180:BI180", T) - True_AISSLR = read_params(f, "B181:BI181", T) - True_TotSLR = read_params(f, "B182:BI182", T) - -# Test that the values are the same + # Extract the generated variables + ThermSLR = m[:sealevelrise, :ThermSLR] + GSICSLR = m[:sealevelrise, :GSICSLR] + GISSLR = m[:sealevelrise, :GISSLR] + AISSLR = m[:sealevelrise, :AISSLR] + TotSLR = m[:sealevelrise, :TotSLR] + + # Extract the true values + True_ThermSLR = read_params(f, "B178:BI178", T) + True_GSICSLR = read_params(f, "B179:BI179", T) + True_GISSLR = read_params(f, "B180:BI180", T) + True_AISSLR = read_params(f, "B181:BI181", T) + True_TotSLR = read_params(f, "B182:BI182", T) + + # Test that the values are the same @test maximum(abs, ThermSLR .- True_ThermSLR) ≈ 0. atol = Precision @test maximum(abs, GSICSLR .- True_GSICSLR) ≈ 0. atol = Precision @test maximum(abs, GISSLR .- True_GISSLR) ≈ 0. atol = Precision @test maximum(abs, AISSLR .- True_AISSLR) ≈ 0. atol = Precision @test maximum(abs, TotSLR .- True_TotSLR) ≈ 0. atol = Precision -end \ No newline at end of file +end diff --git a/test/test_welfare.jl b/test/test_welfare.jl index 4abdc15..786a786 100644 --- a/test/test_welfare.jl +++ b/test/test_welfare.jl @@ -12,10 +12,10 @@ include("../src/components/welfare_component.jl") add_comp!(m, welfare, :welfare) -# Set the parameters that would normally be internal connection from their Excel values + # Set the parameters that would normally be internal connection from their Excel values set_param!(m, :welfare, :CPC, read_params(f, "B126:BI126", T)) -# Load the rest of the external parameters + # Load the rest of the external parameters p = dice2010_excel_parameters(joinpath(@__DIR__, "..", "data", "DICE2010_082710d.xlsx")) set_param!(m, :welfare, :l, p[:l]) set_param!(m, :welfare, :elasmu, p[:elasmu]) @@ -23,22 +23,22 @@ include("../src/components/welfare_component.jl") set_param!(m, :welfare, :scale1, p[:scale1]) set_param!(m, :welfare, :scale2, p[:scale2]) -# Run the one-component model + # Run the one-component model run(m) -# Extract the generated variables - CEMUTOTPER = m[:welfare, :CEMUTOTPER] - PERIODU = m[:welfare, :PERIODU] - UTILITY = m[:welfare, :UTILITY] + # Extract the generated variables + CEMUTOTPER = m[:welfare, :CEMUTOTPER] + PERIODU = m[:welfare, :PERIODU] + UTILITY = m[:welfare, :UTILITY] -# Extract the true values + # Extract the true values True_CEMUTOTPER = read_params(f, "B129:BI129", T) - True_PERIODU = read_params(f, "B128:BI128", T) - True_UTILITY = read_params(f, "B130") + True_PERIODU = read_params(f, "B128:BI128", T) + True_UTILITY = read_params(f, "B130") -# Test that the values are the same + # Test that the values are the same @test maximum(abs, CEMUTOTPER .- True_CEMUTOTPER) ≈ 0. atol = Precision @test maximum(abs, PERIODU .- True_PERIODU) ≈ 0. atol = Precision @test abs(UTILITY - True_UTILITY) ≈ 0. atol = Precision -end \ No newline at end of file +end