add documentation for the module 15_climate (done by Jessi)

modules/15_climate/box/equations.gms

@@ -8,129 +8,192 @@
 ta2ttot10(ta10, ttot)$((ttot.val lt ta10.val) AND (pm_ttot_val(ttot+1) gt ta10.val)AND (ttot.val ge 2005)) = Yes;
 display ttot2ta10, ta2ttot10;
 
-
-*** taken from box-model.inc ---------------------------
-
-***     carbon cycle
-
-* AL changed towards an IRF Model with 2 or 3 time scales
+***---------------------------------------------------------------------------
+*'  Carbon Cycle. CO2 concentration is calculated using and Impulse-Response-Function Model with 3 time scales
+***---------------------------------------------------------------------------
 q15_cc(ta10)$(ord(ta10) ge 1)..
-                   v15_conc(ta10,'CO2') =e= s15_c0 + ((s15_c2000-s15_c0)/(s15_ca0*s15_cq0+s15_ca1*s15_cq1+s15_ca2*s15_cq2+s15_ca3*s15_cq3))*
-                                                                 (s15_ca0*s15_cq0 + s15_ca1*s15_cq1*exp(-(ord(ta10)-1)/s15_ctau1) + s15_ca2*s15_cq2*exp(-(ord(ta10)-1)/s15_ctau2) + s15_ca3*s15_cq3*exp(-(ord(ta10)-1)/s15_ctau3))+
-                                                                                                 s15_cconvi*sum((tx)$(ord(tx)lt ord(ta10)),v15_emi(tx,'CO2')*p15_epsilon(ta10-ord(tx)));
-
-***     N20 concentration (from ACC2)
-
-q15_concN2OQ(ta10) $(Ord(ta10) GT 1)     .. (v15_conc(ta10,'N2O')-v15_conc(ta10-1,'N2O'))/s15_DELTAT =E= 0.5*
-                                 ((1/s15_CNVN2O*(v15_emi(ta10,'N2O')+ s15_NATN2O)-(1/s15_TAUN2O*(v15_conc(ta10,'N2O')/s15_CONN2O2000R)**(-s15_SENTAUN2O)*v15_conc(ta10,'N2O')))+
-                                 (1/s15_CNVN2O*(v15_emi(ta10-1,'N2O')+ s15_NATN2O)-(1/s15_TAUN2O*(v15_conc(ta10-1,'N2O')/s15_CONN2O2000R)**(-s15_SENTAUN2O)*v15_conc(ta10-1,'N2O'))));
-*** end of stuff taken from box-model.inc ---------------------------
-
-*** taken from box-model.inc ---------------------------
-q15_concCH4Q(ta10) $(Ord(ta10) GT 1)     .. (v15_conc(ta10,'CH4')-v15_conc(ta10-1,'CH4'))/s15_DELTAT =E= 0.5*
-                                 ((1/s15_CNVCH4*(v15_emi(ta10,'CH4')+ s15_NATCH4)-
-                                 (p15_conroh(ta10)/s15_TAUCH4OH+1/s15_TAUCH4SS)*v15_conc(ta10,'CH4'))+
-                                 (1/s15_CNVCH4*(v15_emi(ta10-1,'CH4')+ s15_NATCH4)-
-                                 (p15_conroh(ta10-1)/s15_TAUCH4OH+1/s15_TAUCH4SS)*v15_conc(ta10-1,'CH4')));
-
-***    CH4 radiative forcing (from ACC2)
-
-q15_forcCH4Q(ta10)                    .. v15_forcComp(ta10,'CH4') =E= s15_RHOCH4*(SQRT(v15_conc(ta10,'CH4'))-SQRT(s15_CONCH4PRE))
-                                 -s15_OVERLFAC1*LOG(1+s15_OVERLFAC2*(v15_conc(ta10,'CH4')*s15_CONN2OPRE)**s15_OVERLEXP1+
-                                 s15_OVERLFAC3*v15_conc(ta10,'CH4')*(v15_conc(ta10,'CH4')*s15_CONN2OPRE)**s15_OVERLEXP2)
-                                 +s15_OVERLFAC1*LOG(1+s15_OVERLFAC2*(s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP1+
-                                 s15_OVERLFAC3*s15_CONCH4PRE*(s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP2);
-***    N2O radiative forcing (from ACC2)
-q15_forcN2OQ(ta10)                    .. v15_forcComp(ta10,'N2O') =E= s15_RHON2O*(SQRT(v15_conc(ta10,'N2O'))-SQRT(s15_CONN2OPRE))
-                                 -s15_OVERLFAC1*LOG(1+s15_OVERLFAC2*(s15_CONCH4PRE*v15_conc(ta10,'N2O'))**s15_OVERLEXP1+
-                                 s15_OVERLFAC3*s15_CONCH4PRE*(s15_CONCH4PRE*v15_conc(ta10,'N2O'))**s15_OVERLEXP2)
-                                 +s15_OVERLFAC1*LOG(1+s15_OVERLFAC2*(s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP1+
-                                 s15_OVERLFAC3*s15_CONCH4PRE*(s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP2);
-
-*AL* changed towards the ACC2 approach to use the foghg
-***     so2 radiative forcing module
-
-    q15_forcso2(ta10)..
-                    v15_forcComp(ta10,'SO2') =e= s15_dso1990 * v15_emi(ta10,'SO2') / s15_so1990 +
-                             s15_iso1990 * log(1 + v15_emi(ta10,'SO2')/s15_enatso2) /
-                                       log(1 + s15_so1990/s15_enatso2) ;
+    v15_conc(ta10,'CO2')
+	=e=
+	s15_c0 + ((s15_c2000-s15_c0)/(s15_ca0*s15_cq0+s15_ca1*s15_cq1+s15_ca2*s15_cq2+s15_ca3*s15_cq3))
+	         * (s15_ca0*s15_cq0 + s15_ca1*s15_cq1*exp(-(ord(ta10)-1)/s15_ctau1) + s15_ca2*s15_cq2*exp(-(ord(ta10)-1)/s15_ctau2) + s15_ca3*s15_cq3*exp(-(ord(ta10)-1)/s15_ctau3))
+		   + s15_cconvi * sum((tx)$(ord(tx)lt ord(ta10)),v15_emi(tx,'CO2')*p15_epsilon(ta10-ord(tx)));
+
+***---------------------------------------------------------------------------
+*'  CH4 concentration
+***---------------------------------------------------------------------------
+q15_concCH4Q(ta10)$(Ord(ta10) GT 1)..
+    (v15_conc(ta10,'CH4')-v15_conc(ta10-1,'CH4'))/s15_DELTAT
+	=e=
+	0.5 * (
+	      (1/s15_CNVCH4 * (v15_emi(ta10,'CH4')+ s15_NATCH4) - (p15_conroh(ta10) / s15_TAUCH4OH + 1 / s15_TAUCH4SS) * v15_conc(ta10,'CH4'))
+		  + (1/s15_CNVCH4*(v15_emi(ta10-1,'CH4') + s15_NATCH4) - (p15_conroh(ta10-1) / s15_TAUCH4OH + 1 / s15_TAUCH4SS) * v15_conc(ta10-1,'CH4'))
+		  );
+
+***---------------------------------------------------------------------------
+*'  N20 concentration (from ACC2)
+***---------------------------------------------------------------------------
+q15_concN2OQ(ta10)$(Ord(ta10) GT 1)..
+    (v15_conc(ta10,'N2O')-v15_conc(ta10-1,'N2O'))/s15_DELTAT
+	=e=
+	0.5 * (
+	(1/s15_CNVN2O * (v15_emi(ta10,'N2O')+ s15_NATN2O) - (1/s15_TAUN2O*(v15_conc(ta10,'N2O')/s15_CONN2O2000R)**(-s15_SENTAUN2O)*v15_conc(ta10,'N2O')))
+	+ (1/s15_CNVN2O * (v15_emi(ta10-1,'N2O') + s15_NATN2O) - (1/s15_TAUN2O*(v15_conc(ta10-1,'N2O')/s15_CONN2O2000R)**(-s15_SENTAUN2O)*v15_conc(ta10-1,'N2O')))
+	);
+
+***---------------------------------------------------------------------------
+*'    CO2 radiative forcing
+***---------------------------------------------------------------------------
+q15_forcco2(ta10)..
+    v15_forcComp(ta10,'CO2')
+	=e=
+	s15_fcodb * log(v15_conc(ta10,'CO2')/s15_c0) / log(2) ;
+
+***---------------------------------------------------------------------------
+*'    CH4 radiative forcing (from ACC2)
+***---------------------------------------------------------------------------
+q15_forcCH4Q(ta10)..
+    v15_forcComp(ta10,'CH4')
+	=e=
+	s15_RHOCH4 * (SQRT(v15_conc(ta10,'CH4')) - SQRT(s15_CONCH4PRE))
+    - s15_OVERLFAC1 * LOG(1 + s15_OVERLFAC2 *(v15_conc(ta10,'CH4') * s15_CONN2OPRE)**s15_OVERLEXP1
+	                        + s15_OVERLFAC3 * v15_conc(ta10,'CH4') * (v15_conc(ta10,'CH4')*s15_CONN2OPRE)**s15_OVERLEXP2
+							)
+    + s15_OVERLFAC1 * LOG(1 + s15_OVERLFAC2 * (s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP1
+	                        + s15_OVERLFAC3 * s15_CONCH4PRE * (s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP2
+							);
+							
+***---------------------------------------------------------------------------
+*'    N2O radiative forcing (from ACC2)
+***---------------------------------------------------------------------------
+q15_forcN2OQ(ta10)..
+    v15_forcComp(ta10,'N2O')
+	=e=
+	s15_RHON2O * (SQRT(v15_conc(ta10,'N2O'))-SQRT(s15_CONN2OPRE))
+    - s15_OVERLFAC1 * LOG(1 + s15_OVERLFAC2 * (s15_CONCH4PRE*v15_conc(ta10,'N2O'))**s15_OVERLEXP1
+	                      + s15_OVERLFAC3 * s15_CONCH4PRE * (s15_CONCH4PRE*v15_conc(ta10,'N2O'))**s15_OVERLEXP2
+						  )
+    + s15_OVERLFAC1 * LOG(1 + s15_OVERLFAC2 * (s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP1
+	                      + s15_OVERLFAC3 * s15_CONCH4PRE * (s15_CONCH4PRE*s15_CONN2OPRE)**s15_OVERLEXP2
+						  );
+
+***---------------------------------------------------------------------------
+*'    SO2 radiative forcing
+***---------------------------------------------------------------------------
+q15_forcso2(ta10)..
+    v15_forcComp(ta10,'SO2')
+	=e=
+	s15_dso1990 * v15_emi(ta10,'SO2') / s15_so1990 
+	+ s15_iso1990 * log(1 + v15_emi(ta10,'SO2')/s15_enatso2) / log(1 + s15_so1990/s15_enatso2) ;
 									   
-*** BC and OC from fossil fuels forcing
-
-	q15_forcbc(ta10)..
-					v15_forcComp(ta10,'BC') =e= (v15_emi(ta10,'BC') / s15_bc2005) * p15_oghgf_ffbc('2005');
+***---------------------------------------------------------------------------
+*'    BC and OC from fossil fuels radiative forcing; scales linear with emissions.
+***---------------------------------------------------------------------------
+q15_forcbc(ta10)..
+	v15_forcComp(ta10,'BC')
+	=e=
+	(v15_emi(ta10,'BC') / s15_bc2005) * p15_oghgf_ffbc('2005');
 					
-	q15_forcoc(ta10)..
-					v15_forcComp(ta10,'OC') =e= (v15_emi(ta10,'OC') / s15_oc2005) * p15_oghgf_ffoc('2005');					
-
-
-
-***     co2 radiative forcing module
-
-  q15_forcco2(ta10)..
-                   v15_forcComp(ta10,'CO2') =e= s15_fcodb * log(v15_conc(ta10,'CO2')/s15_c0) / log(2) ;
-
-***     total radiative forcing (foghg given exogenously)
-
-  q15_forctotal(ta10)..
-                   v15_forcComp(ta10,'TTL') =e= v15_forcComp(ta10,'CO2') +
-                                v15_forcComp(ta10,'SO2') +
-                                v15_forcComp(ta10,'CH4') +
-                                v15_forcComp(ta10,'N2O') +
-                                v15_forcComp(ta10,'BC') + 
-                                v15_forcComp(ta10,'OC') + 
-                                v15_forcComp(ta10,'oghg_nokyo') +
-                                v15_forcComp(ta10,'oghg_kyo');
+q15_forcoc(ta10)..
+	v15_forcComp(ta10,'OC')
+	=e=
+	(v15_emi(ta10,'OC') / s15_oc2005) * p15_oghgf_ffoc('2005');					
+
+***---------------------------------------------------------------------------
+*'    Total radiative forcing (foghg given exogenously)
+***---------------------------------------------------------------------------
+q15_forctotal(ta10)..
+    v15_forcComp(ta10,'TTL')
+	=e=
+	hv15_forcComp(ta10,'CO2') 
+	+ v15_forcComp(ta10,'SO2') 
+	+ v15_forcComp(ta10,'CH4') 
+	+ v15_forcComp(ta10,'N2O') 
+	+ v15_forcComp(ta10,'BC') 
+	+ v15_forcComp(ta10,'OC') 
+	+ v15_forcComp(ta10,'oghg_nokyo') 
+	+ v15_forcComp(ta10,'oghg_kyo');
 								
-*GL: Total forcing of Kyoto gases - required for formulation of policy targets
-  q15_forc_kyo(ta10)..
-                   v15_forcKyo(ta10) =e= v15_forcComp(ta10,'CO2') +
-                                v15_forcComp(ta10,'CH4') +
-                                v15_forcComp(ta10,'N2O') +
-                                 v15_forcComp(ta10,'oghg_kyo');
-
-*JeS* rcp forcing
-	q15_forc_rcp(ta10)..
-					v15_forcRcp(ta10) =e= v15_forcComp(ta10,'CO2') +
-                                v15_forcComp(ta10,'CH4') +
-                                v15_forcComp(ta10,'N2O') +
-								v15_forcComp(ta10,'oghg_kyo') + 
-								v15_forcComp(ta10,'SO2') + 
-								v15_forcComp(ta10,'BC') + 
-								v15_forcComp(ta10,'OC') +
-								v15_forcComp(ta10,'oghg_nokyo_rcp');
+***---------------------------------------------------------------------------
+*'    Total radiative forcing of Kyoto gases
+***---------------------------------------------------------------------------
+q15_forc_kyo(ta10)..
+    v15_forcKyo(ta10) 
+	=e=
+	v15_forcComp(ta10,'CO2') 
+	+ v15_forcComp(ta10,'CH4') 
+	+ v15_forcComp(ta10,'N2O') 
+	+ v15_forcComp(ta10,'oghg_kyo');
+
+***---------------------------------------------------------------------------
+*'    RCP forcing
+***---------------------------------------------------------------------------
+q15_forc_rcp(ta10)..
+	v15_forcRcp(ta10)
+	=e=
+	v15_forcComp(ta10,'CO2') 
+	+ v15_forcComp(ta10,'CH4') 
+	+ v15_forcComp(ta10,'N2O') 
+	+ v15_forcComp(ta10,'oghg_kyo') 
+	+ v15_forcComp(ta10,'SO2') 
+	+ v15_forcComp(ta10,'BC') 
+	+ v15_forcComp(ta10,'OC') 
+	+ v15_forcComp(ta10,'oghg_nokyo_rcp');
 																
-***     new temperature equations
-
-q15_clisys1(ta10) $(Ord(ta10) > 1) .. (v15_tempFast(ta10)-v15_tempFast(ta10-1))/s15_deltat_box =E= 0.5/s15_RPCTT1*
-                                         ((s15_RPCTA1*v15_forcComp(ta10,'TTL')/3.7-v15_tempFast(ta10))+
-                                          (s15_RPCTA1*v15_forcComp(ta10-1,'TTL')/3.7-v15_tempFast(ta10-1)));
-
-q15_clisys2(ta10) $(Ord(ta10) > 1) .. (v15_tempSlow(ta10)-v15_tempSlow(ta10-1))/s15_deltat_box =E= 0.5/s15_RPCTT2*
-                                         ((s15_RPCTA2*v15_forcComp(ta10,'TTL')/3.7-v15_tempSlow(ta10))+
-                                          (s15_RPCTA2*v15_forcComp(ta10-1,'TTL')/3.7-v15_tempSlow(ta10-1)));
-
-q15_clisys(ta10) .. v15_temp(ta10) =E= v15_tempFast(ta10) + v15_tempSlow(ta10);
-
-q15_clisys01(ta10) $(ord(ta10)=1)..
-                    v15_tempFast(ta10) =e= s15_RPCTA1 * s15_temp2000/s15_tsens;
-
-q15_clisys02(ta10) $(ord(ta10)=1)..
-                    v15_tempSlow(ta10) =e= s15_RPCTA2 * s15_temp2000/s15_tsens;
-
-*JeS* forcing overshoot for damage function
+***---------------------------------------------------------------------------
+*'    Temperature equations, consisting of a fast and a slow response function
+***---------------------------------------------------------------------------
+q15_clisys01(ta10)$(ord(ta10)=1)..
+    v15_tempFast(ta10)
+	=e=
+	s15_RPCTA1 * s15_temp2000 / s15_tsens;
+
+q15_clisys02(ta10)$(ord(ta10)=1)..
+    v15_tempSlow(ta10)
+	=e=
+	s15_RPCTA2 * s15_temp2000 / s15_tsens;
+
+q15_clisys1(ta10)$(Ord(ta10) > 1)..
+    (v15_tempFast(ta10) - v15_tempFast(ta10-1)) / s15_deltat_box
+    =e=
+    0.5 / s15_RPCTT1 * ( (s15_RPCTA1 * v15_forcComp(ta10,'TTL') / 3.7 - v15_tempFast(ta10))
+	                   + (s15_RPCTA1 * v15_forcComp(ta10-1,'TTL') / 3.7 - v15_tempFast(ta10-1))
+					   );
+
+q15_clisys2(ta10)$(Ord(ta10) > 1)..
+    (v15_tempSlow(ta10) - v15_tempSlow(ta10-1)) / s15_deltat_box
+	=e=
+	0.5 / s15_RPCTT2 * ( (s15_RPCTA2 * v15_forcComp(ta10,'TTL') / 3.7 - v15_tempSlow(ta10))
+	                   + (s15_RPCTA2 * v15_forcComp(ta10-1,'TTL') / 3.7 - v15_tempSlow(ta10-1))
+					   );
+
+q15_clisys(ta10)..
+    v15_temp(ta10)
+	=e=
+	v15_tempFast(ta10) + v15_tempSlow(ta10);
+
+***---------------------------------------------------------------------------
+*'    Forcing overshoot for damage function
+***---------------------------------------------------------------------------
 q15_forc_os(t)..
-        vm_forcOs(t) =E= v15_forcKyo(t)-s15_gr_forc_kyo+v15_slackForc(t);
-*** end of stuff taken from box-model.inc ---------------------------
-
-*-----------------------------------------------------------------------------------------
-*** ------------------------link to core--------------------------------------------------
-*-----------------------------------------------------------------------------------------
-
-q15_linkEMI(ttot2ta10(ttot, ta10),emis2climate10(enty,FOB10)).. vm_emiAllGlob(ttot,enty) =e= v15_emi(ta10,FOB10);
+    vm_forcOs(t)
+	=e=
+	v15_forcKyo(t) - s15_gr_forc_kyo + v15_slackForc(t);
+
+***-----------------------------------------------------------------------------------------
+*'     link to core
+***-----------------------------------------------------------------------------------------
+q15_linkEMI(ttot2ta10(ttot, ta10),emis2climate10(enty,FOB10))..
+    vm_emiAllGlob(ttot,enty)
+	=e=
+	v15_emi(ta10,FOB10);
 $IF %cm_so2_out_of_opt% == "on" q15_linkEMI_aer(ttot2ta10(ttot, ta10),emiaer2climate10(emiaer,FOB10)).. p15_so2emi(ttot,emiaer) =e= v15_emi(ta10,FOB10);
 
-*---interpolation for linking (annual resolution of climate modules)
-q15_interEMI(ta2ttot10(ta10, ttot),FOBEMI(FOB10)).. v15_emi(ta10,FOB10) =e= (1-p15_interpol(ta10)) * v15_emi(ttot,FOB10) + p15_interpol(ta10) *  v15_emi(ttot+1,FOB10);
+***-----------------------------------------------------------------------------------------
+*'     interpolation for linking (annual resolution of climate modules)
+***-----------------------------------------------------------------------------------------
+q15_interEMI(ta2ttot10(ta10, ttot),FOBEMI(FOB10))..
+    v15_emi(ta10,FOB10)
+	=e=
+	(1-p15_interpol(ta10)) * v15_emi(ttot,FOB10) + p15_interpol(ta10) *  v15_emi(ttot+1,FOB10);
 
 *** EOF ./modules/15_climate/box/equations.gms

modules/15_climate/box/realization.gms

@@ -6,6 +6,9 @@
 *** |  Contact: remind@pik-potsdam.de
 *** SOF ./modules/15_climate/box.gms
 
+*' @description 
+*' In this realization, concentration, forcing, and temperature values are calculated using a simple model that can be used within the optimization routine.
+
 *####################### R SECTION START (PHASES) ##############################
 $Ifi "%phase%" == "sets" $include "./modules/15_climate/box/sets.gms"
 $Ifi "%phase%" == "declarations" $include "./modules/15_climate/box/declarations.gms"

modules/15_climate/magicc/postsolve.gms

@@ -6,45 +6,50 @@
 *** |  Contact: remind@pik-potsdam.de
 *** SOF ./modules/15_climate/magicc/postsolve.gms
 
-*cb 20140218 HERE goes the code for the iterative adjustment of the emission budget for SSP runs
-*** emission budgets are adjusted, such that a predefined forcing target in 2100 is met
-*** the actual 2100 forcing after each iteration is calculated by a magicc run started from GAMS
+***---------------------------------------------------------------------------
+*' HERE goes the code for the iterative adjustment of the emission budget for SSP runs
+*' emission budgets are adjusted, such that a predefined forcing target in 2100 is met
+*' the actual 2100 forcing after each iteration is calculated by a magicc run started from GAMS
+***---------------------------------------------------------------------------
 
 *** Generate MAGICC scenario file
 $include "./core/magicc.gms";
-
-* execute MAGICC (this is cheap enough, ~2s)
+*** execute MAGICC (this is cheap enough, ~2s)
 Execute "Rscript run_magicc.R";
-* read in results
+*** read in results
 Execute "Rscript read_DAT_TOTAL_ANTHRO_RF.R";
 Execute_Loadpoint 'p15_forc_magicc'  p15_forc_magicc;
 Execute "Rscript read_DAT_SURFACE_TEMP.R";
 Execute_Loadpoint 'p15_magicc_temp' pm_globalMeanTemperature = pm_globalMeanTemperature;
-* MAGICC only reports unitl 2300:
+*** MAGICC only reports unitl 2300:
 pm_globalMeanTemperature(tall)$(tall.val gt 2300) = 0;
 
-
+***---------------------------------------------------------------------------
+*' calibrate temperature (GMT anomaly) to match HADCRUT4 in 2000. 
+*' This ensures that different MAGICC configurations start at the same observed temperature.
+***---------------------------------------------------------------------------
 $ifthen.cm_magicc_calibrateTemperature2000 %cm_magicc_calibrateTemperature2000% == "HADCRUT4"
-*AJS* calibrate temperature (GMT anomaly) to match HADCRUT4 in 2000. This ensures that different MAGICC configurations start at the same observed temperature.
-*OLD* The HADCRUT4 offset from 1861-1880 (the SR1.5 reference period) to 1985-2015 is  0.67 degree Celsius (median).
-*UPDATE* Use the 2010 offset relative to 2006-2015 reference period (AR6SR Table 2.2, footnote 1; email GL from 19122018)
+
+***---------------------------------------------------------------------------
+*' Calibrate temperature such that anomaly in 2006-2015 reference period is 0.97 (SR1.5 Table 2.2, footnote 1)
+***---------------------------------------------------------------------------
 s15_tempOffset2010 = sum(tall$(tall.val gt 2005 and tall.val le 2015),pm_globalMeanTemperature(tall))/10; 
 display s15_tempOffset2010;
 pm_globalMeanTemperature(tall) = pm_globalMeanTemperature(tall) - s15_tempOffset2010 + 0.97;
 display pm_globalMeanTemperature;
 
-*temperature convergence indicator
+*** temperature convergence indicator
 p15_gmt_conv = 100*smax(t,abs(pm_globalMeanTemperature(t)/max(p15_gmt0(t),1e-8) -1));
 display p15_gmt_conv;
-*save temp from last iteration
+*** save temp from last iteration
 p15_gmt0(tall) = pm_globalMeanTemperature(tall);
 
 $endif.cm_magicc_calibrateTemperature2000
 
-* offset from HADCRUT4 to zero temperature in 1900, instead of the default 1870 (20 year averages each).
+*** offset from HADCRUT4 to zero temperature in 1900, instead of the default 1870 (20 year averages each).
 pm_globalMeanTemperatureZeroed1900(tall)  = pm_globalMeanTemperature(tall) + 0.092; 
 
-* derive temperature impulse response (TIRF) from MAGICC pulse scenarios
+*** derive temperature impulse response (TIRF) from MAGICC pulse scenarios
 $ifthen.cm_magicc_tirf "%cm_magicc_temperatureImpulseResponse%" == "on"
 * the TIRF does not change much with emissions profile (see, e.g., figure in Schultes et al. 2017); 
 * thus only compute TIRF after each of the first 10 iterations, then only every fifth iteration. 
@@ -56,16 +61,22 @@ if( ((iteration.val le 10) or ( mod(iteration.val,5 ) eq 0)) ,
 *NOTE the MAGICC results (*.OUT files) are from  the last pulse experiment now, so take care if reading them in after this point.
 $endif.cm_magicc_tirf
 
-*** Iterative adjustment based on 
+***---------------------------------------------------------------------------
+*' Iterative adjustment of budgets or carbon taxes to meet forcing target 
+***---------------------------------------------------------------------------
 if (cm_iterative_target_adj eq 2, !! otherwise adjustment happens in core/postsolve.gms 
   
-*** Iterative adjustment for budget runs  
+***---------------------------------------------------------------------------
+*' Iterative adjustment for budget runs: scale current budget with the ratio of target forcing s15_gr_forc_os to current forcing p15_forc_magicc.
+*' The offset is only there to increase the speed of convergence, the values have no physical meaning.
+*' For low stabilization targets (rcp2.0, rcp2.6, rcp3.7) the target is the 2100 forcing target (s15_rcpCluster eq 1),
+*' for lower targets the forcing target is valid during the full century (s15_rcpCluster eq 0).
+***---------------------------------------------------------------------------
   if ((cm_emiscen eq 6),
    
       display sm_budgetCO2eqGlob, s15_gr_forc_os, p15_forc_magicc;
    
       if (s15_rcpCluster eq 1,
-display ' Liebe Lavinia ';
         sm_budgetCO2eqGlob 
         = 
           sm_budgetCO2eqGlob 
@@ -98,7 +109,12 @@ display ' Liebe Lavinia ';
      display sm_budgetCO2eqGlob;
      );
 
-*** Iterative adjustments for tax runs
+***---------------------------------------------------------------------------
+*' Iterative adjustment for carbon tax runs: scale current tax pathway with the ratio of target forcing s15_gr_forc_os to current forcing p15_forc_magicc.
+*' The offset is only there to increase the speed of convergence, the values have no physical meaning.
+*' For low stabilization targets (rcp2.0, rcp2.6, rcp3.7) the target is the 2100 forcing target (s15_rcpCluster eq 1),
+*' for lower targets the forcing target is valid during the full century (s15_rcpCluster eq 0).
+***---------------------------------------------------------------------------
   if (cm_emiscen eq 9, 
   
     display pm_taxCO2eq, s15_gr_forc_os, p15_forc_magicc;

modules/15_climate/magicc/realization.gms

@@ -6,6 +6,10 @@
 *** |  Contact: remind@pik-potsdam.de
 *** SOF ./modules/15_climate/magicc.gms
 
+*' @description 
+*' In this realization, concentration, forcing, and temperature values are calculated using a MAGICC6.4. 
+*' MAGICC is run in between iterations and can be used to adapt carbon tax pathways and budgets to meet a give climate target.
+
 *####################### R SECTION START (PHASES) ##############################
 $Ifi "%phase%" == "sets" $include "./modules/15_climate/magicc/sets.gms"
 $Ifi "%phase%" == "declarations" $include "./modules/15_climate/magicc/declarations.gms"

modules/15_climate/module.gms

@@ -6,6 +6,12 @@
 *** |  Contact: remind@pik-potsdam.de
 *** SOF ./modules/15_climate/15_climate.gms
 
+*' @title climate
+*'
+*' @description  The 15_climate module calculates the resulting climate variables using either MAGICC6.4 or a stylized box model that can be used within the optimization routine.
+*'
+*' @authors Jessica Strefler, Michaja Pehl, Christoph Bertram
+
 *###################### R SECTION START (MODULETYPES) ##########################
 $Ifi "%climate%" == "box" $include "./modules/15_climate/box/realization.gms"
 $Ifi "%climate%" == "magicc" $include "./modules/15_climate/magicc/realization.gms"