Division by zero error

Problems with modeling
Post Reply
m_tanawat
User
User
Posts: 13
Joined: 6 years ago

Division by zero error

Post by m_tanawat »

Dear admin,

Could you please check my program?
It showed division by zero error.
I try solving it. But I cannot.

I used to study other posts.
You suggested to assign syntax $xx ,<0>
For my model, can i do it too?
If so, how can i do?


I am sorry I didn't provide input file since this platform could not attach file.
I believe that you can check only formula you can solve it.

Best regards,
TItech


$Ontext
Sample CGE for MPSGE
2 goods, 2 energy, 2 final demand and 2 factors model
$Offtext


set

sec sector and goods /agr,cha,coa,pet,min,foo,tex,pap,pec,per,opp,nom,mem,ret,woo,ele,gas,wat,ctr,who,rai,pub,roa,waf,air,sil,pos,ser,ipa/
i goods /agr,cha,coa,pet,min,foo,tex,pap,pec,per,opp,nom,mem,ret,woo,ele,gas,wat,ctr,who,rai,pub,roa,waf,air,sil,pos,ser,ipa/
en(i) energy goods / cha,coa,pet,min ,per ,ele,gas/
ne(i) non energy goods /agr, foo,tex,pap,pec ,opp,nom,mem,ret,woo, wat,ctr,who,rai,pub,roa,waf,air,sil,pos,ser,ipa/
qm(ne) equipment / mem /
j(sec) sector /agr,cha,coa,pet,min,foo,tex,pap,pec,per,opp,nom,mem,ret,woo,ele,gas,wat,ctr,who,rai,pub,roa,waf,air,sil,pos,ser,ipa/
ej(j) energy sector /cha,coa,pet,min,per,ele,gas/
nj(j) non energy sector /agr,foo,tex,pap,pec,opp,nom,mem,ret,woo,wat,ctr,who,rai,pub,roa,waf,air,sil,pos,ser,ipa/
fd final demand /hou,gov,inv,stc,exp,imp,tx_m/
tec technology /tec0*tec2/
sv household energy service /sv1/
;
* 23 sectors and 23 goods ( no seperation of commodity and activitiy)
* 9 energy goods, 7 final demands

$ontext
set
sec sector and goods /2nd,eqm,3rd,en1,en2/
i(sec) goods /2nd,eqm,3rd,en1,en2/
en(i) energy goods / en1,en2/
ne(i) non energy goods /2nd,eqm,3rd /
qm(ne) equipment / eqm /
j(sec) sector /2nd,eqm,3rd,en1,en2/
ej(j) energy sector / en1,en2/
nj(j) non energy sector /2nd,eqm,3rd /
fd final demand /h_c, h_i, g_c, g_i, stc, exp, imp, tx_m/
tec technology /tec0*tec2/
sv household energy service /sv1/
;
$offtext

alias (jj,j) ;
alias (ii,i) ;
alias (en,en2) ;
alias (ne,ne2) ;
alias (tec,tc2) ;

parameter
el_lke(j) elasticity of substitution between labor and energy-capital
el_ke(j) elasticity of substitution between energy and capital
el_e(j) elasticity of substitution among energies
el_ex(i) elasticity of transformation
el_im(i) elasticity of substitution

kle0(j) total value added in sector j
ke0(j) total energy caital composite in sector j
ene0(j) total energy composite in sector j
e0(en,j) energy input with CO2 in sector j
k0(j) capital input in sector j
l0(j) labor input in sector j
x0(i,j) intermediate input in sector j
em0(en,j) energy input as material

q0(j) total produced goods
s0(i) total domestic supply
d0(i) total supply
o0(i) total domestic output

c0(i) final consumption by household
i0(i) investment(fixed capital formation) by household
cg0(i) final consumption by government
ig0(i) investment(fixed capital formation) by government
sc0(i) stock change
ex0(i) export
im0(i) import
utl0 total household consumption
ivh0 total household investment
gvc0 total government consumption
gvi0 total government investment

eh0(en) energy with CO2 by household
wc0(i) waste generation from household consumption
wi0(i) waste generation from household investment

mh0 household income
mg0 government income
g_h transfer between household and government

cmb_r(en,j) combustion rate of energy

fs foreign saving
pintx(i) international price of exported goods
pintm(i) international price of imported goods

ef_co2(en) CO2 emission factor
efa_co2(j) CO2 emission factor
CO2e0(en,j) energy related CO2 emission from sector j
CO2a0(j) activity related CO2 emission from sector j
CO2f0(en) CO2 emission from household

end_l labor endowment
end_k(j,tec) capital endowment
end_kn capital endowment
end_CO2 total endowment

tx_m(i) import tax
tx_o(i) output tax
tx_k(j) capital tax
tx_l(j) labor tax

ql0(j) quantity of labor
qk0(j) quantity of capital
qim0(i) quantity of import
qd0(i) quantity of domestic demand
qo0(i)

invh household investment
invg government investment
sc(i) stock change
wc(i) waste from household consumption
wi(i) waste from household investment

acnt(j) activity constraint
q_up(j) upple limit of sector

ad0(qm,sv,tec)

depr depreciation rate
dep_h depreciation rate (household energy equipment)

c_tax carbon tax policy / 0 /
* 0: no-policy
* 1: introduction
;


*$ontext
parameter your_use(*,*) ;
$libinclude xlimport your_use thailandiodynamic_transport.xlsx sheet1!A1:AL34


parameter your_make(*,*) ;
$libinclude xlimport your_make thailandiodynamic_transport.xlsx sheet4!A1:AE31

;

$ontext
parameter your_tax(*,*) ;
$libinclude xlimport your_tax thailandiodynamic.xlsx sheet3!A1:X4
$offtext

parameter use(*,*) use matrix ;
parameter out(*,*) make matrix ;

$ontext
parameter tax (*,*) tax matrix ;
$offtext

* ---
use(i,j) = your_use(i,j) ;
use("cap",j) = your_use("cap",j) ;
use("lab",j) = your_use("lab",j) ;
use(i,fd) = your_use(i,fd) ;

out(j,i) = your_make(i,j) ;

* ---
*$offtext

ef_co2("cha") = 0 ;
ef_co2("coa") = 0.84;
ef_co2("pet") = 0.2 ;
ef_co2("per") = 0.01 ;
ef_co2("ele") = 0.01 ;
efa_co2(j) = 0 ;
efa_co2("coa") = 0.12

* primary energy conversion method to get emission factor (monetary unit)
*Data derive from IO table and energy balance table
;

$ontext
ef_co2("en1") = 30 ;
ef_co2("en2") = 3 ;
efa_co2(j) = 0 ;
efa_co2("2nd") = 0.01 ;
$offtext

$ontext
table use(*,*) use matrix
2nd eqm 3rd en1 en2 h_c h_i g_c g_i stc exp imp tx_m tot
2nd 110 30 40 8 3 60 10 10 4 5 30 -10 -3 297
eqm 19 4 9 2 2 20 28 10 6 5 10 -5 -2 108
3rd 50 23 137 2 3 148 150 40 40 -20 100 -40 -10 623
en1 7 1 4 1 2 12 0 0 0 2 0 -1 -1 27
en2 1 1 9 2 1 2 0 1 0 1 0 0 0 18
cap 30 20 130 4 2
lab 38 21 194 2 2
tot 255 100 523 21 15
;

table out(*,*) make matrix
2nd eqm 3rd en1 en2 tot
2nd 240 13 2 0 0 255
eqm 5 79 16 0 0 100
3rd 0 0 523 0 0 523
en1 1 1 0 19 0 21
en2 0 0 0 2 13 15
h_c 43 0 40 4 4 91
h_i 3 0 2 0 0 5
tot 292 93 583 25 17
;
$offtext

table tax(*,*) tax matrix
agr cha coa pet min foo tex pap pec per opp nom mem ret woo ele gas wat ctr who tra ser ipa
tx_o -1804710 2030 24293 47530797 1351847 201219215 8471879 2944483 12116089 165285232 4286503 1052702 132096369 14161674 1916579 10800716 12277633 1328589 10060544 77539508 30777620 129832403 4180882
tx_k 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
tx_l 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
;

tax("tx_o",i) = your_use("tx_o",i) ;


el_ke(j) = 0.8 ;
el_lke(j) = 0.5 ;
el_e(j) = 1 ;
el_ex(i) = 2 ;
el_im(i) = 2 ;

depr = 0.05 ;
dep_h = 0.1 ;

cmb_r(en,j) = 1 ;
cmb_r("coa","coa") = 1 ;

x0(i,j) = use(i,j) ;
c0(i) = use(i,"hou") ;
i0(i) = use(i,"inv") ;
cg0(i) = use(i,"gov") ;
ig0(i) = 0 ;
sc0(i) = use(i,"stc") ;
im0(i) = -use(i,"imp")- use(i,"tx_m") ;
ex0(i) = use(i,"exp") ;
k0(j) = use("cap",j) ;
l0(j) = use("lab",j) ;
q0(j) = sum(i, x0(i,j)) + k0(j) + l0(j) ;
s0(i) = sum(j, x0(i,j)) + c0(i) + i0(i) + cg0(i) + ig0(i) + sc0(i) ;
wc0(i) = 0 ;
wi0(i) = 0 ;
o0(i) = sum(j, out(j,i)) + wc0(i) + wi0(i) ;
d0(i) = o0(i) + tax("tx_o",i) - ex0(i) ;

*utl0 = sum(i, c0(i)) ;
utl0 = sum(ne$(not qm(ne)), c0(ne)) ;
ivh0 = sum(i, i0(i)) ;
gvc0 = sum(i, cg0(i)) ;
gvi0 = sum(i, ig0(i)) ;


tx_m(i)$(im0(i)) = -use(i,"tx_m")/(-use(i,"imp")) ;
tx_o(i) = tax("tx_o",i)/(o0(i)) ;
tx_k(j) = tax("tx_k",j)/(k0(j) - tax("tx_k",j)) ;
tx_l(j) = tax("tx_l",j)/(l0(j) - tax("tx_l",j)) ;

ql0(j) = l0(j) /(1+tx_l(j)) ;
qk0(j) = k0(j) /(1+tx_k(j)) ;
qim0(i) = im0(i)/(1+tx_m(i)) ;
qo0(i) = o0(i) /(1+tx_o(i)) ;

g_h = sum(i, tx_m(i)*qim0(i))
+ sum(i, tx_o(i)*o0(i))
+ sum(j, tx_k(j)*qk0(j) + tx_l(j)*ql0(j))
- sum(i, cg0(i) + ig0(i)) ;


e0(en,j) = x0(en,j)*cmb_r(en,j) ;
ene0(j) = sum(en, e0(en,j)) ;
ke0(j) = k0(j) + ene0(j) ;
kle0(j) = l0(j) + ke0(j) ;

em0(en,j) = x0(en,j)*(1-cmb_r(en,j)) ;

eh0(en) = c0(en) ;
invh = ivh0 ;
invg = gvi0 ;
sc(i) = sc0(i) ;
wc(i) = wc0(i) ;
wi(i) = wi0(i) ;

CO2e0(en,j) = ef_co2(en)*e0(en,j) ;
CO2a0(j) = efa_co2(j)*q0(j) ;
CO2f0(en) = ef_co2(en)*c0(en) ;

end_k(j,tec) = 0 ;
end_k(j,"tec0") = qk0(j) ;
end_kn = 0 ;
end_l = sum(j, ql0(j)) ;
end_CO2 = sum(en, sum(j, CO2e0(en,j)) + CO2f0(en)) + sum(j, CO2a0(j)) ;

fs = sum(i, qim0(i) - ex0(i)) ;
pintx(i) = 1 ;
pintm(i) = 1 ;

mh0 = sum(j, k0(j)*(1-tx_k(j)) + l0(j)*(1-tx_l(j))) + fs + sum(i, sc0(i) + wc0(i) + wi0(i)) + g_h ;
mg0 = sum(i, tax("tx_o",i) - use(i,"tx_m")) + sum(j, tax("tx_k",j) + tax("tx_l",j)) - g_h ;


acnt(j) = 0 ;
q_up(j) = q0(j) ;

ad0(qm,sv,tec) = 0 ;


parameter
grs_out(j,i,tec) change of output to activity
grs_ne(ne,j,tec) change of non-energy inputs productivity
grs_co2a(j,tec) change of co2 emissions from activity
grs_l(j,tec) change of labor productivity
grs_k(j,tec) change of capital productivity
grs_en(en,j,tec) change of energy productivity
grs_co2e(en,j,tec) change of co2 emissions from energy

grn_out(j,i,tec) change of output to activity
grn_ne(ne,j,tec) change of non-energy inputs productivity
grn_co2a(j,tec) change of co2 emissions from activity
grn_l(j,tec) change of labor productivity
grn_k(j,tec) change of capital productivity
grn_en(en,j,tec) change of energy productivity
grn_co2e(en,j,tec) change of co2 emissions from energy

avt(*,tec) available technoloigy
* 1: available
* 0: not available
ad_inv(ne,j,tec) additional investment

gr_c(ne) change of final consumption
gr_ehn(sv,en,tec) change of energy consumption
gr_ehs(sv,en,tec) change of energy consumption
;

grs_out(j,i,tec) = 1 ;
grs_ne(ne,j,tec) = 1 ;
grs_co2a(j,tec) = 1 ;
grs_l(j,tec) = 1 ;
grs_k(j,tec) = 1 ;
grs_en(en,j,tec) = 1 ;
grs_co2e(en,j,tec) = 1 ;

grn_out(j,i,"tec0") = 1 ;
grn_ne(ne,j,"tec0") = 1 ;
grn_co2a(j,"tec0") = 1 ;
grn_l(j,"tec0") = 0.9 ;
grn_k(j,"tec0") = 0.9 ;
grn_en(en,j,"tec0") = 0.9 ;
grn_co2e(en,j,"tec0") = 1 ;

grn_out(j,i,"tec1") = 1 ;
grn_ne(ne,j,"tec1") = 1 ;
grn_co2a(j,"tec1") = 1 ;
grn_l(j,"tec1") = 0.9 ;
grn_k(j,"tec1") = 0.9 ;
grn_en(en,j,"tec1") = 0.8 ;
grn_co2e(en,j,"tec1") = 1 ;

grn_out(j,i,"tec2") = 1 ;
grn_ne(ne,j,"tec2") = 1 ;
grn_co2a(j,"tec2") = 1 ;
grn_l(j,"tec2") = 0.9 ;
grn_k(j,"tec2") = 0.9 ;
grn_en(en,j,"tec2") = 0.7 ;
grn_co2e(en,j,"tec2") = 1 ;


gr_c(ne) = 1 ;
gr_ehn(sv,en,"tec0") = 1 ;
gr_ehn(sv,en,"tec1") = 0.9 ;
gr_ehn(sv,en,"tec2") = 0.8 ;
gr_ehs(sv,en,"tec0") = 1.1 ;
gr_ehs(sv,en,"tec1") = 0.9 ;
gr_ehs(sv,en,"tec2") = 0.8 ;

avt(j,tec) = 0 ;
avt(j,"tec0") = 1 ;
avt("hus",tec) = 0 ;
avt("hus","tec0") = 1 ;

ad_inv(ne,j,tec) = 0 ;
ad_inv(ne,j,"tec1") = qk0(j)*0.01 ;
ad_inv(ne,j,"tec2") = qk0(j)*0.02 ;


parameter
es(sv) energy service demand
es0(sv) energy service demand
seh0(en,sv) energy demand by energy service demand

cqs(qm,sv,tec)
cqs0(qm,sv,tec)
cqn0(qm,sv,tec)

sehn0(en,sv)
sehs0(en,sv)

esn0(sv,tec)
ess0(sv,tec)
;


cqn0(qm,sv,tec) = c0(qm) ;
cqs0(qm,sv,tec) = c0(qm)/dep_h ;
cqs(qm,sv,tec) = 0 ;
cqs(qm,sv,"tec0") = cqs0(qm,sv,"tec0") ;

sehn0(en,sv) = eh0(en)*dep_h /gr_ehn(sv,en,"tec0") ;
sehs0(en,sv) = eh0(en)*(1-dep_h)/gr_ehs(sv,en,"tec0") ;

esn0(sv,tec) = sum(en, sehn0(en,sv)*gr_ehn(sv,en,tec)) + sum(qm, cqn0(qm,sv,tec) + ad0(qm,sv,tec)) ;
ess0(sv,tec) = sum(en, sehs0(en,sv)*gr_ehs(sv,en,tec)) + sum(qm, cqs0(qm,sv,tec)) ;

es(sv) = esn0(sv,"tec0") + ess0(sv,"tec0") ;
es0(sv) = es(sv) ;
seh0(en,sv) = eh0(en) ;



$ontext
$MODEL:CGE4_mpsge
*$FUNLOG:.TRUE
$PEPS:0

$SECTORS:
ACTS(j,tec)$(q0(j) and end_k(j,tec)) ! activity
ACTN(j,tec)$(q0(j) and end_kn and avt(j,tec)) ! activity
COH$utl0 ! household consumption
ESHN(sv,tec)$(avt("hus",tec))
ESHS(sv,tec)$(sum(qm, cqs(qm,sv,tec)) and avt("hus",tec))
COH_E(en)$eh0(en) ! household energy
IVH$ivh0 ! household investment
COG$gvc0 ! government consumption
IVG$gvi0 ! government investment

export(i) ! distribution between export and domestic
import(i) ! composition of domestic and import
pr_ex(i)$ex0(i) ! export
pr_im(i)$im0(i) ! import

$COMMODITIES:
PO(i) ! price of activity in sector j
PES(sv)
PEH(en)$eh0(en) !
PZ(qm,sv,tec)$(avt("hus",tec) and cqs(qm,sv,tec))

PK(j,tec)$end_k(j,tec) ! rent (price of existing capital)
PKN$end_kn ! rent (price of new capital)
PL ! wage (price of labor)
PCO2 ! price of CO2
PACT(j)$(acnt(j)) ! upper limit

PU$utl0 ! price of aggregated household consumption
PIVH$ivh0 ! price of aggregated household investment
PGC$gvc0 ! price of aggregated government consumption
PGI$gvi0 ! price of aggregated government investment

PD(i) ! price of produced domstic goods i
PS(i) ! price of supplied goods i
PEX(i)$(ex0(i)) ! price of exported goods i
PIM(i)$(im0(i)) ! price of imported goods i
EXR ! exchange rate


$CONSUMERS:
MH ! income
GOV ! government


$AUXILIARY:
CO2_PRI$(c_tax=1) ! carbon policy


$PROD:ACTS(j,tec)$(q0(j) and end_k(j,tec)) t:0 s:0 en.tl:0
o:PO(i) q:(out(j,i)*grs_out(j,i,tec)) p:1
i:PS(ne) q:(x0(ne,j)*grs_ne(ne,j,tec)) p:1
i:PCO2 q:(CO2a0(j)*grs_co2a(j,tec)) p:0
i:PS(en) q:(em0(en,j)) p:1
i:PACT(j)$acnt(j) q:q0(j) q:0
i:PL q:(ql0(j)*grs_l(j,tec)) p:(1+tx_l(j)) A:GOV T:tx_l(j)
i:PK(j,tec) q:(qk0(j)*grs_k(j,tec)) p:(1+tx_k(j)) A:GOV T:tx_k(j)
i:PS(en) q:(e0(en,j)) p:1 en.tl:
i:PCO2#(en) q:(co2e0(en,j)*grs_co2e(en,j,tec)) p:0 en.tl:


$PROD:ACTN(j,tec)$(q0(j) and end_kn and avt(j,tec)) t:0 s:0 en.tl:0
o:PO(i) q:(out(j,i)*grn_out(j,i,tec)) p:1
i:PS(ne) q:(x0(ne,j)*grn_ne(ne,j,tec)) p:1
i:PCO2 q:(CO2a0(j)*grn_co2a(j,tec)) p:0
i:PS(en) q:(em0(en,j)) p:1
i:PACT(j)$acnt(j) q:q0(j) q:0
i:PL q:(ql0(j)*grn_l(j,tec)) p:(1+tx_l(j)) A:GOV T:tx_l(j)
i:PKN q:(qk0(j)*grn_k(j,tec)) p:(1+tx_k(j)) A:GOV T:tx_k(j)
i:PS(ne) q:(ad_inv(ne,j,tec)) p:1
i:PS(en) q:(e0(en,j)) p:1 en.tl:
i:PCO2#(en) q:(co2e0(en,j)*grn_co2e(en,j,tec)) p:0 en.tl:


$PROD:export(i) t:el_ex(i)
o:PEX(i) q:ex0(i) p:1
o:PD(i) q:d0(i) p:1
i:PO(i) q:o0(i) p:(1+tx_o(i)) A:GOV T:tx_o(i)

$PROD:import(i) s:el_im(i)
o:PS(i) q:s0(i) P:1
i:PIM(i) q:qim0(i) p:(1+tx_m(i)) A:GOV T:tx_m(i)
i:PD(i) q:d0(i) p:1

$PROD:pr_ex(i)$ex0(i)
o:EXR q:(ex0(i)*pintx(i))
i:PEX(i) q:(ex0(i))

$PROD:pr_im(i)$im0(i)
o:PIM(i) q:(qim0(i))
i:EXR q:(qim0(i)*pintm(i))


$prod:COH$utl0 s:1 b:0
O:PU$utl0 q:utl0 p:1
I:PS(ne)$(not qm(ne)) q:(c0(ne)*gr_c(ne)) p:1

$prod:ESHN(sv,tec)$avt("hus",tec) s:0
O:PES(sv) q:esn0(sv,tec) p:1
I:PEH(en) q:(sehn0(en,sv)*gr_ehn(sv,en,tec)) p:1
I:PS(qm) q:(cqn0(qm,sv,tec)) p:1
I:PS(qm) q:(ad0(qm,sv,tec)) p:1

$prod:ESHS(sv,tec)$(sum(qm, cqs(qm,sv,tec)) and avt("hus",tec)) s:0
O:PES(sv) q:ess0(sv,tec) p:1
I:PEH(en) q:(sehs0(en,sv)*gr_ehs(sv,en,tec)) p:1
I:PZ(qm,sv,tec) q:(cqs0(qm,sv,tec)) p:1

$prod:COH_E(en)$eh0(en) s:0
O:PEH(en) q:eh0(en) p:1
I:PS(en) q:c0(en) p:1
I:PCO2 q:co2f0(en) p:0


$prod:IVH$ivh0 s:0
O:PIVH$ivh0 q:ivh0 p:1
I:PS(ne) q:i0(ne) p:1

$prod:COG$gvc0 s:0
O:PGC$gvc0 q:gvc0 p:1
I:PS(i) q:cg0(i) p:1

$prod:IVG$gvi0 s:0
O:PGI$gvi0 q:gvi0 p:1
I:PS(i) q:ig0(i) p:1


$demand:MH
D:PU q:utl0
E:PES(sv) q:(-es(sv))
E:PIVH q:(-invh)
E:PK(j,tec)$end_k(j,tec) q:end_k(j,tec)
E:PKN$(end_kn) q:end_kn
E:PL q:end_l
E:EXR q:(fs)
E:PS(i) q:(-sc(i))
E:PO(i) q:wc(i)
E:PO(i) q:wi(i)
E:PL q:g_h
E:PACT(j)$acnt(j) q:q_up(j)
E:PZ(qm,sv,tec) q:cqs(qm,sv,tec)

$demand:GOV s:0
D:PGC q:gvc0
E:PGI q:(-invg)
E:PCO2$(c_tax=0) q:(end_co2)
E:PCO2$(c_tax=1) q:(end_co2) r:CO2_PRI
E:PL q:(-g_h)


$report:
V:QKS(j,tec)$(end_k(j,tec) and ke0(j)) i:PK(j,tec) PROD:ACTS(j,tec)
V:QLS(j,tec)$(end_k(j,tec) and kle0(j)) i:PL PROD:ACTS(j,tec)
V:QXS(ne,j,tec)$(end_k(j,tec) and ene0(j)) i:PS(ne) PROD:ACTS(j,tec)
V:QES(en,j,tec)$(end_k(j,tec) and e0(en,j)) i:PS(en) PROD:ACTS(j,tec)
V:QKN(j,tec)$(end_kn and avt(j,tec)) i:PKN PROD:ACTN(j,tec)
V:QLN(j,tec)$(end_kn and avt(j,tec)) i:PL PROD:ACTN(j,tec)
V:QXN(ne,j,tec)$(end_kn and avt(j,tec)) i:PS(ne) PROD:ACTN(j,tec)
V:QEN(en,j,tec)$(end_kn and avt(j,tec)) i:PS(en) PROD:ACTN(j,tec)
V:QCN(ne) i:PS(ne) PROD:COH
V:QCE(en) i:PS(en) PROD:COH_E(en)


$constraint:CO2_PRI$(c_tax=1)
end_co2*CO2_PRI =e= sum((j,tec)$end_k(j,tec), ACTS(j,tec)*CO2a0(j)*grs_co2a(j,tec) + sum(en, ACTS(j,tec)*co2e0(en,j)*grs_co2e(en,j,tec)) )
+ sum((j,tec)$avt(j,tec), ACTN(j,tec)*CO2a0(j)*grn_co2a(j,tec) + sum(en, ACTN(j,tec)*co2e0(en,j)*grn_co2e(en,j,tec)) )
+ sum(en, COH_E(en)*co2f0(en)) ;


$offtext

$SYSINCLUDE MPSGESET CGE4_mpsge


ACTS.l(j,tec)$(avt(j,tec) = 0) = 0 ;
ACTN.l(j,tec) = 0 ;

ESHN.l(sv,tec)$(avt("hus",tec) = 0) = 0 ;
ESHS.l(sv,tec)$(avt("hus",tec) = 0) = 0 ;


PL.fx = 1 ;
PCO2.l = 0 ;
PACT.l(j) = 0 ;


option iterlim = 0 ;
OPTION SYSOUT = ON;
OPTION SOLPRINT = ON ;

$INCLUDE CGE4_mpsge.GEN
SOLVE CGE4_mpsge USING MCP;

CGE4_mpsge.OPTFILE = 1 ;

*$exit

set
t time perid /2010*2031/
;

parameter
yr(t) year
cap_stk_t total capital stock
cap_stk(j,tec) capital stock
cap_cnv
inv_t total investment

output
outputK
outputY

;


table scen_main(t,*) future scenario
gdp pop ghg import export
2010 10802 63878 345 6768 7136
2011 11300 64076 359 7005 7386
2012 12349 64456 373 7250 7644
2013 12901 64758 388 7504 7912
2014 13132 65124 404 7766 8189
2015 13533 65729 420 8038 8475
2016 14034 65931 428 8320 8772
2017 14525 66261 437 8611 9079
2018 15034 66592 445 8912 9397
2019 15560 66925 454 9224 9726
2020 16104 67260 463 9547 10066
2021 16910 67327 473 10024 10569
2022 17755 67394 482 10525 11098
2023 18643 67462 492 11052 11653
2024 19575 67529 502 11604 12235
2025 20554 67597 512 12185 12847
2026 21581 67664 522 12794 13489
2027 22660 67732 532 13433 14164
2028 23793 67800 543 14105 14872
2029 24983 67867 554 14810 15616
2030 26232 67935 555 15551 16397
2031 27544 68003 559 16328 17216
;

parameter grow(t,*) growth rate ;

grow(t,"GDP")$(ord(t) lt card(t)) = (scen_main(t+1,"GDP")/scen_main(t,"GDP")-1)*100 ;
grow(t,"POP")$(ord(t) lt card(t)) = (scen_main(t+1,"POP")/scen_main(t,"POP")-1)*100 ;
grow(t,"GHG")$(ord(t) lt card(t)) = (scen_main(t+1,"GHG")/scen_main(t,"GHG")-1)*100 ;


yr(t) = 2009 + ord(t) ;

inv_t = sum(i, i0(i)+ig0(i)) ;
cap_stk_t = inv_t/(grow("2010","gdp")/100 + depr) ;
cap_stk(j,tec) = cap_stk_t*end_k(j,tec)/sum((jj,tc2), end_k(jj,tc2)) ;
cap_cnv = sum((j,tec), end_k(j,tec))/cap_stk_t ;



loop(t$(yr(t) le 2030),
*loop(t$(yr(t) le 2013),


option iterlim = 100000 ;

$INCLUDE CGE4_mpsge.GEN
SOLVE CGE4_mpsge USING MCP;


output(t,"C_H") = sum(ne, QCN.l(ne)) + sum(en, QCE.l(en)) ;
output(t,"C_G") = sum(i, COG.l*cg0(i)) ;
output(t,"I_H") = sum(ne, IVH.l*i0(ne)) ;
output(t,"I_G") = sum(i, IVG.l*ig0(i)) ;
output(t,"EXP") = sum(i, pr_ex.l(i)*ex0(i)) ;
output(t,"IMP") = sum(i, pr_im.l(i)*qim0(i)) ;
output(t,"TXM") = sum(i, pr_im.l(i)*qim0(i)*tx_m(i)) ;
output(t,"STC") = sum(i, sc(i)) ;
output(t,"GDP_cal") = output(t,"C_H") + output(t,"C_G") + output(t,"I_H") + output(t,"I_G")
+ output(t,"EXP") - output(t,"IMP") - output(t,"TXM") + output(t,"STC") ;
output(t,"GDP_exp")$(yr(t) = 2010) = output(t,"GDP_cal") ;
output(t,"GDP_exp")$(yr(t) > 2010) = output(t-1,"GDP_exp")*(1+grow(t,"GDP")/100) ;

output(t,"CO2") = sum((j,tec), ACTS.l(j,tec)*CO2a0(j)*grs_co2a(j,tec) + sum(en, ACTS.l(j,tec)*co2e0(en,j)*grs_co2e(en,j,tec))
+ ACTN.l(j,tec)*CO2a0(j)*grn_co2a(j,tec) + sum(en, ACTN.l(j,tec)*co2e0(en,j)*grn_co2e(en,j,tec)))
+ sum(en, COH_E.l(en)*co2f0(en)) ;
output(t,"PCO2") = PCO2.l ;
output(t,"cap_st") = cap_stk_t ;

outputK(j,"S",tec,t) = QKS.l(j,tec) ;
outputK(j,"N",tec,t) = QKN.l(j,tec) ;

outputY(j,"S",tec,t) = ACTS.l(j,tec)*sum(i, out(j,i)*grs_out(j,i,tec)) ;
outputY(j,"N",tec,t) = ACTN.l(j,tec)*sum(i, out(j,i)*grn_out(j,i,tec)) ;
outputY(j,"tot",tec,t) = outputY(j,"S",tec,t) + outputY(j,"N",tec,t) ;
outputY("tot","S",tec,t) = sum(j, outputY(j,"S",tec,t)) ;
outputY("tot","N",tec,t) = sum(j, outputY(j,"N",tec,t)) ;
outputY("tot","tot",tec,t) = sum(j, outputY(j,"tot",tec,t)) ;



* update of efficiency
grs_out(j,i,tec)$(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) = (grs_out(j,i,tec) *outputK(j,"S",tec,t) + grn_out(j,i,tec) *outputK(j,"N",tec,t))
/(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) ;
grs_ne(ne,j,tec)$(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) = (grs_ne(ne,j,tec) *outputK(j,"S",tec,t) + grn_ne(ne,j,tec) *outputK(j,"N",tec,t))
/(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) ;
grs_co2a(j,tec)$(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) = (grs_co2a(j,tec) *outputK(j,"S",tec,t) + grn_co2a(j,tec) *outputK(j,"N",tec,t))
/(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) ;
grs_l(j,tec)$(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) = (grs_l(j,tec) *outputK(j,"S",tec,t) + grn_l(j,tec) *outputK(j,"N",tec,t))
/(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) ;
grs_k(j,tec)$(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) = (grs_k(j,tec) *outputK(j,"S",tec,t) + grn_k(j,tec) *outputK(j,"N",tec,t))
/(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) ;
grs_en(en,j,tec)$(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) = (grs_en(en,j,tec) *outputK(j,"S",tec,t) + grn_en(en,j,tec) *outputK(j,"N",tec,t))
/(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) ;
grs_co2e(en,j,tec)$(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) = (grs_co2e(en,j,tec)*outputK(j,"S",tec,t) + grn_co2e(en,j,tec)*outputK(j,"N",tec,t))
/(outputK(j,"S",tec,t) + outputK(j,"N",tec,t)) ;

gr_ehs(sv,en,tec)$(ESHS.l(sv,tec)*sum(qm, cqs0(qm,sv,tec)) + ESHN.l(sv,tec)*sum(qm, cqn0(qm,sv,tec)))
= (gr_ehs(sv,en,tec)*ESHS.l(sv,tec)*sum(qm, cqs0(qm,sv,tec)) + gr_ehn(sv,en,tec)*ESHN.l(sv,tec)*sum(qm, cqn0(qm,sv,tec)))
/(ESHS.l(sv,tec)*sum(qm, cqs0(qm,sv,tec)) + ESHN.l(sv,tec)*sum(qm, cqn0(qm,sv,tec))) ;

grn_out(j,i,tec) = grn_out(j,i,tec) *0.98 ;
grn_ne(ne,j,tec) = grn_ne(ne,j,tec) *1 ;
grn_co2a(j,tec) = grn_co2a(j,tec) *1 ;
grn_l(j,tec) = grn_l(j,tec) *(1- grow(t,"GDP")/100) ;
grn_k(j,tec) = grn_k(j,tec) *(1- grow(t,"GDP")/100) ;
grn_en(en,j,tec) = grn_en(en,j,tec) *0.98 ;
grn_co2e(en,j,tec) = grn_co2e(en,j,tec)*1 ;

gr_ehn(sv,en,tec) = gr_ehn(sv,en,tec)*1 ;
*gr_c("3rd") = gr_c("3rd") *1.01 ;


* update of endowment
end_l = end_l*(1+grow(t,"pop")/100) ;
end_co2 = end_co2*1.1 ;
*end_co2 = end_co2*(1+grow(t,"GHG")/100) ;

if((yr(t) ge 2020),
*c_tax = 1 ;
*PCO2.fx = 0.001 ;
) ;

end_k(j,tec) = (outputK(j,"S",tec,t) + outputK(j,"N",tec,t))*(1-depr) ;
end_kn = inv_t*cap_cnv ;

cqs(qm,sv,tec) = (ESHS.l(sv,tec)*cqs0(qm,sv,tec) + ESHN.l(sv,tec)*cqn0(qm,sv,tec))*(1-dep_h) ;

cap_stk_t = (sum((j,tec), end_k(j,tec)) + end_kn)/(cap_cnv) ;
inv_t = cap_stk_t*(grow(t+1,"gdp")/100 + depr) ;
invh = inv_t*ivh0/(ivh0+gvi0) ;
invg = inv_t*gvi0/(ivh0+gvi0) ;

fs = sum(i, qim0(i)*scen_main(t+1,"import")/scen_main("2010","import")
- ex0(i) *scen_main(t+1,"export")/scen_main("2010","export")) ;

sc(i)$(sc(i) < 0) = 0 ;
sc(i)$(sc(i) > 0) = -sc(i) ;

*avt(nj,"tec1")$(yr(t) ge 2020) = 1 ;
*avt(nj,"tec2")$(yr(t) ge 2025) = 1 ;

q_up(j) = q_up(j)*(1+grow(t,"gdp")/100) ;


ACTN.up(j,tec) = sum(tc2, ACTS.l(j,tc2)+ACTN.l(j,tc2))*0.2 ;


);

display outputY, output ;
User avatar
Renger
Posts: 639
Joined: 7 years ago

Re: Division by zero error

Post by Renger »

Hi

You can add files: see below the window where you write your question. There is an "Attachment"-Tab.
Furthermore, if you do post code, put it between code tags (the </> symbol will insert those tags on top of the writing window), so the format is saved.

Cheers
Renger

PS. The administrator is usually not the one who answers questions. This is a forum, where GAMS-users help other GAMS-users
____________________________________
Enjoy modeling even more: Read my blog on modeling at The lazy economist
m_tanawat
User
User
Posts: 13
Joined: 6 years ago

Re: Division by zero error

Post by m_tanawat »

Hi,

Thank you very much for your suggestion.
I send my attached file here.

1. gams file
2. input file

There is only division by zero error remaining in the program.
If you can help solve it, it will run normal completion.

If would be appreciated if you can advise how you solve it.

Best regards,
Titech
Attachments
thailandiodynamic_transport.xlsx
(144.87 KiB) Downloaded 301 times
sample_4_Thai_dynamic 12Mar.gms
(27.41 KiB) Downloaded 287 times
User avatar
Renger
Posts: 639
Joined: 7 years ago

Re: Division by zero error

Post by Renger »

Hi
Gams points you to the line 326:

tx_o(i) = tax("tx_o",i)/(o0(i)) ;

It therefore looks like you have a o0(i) that is for one or more i zero.
I added the following code to find the zero values

o0(i) = o0(i) + EPS;
display i, o0;

(I added EPS, as GAMS doesn't report 0. EPS is a very small value).
And I got

Code: Select all

----    327 PARAMETER o0  total domestic output

agr 1.783453E+9,    cha 5709282.000,    coa 1.739714E+7,    pet 3.392892E+8
min 7.067952E+7,    foo 2.180018E+9,    tex 7.632852E+8,    pap 2.532224E+8
pec 1.566716E+9,    per 1.045041E+9,    opp 1.110240E+8,    nom 3.985297E+8
mem 5.743578E+9,    ret 1.036455E+9,    woo 1.973150E+8,    ele 6.949323E+8
gas 3.907113E+8,    wat 5.717528E+7,    ctr 9.049838E+8,    who 2.859890E+9
rai         EPS,    pub 2.981214E+8,    roa 4.161443E+8,    waf 2.222584E+8
air 2.629094E+8,    sil 1.262934E+8,    pos 3.039277E+8,    ser 4.375924E+9
ipa 2.104718E+8
It looks like there is no domestic output for the sector "rai". You should check this, and if this is correct, you can rewrite your tax line as follows:

Code: Select all

tx_o(i)$o0(i)          = tax("tx_o",i)/(o0(i)) ;
Gams would skip this line (setting tx_o to zero) if there is no domestic output (don't forget to delete the the line with EPS).

Further advise:
  • Your data is badly scaled. A rule of thumb is to scale your values around 1 (in your case it looks like if you divide everything by 1E9, this would be a good choice.
    Why should you scale? If you scale the solver might more easily solve as you reduce the precision. The model is solved if there is no difference smaller than 1E7. If you don't scale this might mean that you want to solve your model to the precision of a very small money amount (e.g. 0.0001 Baht) which is certainly smaller than the precision of your data.

    It looks like you adjust your energy technology but forgot to calibrate it properly with an additional loop.
Cheers

Renger
____________________________________
Enjoy modeling even more: Read my blog on modeling at The lazy economist
Post Reply