French Forest Sector Model (FFSM++)

Laboratoire d'Economie Forestière - LEF - Nancy, France

User Tools

Site Tools


en:dev:initialisation_mkt

This is an old revision of the document!


This is the workflow of the model in the first few years.

p → products; pp → primary products; tp → transformed products; r → region; t → time

$pl_{r,tp,t=2}= \sum_{pp} {pl_{r,pp,t=2} * a_{r,pp,t=2,tp} } + (m_{r,tp,t=2} * pol\_trSub_{r,tp,t=2}) $ $dl_{r,pp,t=2}= \sum_{pt} {sl_{r,tp,t=2} * a_{r,pp,t=2,tp} } $ $st_{r,p,t=2} = sl_{r,p,t=2} r,p,t=2} + sa_{r,p,t=2} r,p,t=2} $ $dt_{r,p,t=2} = dl_{r,p,t=2} r,p,t=2} + da_{r,p,t=2} r,p,t=2} $

ModelCoreSpatial::initMarketModule(){

msgOut(MSG_INFO, "Starting market module (init stage)..");
for(uint i=0;i<regIds2.size();i++){
  int r2 = regIds2[i];

RPAR('st',i,prd,t-1) = RPAR('sl',i,prd,t-1)+RPAR('sa',i,prd,t-1); RPAR('dt',i,prd,t-1) = RPAR('dl',i,prd,t-1)+RPAR('da',i,prd,t-1);

  for (uint ap=0;ap<allProducts.size();ap++){
    //double debug = gpd("dl",r2,allProducts[ap],secondYear);
    double stvalue =   gpd("sl",r2,allProducts[ap],secondYear)
                     + gpd("sa",r2,allProducts[ap],secondYear);
    double dtvalue =   gpd("dl",r2,allProducts[ap],secondYear)
                     + gpd("da",r2,allProducts[ap],secondYear);
    spd(stvalue,"st",r2,allProducts[ap],secondYear,true);
    spd(stvalue,"stFromHarvesting",r2,allProducts[ap],secondYear,true);
    spd(dtvalue,"dt",r2,allProducts[ap],secondYear,true);
  }
  // q1(i,p_tr)  =  1/(1+((RPAR('dl',i,p_tr,t-1)/RPAR('da',i,p_tr,t-1))**(1/psi(i,p_tr)))*(RPAR('pl',i,p_tr,t-1)/PT(p_tr,t-1)));
  // p1(i,p_tr)              =  1-q1(i,p_tr);
  // RPAR('dc',i,p_tr,t-1)   =  (q1(i,p_tr)*RPAR('da',i,p_tr,t-1)**((psi(i,p_tr)-1)/psi(i,p_tr))+ p1(i,p_tr)*RPAR('dl',i,p_tr,t-1)**((psi(i,p_tr)-1)/psi(i,p_tr)))**(psi(i,p_tr)/(psi(i,p_tr)-1));
  // RPAR('pc',i,p_tr,t-1)   =  (RPAR('da',i,p_tr,t-1)/RPAR('dc',i,p_tr,t-1))*PT(p_tr,t-1)+(RPAR('dl',i,p_tr,t-1)/RPAR('dc',i,p_tr,t-1))*RPAR('pl',i,p_tr,t-1);
  // RPAR('pc',i,p_pr,t-1)   =  (RPAR('sa',i,p_pr,t-1)/RPAR('sc',i,p_pr,t-1))*PT(p_pr,t-1)+(RPAR('sl',i,p_pr,t-1)/RPAR('sc',i,p_pr,t-1))*RPAR('pl',i,p_pr,t-1);
  // RPAR('pw',i,p_tr,t-1)   =  (RPAR('dl',i,p_tr,t-1)*RPAR('pl',i,p_tr,t-1)+RPAR('da',i,p_tr,t-1)*PT(p_tr,t-1))/RPAR('dt',i,p_tr,t-1) ; //changed 20120419
  // K(i,p_tr,t-1)           =  k1(i,p_tr)*RPAR('sl',i,p_tr,t-1);
  for(uint sp=0;sp<secProducts.size();sp++){
    double psi = gpd("psi",r2,secProducts[sp],secondYear);
    double dl  = gpd("dl",r2,secProducts[sp],secondYear);
    double da  = gpd("da",r2,secProducts[sp],secondYear);
    double pl  = gpd("pl",r2,secProducts[sp],secondYear);
    double sl  = gpd("sl",r2,secProducts[sp],secondYear);
    double k1  = gpd("k1",r2,secProducts[sp],secondYear);
    double pWo = gpd("pl",WL2,secProducts[sp],secondYear); // World price (local price for region 99999)
    double q1  = 1/ ( 1+pow(dl/da,1/psi)*(pl/pWo) );
    double p1  = 1-q1;
    double dc  = pow(
            q1*pow(da,(psi-1)/psi) + p1*pow(dl,(psi-1)/psi)
           ,
            psi/(psi-1)
           );
    double pc = (da/dc)*pWo
           +(dl/dc)*pl;
    double pw = (dl*pl+da*pWo)/(dl+da);
    double k = k1*sl;
    spd(q1,"q1",r2,secProducts[sp],firstYear,true);
    //spd(p1,"p1",r2,secProducts[sp],firstYear,true);
    spd(dc,"dc",r2,secProducts[sp],secondYear,true);
    spd(pc,"pc",r2,secProducts[sp],secondYear,true);
    spd(pw,"pw",r2,secProducts[sp],secondYear,true);
    spd(k,"k",r2,secProducts[sp],secondYear,true);
  }
  // t1(i,p_pr)              =  1/(1+((RPAR('sl',i,p_pr,t-1)/RPAR('sa',i,p_pr,t-1))**(1/eta(i,p_pr)))*(RPAR('pl',i,p_pr,t-1)/PT(p_pr,t-1)));
  // r1(i,p_pr)              =  1-t1(i,p_pr);
  // RPAR('sc',i,p_pr,t-1)   =  (t1(i,p_pr)*RPAR('sa',i,p_pr,t-1)**((eta(i,p_pr)-1)/eta(i,p_pr))+ r1(i,p_pr)*RPAR('sl',i,p_pr,t-1)**((eta(i,p_pr)-1)/eta(i,p_pr)))**(eta(i,p_pr)/(eta(i,p_pr)-1))
  // RPAR('pc',i,p_pr,t-1)   =  (RPAR('sa',i,p_pr,t-1)/RPAR('sc',i,p_pr,t-1))*PT(p_pr,t-1)+(RPAR('sl',i,p_pr,t-1)/RPAR('sc',i,p_pr,t-1))*RPAR('pl',i,p_pr,t-1);
  // RPAR('pw',i,p_pr,t-1)   =  (RPAR('sl',i,p_pr,t-1)*RPAR('pl',i,p_pr,t-1)+RPAR('sa',i,p_pr,t-1)*PT(p_pr,t-1))/RPAR('st',i,p_pr,t-1) ; //changed 20120419
  for(uint pp=0;pp<priProducts.size();pp++){
    double sl  = gpd("sl",r2,priProducts[pp],secondYear);
    double sa  = gpd("sa",r2,priProducts[pp],secondYear);
    double eta = gpd("eta",r2,priProducts[pp],secondYear);
    double pl  = gpd("pl",r2,priProducts[pp],secondYear);
    double pWo = gpd("pl",WL2,priProducts[pp],secondYear); // World price (local price for region 99999)
    double t1 = 1/ ( 1+(pow(sl/sa,1/eta))*(pl/pWo) );
    double r1 = 1-t1;
    double sc = pow(
             t1*pow(sa,(eta-1)/eta) + r1*pow(sl,(eta-1)/eta)
          ,
             eta/(eta-1)
          );
    double pc = (sa/sc)*pWo+(sl/sc)*pl;
    double pw = (sl*pl+sa*pWo)/(sl+sa);
    spd(t1,"t1",r2,priProducts[pp],firstYear,true);
   //spd(r1,"r1",r2,priProducts[pp],firstYear,true);
    spd(sc,"sc",r2,priProducts[pp],secondYear,true); 
    spd(pc,"pc",r2,priProducts[pp],secondYear,true);
    spd(pw,"pw",r2,priProducts[pp],secondYear,true);
  }
  // up to here tested with gams output on 20120628, that's fine !!
} // end for each region in level 2
// initializing the exports to zero quantities
// initializing of the transport cost for the same region to one and distance to zero
for(uint r1=0;r1<l2r.size();r1++){
  for(uint r2=0;r2<l2r[r1].size();r2++){
    for(uint p=0;p<allProducts.size();p++){
      for(uint r2To=0;r2To<l2r[r1].size();r2To++){
        spd(0,"rt",l2r[r1][r2],allProducts[p],secondYear,true,i2s(l2r[r1][r2To]));  // regional trade, it was exp in gams
        if(l2r[r1][r2] == l2r[r1][r2To]){
          spd(1,"ct",l2r[r1][r2],allProducts[p],firstYear,true,i2s(l2r[r1][r2To])); // as long this value is higher than zero, rt within the same region is not choosen by the solver, so the value doesn't really matters. If it is zero, the solver still works and results are the same, but reported rt within the region are crazy high (100000)
        }
      }
    } // end each product
    for(uint r2To=0;r2To<l2r[r1].size();r2To++){
      if(l2r[r1][r2] == l2r[r1][r2To]){
        spd(0,"dist",l2r[r1][r2],"",firstYear,true,i2s(l2r[r1][r2To])); // setting distance zero in code, so no need to put it in the data
      }
    }
  } // end of r2 regions
} // end of r1 region

}

en/dev/initialisation_mkt.1490087764.txt.gz · Last modified: 2018/06/18 16:44 (external edit)