Auto-slicing of shocks

Convergence of GAMS algorithms depends heavily on the distance between initial values and optimal solution. When the shock is large (there are numerous shocks in MIRAGE: labor supply, total factor productivity and capital are always changed between two periods and, depending on you scenario, you can add policy shocks on some variables like tariff rates, protection in services, taxes or subsidies), model's convergence is at risk, because the new equilibrium is far from the previous one. It is easy to reduce a large shock in a succession of small shocks by cutting the shock in small parts. This work can be tedious and is very contingent to the problem. The file NextSlice.gms takes automatically care of this issue.

A simple logic underlies it. It tries to solve the problem (by default the first iteration contains the whole shock). If it fails, it restarts from previous values but slices the shock in smaller parts. It will continue to divide the shock until it achieves the convergence or reaches the limit of a division in 1000 small shocks. The test of convergence is based on conopt solver status. Other solvers may require minor changes.

At each step, the exogeneous variables are incremented to reach their final values, according to the following equation:

 Var(...,Temps,sim)$(ord(Temps)=t1) = step*Varend(...,Temps,sim) + (1-step)*Var(...,Temps-1,sim);

If step equals 1 (usual for the first iteration), then the variable equals its target Varend. For other value of step between 0 and 1, the variable is set at a value between its previous value and its target. The smaller the step, the closer we are from the previous optimal solution. Hence, when a shock is large, the step will be small in order to stay close to the initial solution.

This function involves two specific gams files and some code in Simul.gms (Outdated). In MSD.gms (Outdated), Slicing.gms defines and initializes, for the baseline, all the parameter Varend, and the various parameters used after. NextSlice.gms substitutes for both NextBL.gms and Next.gms. It initializes and solves each year but cuts the shocks in small parts when they are too big. Finally some lines, at the end of Simul.gms (Outdated) initialize the Varend for the simulations.

Attention: If you add to the model other exogeneous variables than those already present in Slicing.gms, you should add the corresponding information in all the places described above. If not, the auto-slicing feature will not be activated for your new variables and will thus be of little interest. You can test the correct insertion of all exogeneous variables by solving the model with no increment. If step=0 then all variables are initialized to their previous year values, and the model must be at equilibrium. This test can be done by changing the fourteenth line in NEXTSlice.gms by:

    step  = step + (1/sum(idiv$(ord(idiv) = iter),div(idiv)))$(period>1);
 t1										= t1+1;
 t(Temps) 									= no;
 t(Temps) 									= yes$(ord(Temps)=t1);
 if(BL=1,
 	GDPVOLend(r,Temps,sim)$(ord(Temps)=t1) 		= GDPVOL.l(r,Temps-1,sim)*G_GDP(r,Temps);
 	Lbarend(r,Temps,Simul)$(ord(Temps)=t1) 		= Lbar.l(r,Temps-1,'ref')*G_P1565(r,Temps);
 	Hbarend(r,Temps,Simul)$(ord(Temps)=t1)		= Hbar.l(r,Temps-1,'ref')*G_P1565(r,Temps);
 ); 
 new_k(i,r,s)           			= sum((sim,Temps)$(ord(Temps)=t1),K.l(i,r,s,Temps-1,sim));
 if(t1=2,
   	old_k(i,r,s)         		= sum(Temps$(ord(Temps)=t1),K.l(i,r,s,Temps-1,'ref')-INV.l(i,r,s,Temps-1,'ref'))/(1-delta(r));
  	else
  		old_k(i,r,s)         	= K_(i,r,s);
 );
 maquette.solvestat 			= 0;
 maquette.modelstat 			= 0;
 iter               			= max(0,iter-2);
 while((maquette.solvestat ne 1) or (maquette.modelstat lt 15),
 $batinclude VARIABLES.gms  UPL
  LEON.l(Temps,sim)$(ord(Temps)=t1) 		= LEON.l(Temps-1,sim);
  step                               		= 0;
  period                            		= 0;
  maquette.solvestat               	 	= 1;
  maquette.modelstat                		= 16;
  error                            		= 0;
  iter                              		= iter+1;
  while(((1-step)>1E-15) and (error ne 3) and ((period gt 1) or ((maquette.solvestat=1) and (maquette.modelstat ge 15))),
 $include LB.gms
    period									= period+1;
    step									= step+(1/sum(idiv$(ord(idiv)=iter),div(idiv)))$(period>=1);
    put log1;
    loop(Temps$(ord(Temps)=t1), put 'Year ', Temps.tl:4:0 /);
    put 'Iteration ', period:2:0, '/', sum(idiv$(ord(idiv)=iter),div(idiv)):4:0 //;
    putclose;
    if(BL=1,GDPVOL.fx(r,Temps,sim)$(ord(Temps)=t1) 	= step*GDPVOLend(r,Temps,sim)    	+(1-step)*GDPVOL.l(r,Temps-1,sim);
 	else PGF.fx(r,Temps,sim)$(ord(Temps)=t1)  	= step*PGFend(r,Temps,sim)     		+(1-step)*PGF.l(r,Temps-1,sim);
    );
    Lbar.fx(r,Temps,sim)$(ord(Temps)=t1)   	= step*Lbarend(r,Temps,sim)      	+(1-step)*Lbar.l(r,Temps-1,sim);
    Hbar.fx(r,Temps,sim)$(ord(Temps)=t1)   	= step*Hbarend(r,Temps,sim)      	+(1-step)*Hbar.l(r,Temps-1,sim);
    K_(i,r,s)                              	= step*new_k(i,r,s)              	+(1-step)*old_k(i,r,s);
        taxP(i,r,Temps,sim)$(ord(Temps)=t1)    	= step*taxPend(i,r,Temps,sim)    	+(1-step)*TAXP(i,r,Temps-1,sim);
    tsubK(i,r,Temps,sim)$(ord(Temps)=t1)   	= step*tsubKend(i,r,Temps,sim)   	+(1-step)*tsubK(i,r,Temps-1,sim);
    tsubTE(i,r,Temps,sim)$(ord(Temps)=t1)  	= step*tsubTEend(i,r,Temps,sim)  	+(1-step)*tsubTE(i,r,Temps-1,sim);
    DD(i,r,s,Temps,sim)$(ord(Temps)=t1)    	= step*DDend(i,r,s,Temps,sim)    	+(1-step)*DD(i,r,s,Temps-1,sim);
    rente(i,r,s,Temps,sim)$(ord(Temps)=t1) 	= step*renteend(i,r,s,Temps,sim) 	+(1-step)*rente(i,r,s,Temps-1,sim);
    taxAMF(i,r,s,Temps,sim)$(ord(Temps)=t1)	= step*taxAMFend(i,r,s,Temps,sim)	+(1-step)*taxAMF(i,r,s,Temps-1,sim);
    solve maquette using cns;
    error 					= (error+1)$((maquette.solvestat ne 1) or (maquette.modelstat lt 15));););

If the model still does not converge, you should think about applying your shocks on successive years, because some policy shocks are too important to pass in only one year (remember that capital is fixed in the short run in dynamic mode, so adjustments take times to make).