* This version (20191010) is an updated version from the original codes developed by Aaheim and Rive 2005
* The required input data of GTAP v.9 can be downloaded from https://www.gtap.agecon.purdue.edu/databases/
* All other input and output data are included along with the codes.
* The aggregated GTAP data "foodpris51.gdx" should be obtained by the aggregation codes in the folder "gtap9gams"
* For each climate change scenario, the simulated restults can be obtained by alter the input file of data on changes in wheat yields in the file "scenarios_declare.gms". The input data will be redefine the values of the parameter "YCh(sr,t,r)".
* For each climate change scenario, the reported data can be obtained by running the codes in the file "resdata.gms" based on simulated results by running the file "core.gms"
* The reported data by scenario are then copied to excel files for further treatment.
* Please contact Taoyuan Wei (taoyuan.wei@cicero.oslo.no) for any comments and suggestions


$if not set ds $set ds foodpris51

$INCLUDE GTAP_input_header
$INCLUDE GTAP_input_main

$INCLUDE Zvalues_header
$INCLUDE Zvalues_main
$INCLUDE scenarios_declare
*$exit
$INCLUDE emissions
$INCLUDE elasticities
$INCLUDE closure
$INCLUDE time_declare
$INCLUDE mpsge_all
*$exit
** check balnace
*$INCLUDE balanchk

** single run to calibrate the basedata and provide initial data for loop
Parameter modstat(SIM)  model statement;
*parameter tem;

$if exist "%ds%_res_BAU.GDX" $GDXIN "%ds%_res_BAU.GDX"
$if exist "%ds%_res_BAU.GDX" $LOAD PERMIT_T
$if exist "%ds%_res_BAU.GDX" $GDXIN

loop(SIM$(ord(SIM)=1),
         Loop(t$(ord(t)=1),
                 calib_toggle = 1 ;
                 calib_KHAT(r) = 0;
*         RORE.FX(r) = sum(rr,ROREZ(rr))/card(rr)+1;
                 KHAT0.L(r) = KHAT(r) ;
*                 INV.LO(r)  = 0.95 ;
*                 INV.FX(r)  = 1 ;
$INCLUDE GRACE.GEN
GLOBTRST.FX= GLOBTRUST;
*                  GRACE.iterlim =2e+9;
                  GRACE.RESLIM = 10000000 ;
                 GRACE.optfile =1;
                  SOLVE GRACE USING MCP;
                  IF (((GRACE.modelstat <> 1) OR (GRACE.solvestat <> 1)),
                         modstat(SIM) = 0;
                         abort$(modstat(SIM)=0) "Couldn't solve!";
                  ELSE modstat(SIM) = 1;
                         );
display KHAT0.L, KHAT, RORE.L;
                 calib_toggle = 1 ;
$INCLUDE GRACE.GEN
GLOBTRST.FX= GLOBTRUST;
                 option savepoint=1;
                  SOLVE GRACE USING MCP;
                  IF (((GRACE.modelstat <> 1) OR (GRACE.solvestat <> 1)),
                         modstat(SIM) = 0;
                         abort$(modstat(SIM)=0) "Couldn't solve!";
                  ELSE modstat(SIM) = 1;
                         );
$ontext
loop(tmax$(sum((f,r)$(not fnat(f)),sqr(EFFTFP.L(f,r)-1))>0.000001),
eff_if(f,a,i,r)$(not fnat(f)) =  eff_if(f,a,i,r)*(1+EFFTFP.L(f,r))/2 ;
$INCLUDE GRACE.GEN
GLOBTRST.FX= GLOBTRUST;
                 option savepoint=1;
                  GRACE.RESLIM = 10000000 ;
                  SOLVE GRACE USING MCP;
                  IF (((GRACE.modelstat <> 1) OR (GRACE.solvestat <> 1)),
                         modstat(SIM) = 0;
                         abort$(modstat(SIM)=0) "Couldn't solve!";
                  ELSE modstat(SIM) = 1;
                         );
     );
$offtext
$INCLUDE time_define
display KHAT0.L, KHAT, RORE.L;
                 );
         IO_ABSZ(y,a,yy,r) = IO_ABS.L(y,a,yy,r) ;
         FD_ABSZ(y,j,r) = FD_ABS.L(y,j,r) ;
         EFFTFPZ(f,r) = EFFTFP.L(f,r) ;
         EFFTFP_SZ(f,y,r) = EFFTFP_S.L(f,y,r) ;
         CAPAZ(a,y,r) = CAPA.L(a,y,r) ;
         TS_PrimenZ(y,r) = TS_Primen.L(y,r) ;
         KHAT(r) = KHAT0.L(r) ;
         calib_KHAT(r) = 0;

execute_unload 'Zvalues.gdx', PFAZ,tFD,PFDZ,PIMZ,PDAZ,CAPAZ
         P_IOZ,PTMZ,dDEP,KHAT,EXRZ,FDZ, FDZ,TRNZ,FTOT,GLOBTRUST,GLOBTRUSTSHARE
         IFZ,IFAZ,IMZ,IMTOTZ,INVZ,IOZ,IOAZ, TRNZ, XDZ,REGINCZ,SAVZ, VKBZ,GDPZ,EMSTOTZ,FOSSZ
         VKEZ,TRANRZ,XZ,XDAZ,XDDZ,SAVINGS,GDPgZ,EFFTFPZ,EFFTFP_SZ,TS_PrimenZ
         PPZ,PERMITZ,FreeF_SH,FreeF_SH_noC,FreeF_SH_noF,NFreeF_SH,C_toggle,calib_toggle,Csub_toggle,calib_capa,calib_primen,calib_GDP,calib_KHAT
         evf,evh,eR,eM,eVA,eKLE,eENE,eNEL,eJTOP,eJENE,eJNEN,eARM,eIMP,ty,ti,tf,tx,tm
         eff_xd,eff,eff_f,eff_if,eghga_io,eghga_hh,IO_ABSZ,FD_ABSZ , ROREZ, RORGFLEX
*        eghga_FD,eghga_f,eghga_xd
*         c_CC,CC_par,CC_pede,tmp,dtmp,pcp,dpcp,permit45z,CO2_ffg
         Prim_eng,Prim_enZ,Ctrad_toggle  ,eRG
         VKE.L,RORE.L,INV.L,PINDEX.L,PINDEXW.L,GDP.L,EMSTOT.L,FOSS.L,EFFTFP.L,EFFTFP_S.L,PERMIT.L,tflnd.L,gene_gz,gene_g_t
         ;
         );
** end of single run
*display PERMIT_T
*$exit

*        Loop Section for multiple run
option solprint=off
option savepoint=0;

scalar tem /1/;
LOOP(SIM$SIMCUR(SIM),
*        eghga_abz_T(r,t)$((ord(SIM)=5) and (ord(t)>7))= eghga_abz_T(r,t)*0.989**(ord(t)-7) ;
         Loop(t$(ord(t)<=1),
*         Loop(t$((ord(SIM)<=5 and ord(t)<=30) or (ord(SIM)>=6 and ord(t)<=22)),
*         Loop(t$(ord(t)<=2 or mod(ord(t),5)=2),
                  If(ord(t)  =1,
                          t_size = 1;
                          execute_load 'Zvalues.gdx', PFAZ,tFD,PFDZ,PIMZ,PDAZ,CAPAZ
                          P_IOZ,PTMZ,dDEP,KHAT,EXRZ,FDZ, FDZ,TRNZ,FTOT,GLOBTRUST,GLOBTRUSTSHARE
                          IFZ,IFAZ,IMZ,IMTOTZ,INVZ,IOZ,IOAZ, TRNZ, XDZ,REGINCZ,SAVZ, VKBZ,GDPZ,EMSTOTZ,FOSSZ
                          VKEZ,TRANRZ,XZ,XDAZ,XDDZ,SAVINGS,GDPgZ,EFFTFPZ,EFFTFP_SZ,TS_PrimenZ ,eRG
                          PPZ,PERMITZ,FreeF_SH,FreeF_SH_noC,FreeF_SH_noF,NFreeF_SH, C_toggle,calib_toggle,Csub_toggle,calib_capa,calib_primen,calib_GDP,calib_KHAT
                          evf,evh,eR,eM,eVA,eKLE,eENE,eNEL,eJTOP,eJENE,eJNEN,eARM,eIMP,ty,ti,tf,tx,tm
                          eff_xd,eff,eff_f,eff_if,eghga_io,eghga_hh,gene_gz,IO_ABSZ,FD_ABSZ,ROREZ, RORGFLEX
*                          eghga_FD,eghga_f,eghga_xd
*                          c_CC,CC_par,tmp,dtmp,pcp,dpcp,Prim_enZ
                          Prim_eng,Prim_enZ,Ctrad_toggle  ,eRG
                          VKE.L,RORE.L,INV.L,PINDEX.L,PINDEXW.L,GDP.L,EMSTOT.L,FOSS.L,EFFTFP.L,EFFTFP_S.L,PERMIT.L,tflnd.L,gene_gz,gene_g_t
                                  ;
                          execute_loadpoint 'GRACE_p.gdx';
$if exist efftfp.gdx calib_toggle = 0;

         eARM(i,r) = eARM(i,r)*0.5;

                  else
                           t_size = 1;
                         );
                 calib_KHAT(r) = 0;
GLOBTRST.FX= GLOBTRUST;

$INCLUDE scenarios_define
$INCLUDE GRACE.GEN
                  GRACE.iterlim =2e+9;
                  GRACE.RESLIM = 10000000 ;
                  GRACE.optfile = 1   ;
                  SOLVE GRACE USING MCP;

                  IF (((GRACE.modelstat <> 1) OR (GRACE.solvestat <> 1)),
                          modstat(SIM) = 0;
                          abort$(modstat(SIM)=0) "Couldn't solve!";
                  ELSE modstat(SIM) = 1;
                          );
$INCLUDE time_define
*$ontext

if ((calib_toggle = 1),
         execute_unload 'efftfp.gdx', f,r,t,EFFTFP_T = efftfpz_T, EFFTFP_S_T = efftfp_sz_T, eff_T,eff_f_T,eff_if_T=eff_if_T_,eghga_abz_T
                                 TS_Primen_T=TS_PrimenZ_T,CAPA_T = CAPAZ_T;
         );
execute_unload '%ds%_res2CCO205.gdx', y, r, f, j,jr, ghga, SIM,calib_toggle
        FD_T,FD_ABS_T,IFAF_T,IFAF_noC_T,IFAF_noF_T,IFANFr_T, IM_T,IO_ABS_T,IF_ABS_T,IF_ABSr_T, TRANS_T,X_T,XD_T, XDA_T
        P_C_T,P_IO_T,P_T,PCAP_T,PD_T,PDA_T,PF_T,PFF_T,PFF_noC_T,PFF_noF_T,PFNF_T,PFA_T,PFD_T,PI_T,PIM_T,PS_T,PT_T,PTRANS_T, PP_T,PI_noC_T,PI_noF_T
        GLOBANK_T,GLOBTRST_T,INC_T,REGINC_T
        INV_T,RORE_T,ROREZ_T,VKE_T,GDP_T,EMSTOT_T,PINDEX_T,PINDEXW_T,FOSS_T, PERMIT_T, tflnd_T,KHAT_T, FreeF_SH_T,FreeF_SH_noC_T,FreeF_SH_noF_T

        FDr_T,FD_ABSr_T,IFAr_T,IMr_T,IO_ABSr_T,TRANSr_T,Xr_T,XDr_T,XDAr_T,PERMITr_T,PERMITrN_T
        P_IOr_T,PFAr_T
        INVr_T,VKEr_T,GDPr_T,EMSTOTr_T,FOSSr_T,FTOT_T

        EXRr_T,VTWRr_T,XDDr_T,TRANRr_T

        SAVr_T,IMc_T,EXc_T,EX_IM_BAL_T                ,eRG  ,calib_capa,calib_primen,calib_GDP, calib_KHAT

        GLOBTRUSTSHARE_T,eff_xd_T,eff_if_T,eff_T,eff_f_T,NFreeF_SH_T,EFFTFP_T,EFFTFP_S_T,CAPA_T,TS_Primen_T
        eghga_xd_T, eghga_hh2_T,eghga_FD_T,eghga_io_T,eghga_f_T,eghga_abz_T
         c_LAB_T,c_CAP_T,eKLE_T,eNEL_T,CO2_ffA_T,gene_g_t,gene_t,eff_ifg
         mapSIM, GDP_I,GDP_E,GDP_FP, CTOP_T
         ;
*         );
*$offtext
display eR, eKLE;
$INCLUDE intertemporal
                 );

         );