Ch-10, Economic Operations of Steam Plants

example 10.2 Page 193

In [3]:
from numpy import arange
from sympy.mpmath import quad
mp=250     
def unit1(p1):
    ic=0.2*p1+30
    return ic
def unit2(p2):
    ic=0.15*p2+40
    return ic
mil=20
ttt=225
def un(ic):
    p1=(ic-30)/0.2
    p2=(ic-40)/0.15
    return [p1, p2]

for x in arange(40,61,5):
    [e,r]=un(x)
    if ttt==e+r:
        print "for the same incremental costs unit1 should supply %dMW and unit 2 shold supply %dMW,for equal sharing each unit should supply %3.1fMW"%(e,r,ttt/2)
        break
    

opo=ttt/2
u1=quad(unit1,[opo,e])
u2=quad(unit2,[r,opo])
uuu=(u1+u2)*8760
print "\nyearly extra cost is (%3.2f-%3.2f)8760 =%dper year"%(u1,u2,uuu)
print "\nthis if the load is equally shared by the two units an extra cost of Rs.%d will be incurred.in other words economic loading would result in saving of Rs.%dper year"%(uuu,uuu)

for the same incremental costs unit1 should supply 125MW and unit 2 shold supply 100MW,for equal sharing each unit should supply 112.0MW

yearly extra cost is (698.10-670.80)8760 =11991564per year

this if the load is equally shared by the two units an extra cost of Rs.11991564 will be incurred.in other words economic loading would result in saving of Rs.11991564per year

example 10.3 Page 198

In [4]:
def unit1(p1):
    ic=0.2*p1+30
    return ic
def unit2(p2):
    ic=0.15*p2+40
    return ic
tol=400
pd=50
u1c=5
u2c=1/0.15#from example10_1
p1pd=u1c/(u1c+u2c)
p2pd=u2c/(u1c+u2c)
pi=p1pd*pd
pt=p2pd*pd
print "p1=%1.5fMW\np2=%1.5fMW"%(pi,pt)
p11=pi+tol/2
p22=pt+tol/2
up1=unit1(p11)
up2=unit2(p22)
print "\nthe total load on 2 units would be %3.2fMW and %3.2fMW respectevily. it is easy to check that incremental cost will be same for two units at these loading.\n incremental cost of unit1 is %3.2fRs.MW,\n incremantal cost of unit 2 is %3.2fRs./MW"%(p11,p22,up1,up2)

p1=21.42857MW
p2=28.57143MW

the total load on 2 units would be 221.43MW and 228.57MW respectevily. it is easy to check that incremental cost will be same for two units at these loading.
 incremental cost of unit1 is 74.29Rs.MW,
 incremantal cost of unit 2 is 74.29Rs./MW

example10.5 Page 200

In [5]:
i1=0.8
i2=1.0
l1=complex(0.04,0.12)
l2=complex(0.03,0.1)
l3=complex(0.03,0.12)
vl=1

i3=i1+i2
v1=vl+i3*(l1)+i1*(l2)
v2=vl+i3*(l1)+i2*(l3)
p1=(i1*v1).real
p2=(i2*v2).real
cos1=(v1).real/abs(v1)
cos2=(v2).real/abs(v2)
b11=abs(((l1).real+(l2).real)/(v1**2*cos1**2))
b22=abs(((l1).real+(l3).real)/(v2**2*cos2**2))
b12=abs(((l1).real/(v1*v2*cos1*cos2)))
pl=(p1**2)*b11+(p2**2)*b22+2*p1*p2*b12
print "i1+i3=%dpu\nv1=%1.3f+%1.3fp.u\nv2=%1.3f+%1.3fp.u\np1=%1.3fp.u\np2=%1.3fp.u\ncos(ph1)=%1.3f\ncos(ph2)=%1.3f\nb11=%1.5fp.u\nb22=%1.5fp.u\nb12=%1.5fp.u\npl=%1.6fp.u"%(i3,v1.real,(v1).imag,v2.real,(v2).imag,p1,p2,cos1,cos2,b11,b22,b12,pl)

i1+i3=1pu
v1=1.096+0.296p.u
v2=1.102+0.336p.u
p1=0.877p.u
p2=1.102p.u
cos(ph1)=0.965
cos(ph2)=0.957
b11=0.05827p.u
b22=0.05764p.u
b12=0.03312p.u
pl=0.178800p.u

example10.7 Page 206

In [6]:
from math import atan, cos
za=complex(0.03,0.09)
zb=complex(0.1,0.3)
zc=complex(0.03,0.09)
zd=complex(0.04,0.12)
ze=complex(0.04,0.12)
ia=complex(1.5,-0.4)
ib=complex(0.5,-0.2)
ic=complex(1,-0.1)
id=complex(1,-0.2)
ie=complex(1.5,-0.3)
il1=.4
il2=.6
na1=1 ;nb1=0.6; nc1=0; nd1=.4; ne1=.6
na2=0 ;nb2=-0.4; nc2=1 ;nd2=.4; ne2=.6
vl=1
#some thing is messed
v1=vl+za*ia
v2=vl-zb*ib+zc*ic
a1=atan((ia).imag/(ia).real)
a2=atan((ic).imag/(ic).real)
cosa=cos(a1-a2)
cosph1=cos(atan((v1).imag/(v1).real)-a1)
cosph2=cos(atan((v2).imag/(v2).real)-a2)
b11=(na1**2*(za).real+nb1**2*(zb).real+nc1**2*(zc).real+nd1**2*(zd).real+ne1**2*(ze).real)/(abs(v1)**2*cosph1)
b22=(na2**2*(za).real+nb2**2*(zb).real+nc2**2*(zc).real+nd2**2*(zd).real+ne2**2*(ze).real)/((abs(v2)**2)*cosph2)
bb12=(abs(v1)*abs(v2)*cosph1*cosph2)
ab12=(na2*na1*(za).real+nb2*nb1*(zb).real+nc1*nc2*(zc).real+nd2*nd1*(zd).real+ne2*ne1*0.03)
b12=cosa*ab12/bb12
print "bus voltages at 2 buses are \nv1=%1.3f+i%1.3f,\nv2=%1.3f+i%1.3f"%((v1).real,(v1).imag,(v2).real,(v2).imag)
print "\nloss coffecients are \nb11=%1.5fp.u\nb22=%1.5fp.u\nb12=%1.5fp.u \n"%(b11,b22,b12)
print "loss coffecients in actual values is \nb11=%eM(W)-1\nb22=%eM(W)-1\nb12=%eM(W)-1\n"%(b11/100,b22/100,b12/100)

bus voltages at 2 buses are 
v1=1.081+i0.123,
v2=0.929+i-0.043

loss coffecients are 
b11=0.07877p.u
b22=0.07735p.u
b12=-0.00714p.u 

loss coffecients in actual values is 
b11=7.877236e-04M(W)-1
b22=7.734557e-04M(W)-1
b12=-7.136298e-05M(W)-1

example 10.8 Page 207

In [7]:
r1=22 ;r2=30 ;q1=0.2 ;q2=0.15
b22=0; b12=0 ;p1=100 ;pl=15#transmission losses are 0
b11=pl/(p1)**2
def power(x): #mathematical computation
    p1=(x-r1)/(q1+2*b11*x)
    p2=(x-r2)/q2
    return [p1, p2]
[a,b]=power(60)
print "l1=1/(1-%.3f*p1)\nl2=[1/(1-0)]=1\ngiven lamda=60\nsince ic1*l1=lamda  ic2*l2=lamda\ntotal load=%dMW"%(b11*2,a+b-(b11*a**2))

l1=1/(1-0.000*p1)
l2=[1/(1-0)]=1
given lamda=60
since ic1*l1=lamda  ic2*l2=lamda
total load=390MW

example 10.9 Page 208

In [8]:
from sympy.mpmath import quad
r1=22 ;r2=30 ;q1=0.2 ;q2=0.15
b22=0; b12=0 ;p1=100 ;pl=15#transmission losses are 0
b11=pl/(p1)**2
def power(x): #mathematical computation
    p1=(x-r1)/(q1+2*b11*x)
    p2=(x-r2)/q2
    return [p1,p2]
[a,b]=power(60)
pt=a+b-(b11*a**2)
z=quad(lambda u:q1*u+r1,[a,161.80])
y=quad(lambda v:q2*v+r2,[b,162.5])
m=z+y
print "net change in cost =Rs.%dper hour"%(m)
print "\nthus scheduling the generation by taking transmission losses into account would mean a saving of Rs.%dper hour in fuel cost"%(m)

net change in cost =Rs.-3757per hour

thus scheduling the generation by taking transmission losses into account would mean a saving of Rs.-3757per hour in fuel cost

example 10.10 Page 208

In [9]:
b11=0.001
b12=-0.0005
b22=0.0024
q1=0.08
r1=16
q2=0.08
r2=12
lamda=20

p2=0
for x in range(1,5):
    p1=(1-(r1/lamda)-(2*p2*b12))/((q1/lamda)+2*b11)
    p2=(1-(r2/lamda)-(2*p1*b12))/((q2/lamda)+2*b22)

pl=b11*p1**2+2*b12*p1*p2+b22*p2**2
pr=p1+p2-pl
print "p1=%2.1fMW,p2=%2.1fMW\npl=%1.1fMW\npower recevied %2.1fMW"%(p1,p2,pl,pr)

p1=189.2MW,p2=135.1MW
pl=54.1MW
power recevied 270.3MW

example 10.11 Page 209

In [10]:
a1=561 ;b1=7.92 ;c1=0.001562
a2=310 ;b2=7.85 ;c2=0.00194
ce=c1*c2/(c1+c2)
print "ce=%e"%(ce)
be=((b1/c1)+(b2/c2))*ce
print "be=%1.4f"%(be)
ae=a1-((b1**2)/4*c1)+a2-((b2**2)/4*c2)+((be**2)/4*ce)
print "ae=%3.3f \ncost characteristics of composite unit for demand pt is \nct=%3.3f+%1.4f*p1+%ep1**2"%(ae,ae,be,ce)

ce=8.652998e-04
be=7.8888
ae=870.959 
cost characteristics of composite unit for demand pt is 
ct=870.959+7.8888*p1+8.652998e-04p1**2

example 10.12 Page 210

In [11]:
a1=7700 ;b1=52.8 ;c1=5.5*10**-3
a2=2500; b2=15 ;c2=0.05#given eqution
plo=200 ;pup=800
ct=1000
l=[500,900,1200,500] ;t=[6, 16 ,20 ,24]#from given graph
def cost(y):
    p1=(2*c2*y-(b1-b2))/(2*(c1+c2))
    p2=y-p1
    return [p1,p2]
ma=max(l)
mi=min(l)
for x in range(0,3):
    [e ,g]=cost(l[x])
    if e<plo or g<plo or e>pup or g>pup:
        if e<plo or g<plo:
            v=min(e,g)
            [u for u, j in enumerate((e,g)) if j == v]
    if u==0:
        e=plo
        g=l(x)-e
    else:
        g=plo
        e=l[x]-g
    
    
print "\np1=%3.2fMW\tp2=%3.2fMW"%(e,g)

p1=1000.00MW	p2=200.00MW

example 10.13 Page 211

In [12]:
a1=2000 ;b1=20 ;c1=0.05; p1=350; p2=550
a2=2750 ;b2=26 ;c2=0.03091
def cost(a,b,c,p):
    co=a+b*p+c*p**2
    return co
print "(a)"
toco=cost(a1,b1,c1,p1)+cost(a2,b2,c2,p2)
print "total cost when each system supplies its own load Rs%.3f per hour"%(toco)
l=p1+p2
p11=(b2-b1+2*c2*l)/(2*(c1+c2))
p22=l-p11
totco=cost(a1,b1,c1,p11)+cost(a2,b2,c2,p22)
sav=toco-totco
tilo=p11-p1
print "(b)"
print "\n total cost when load is supplied in economic load dispatch method Rs%d per hour \n saving %.3f \n tie line load %.3f MW"%(totco,sav,tilo)

(a)
total cost when each system supplies its own load Rs41525.275 per hour
(b)

 total cost when load is supplied in economic load dispatch method Rs41447 per hour 
 saving 77.277 
 tie line load 30.905 MW

example10.14 Page 212

In [13]:
a1=5000 ;b1=450 ;c1=0.5 #for system 1 
e1=0.02 ;e2=-0.02#error
a1c=a1*(1-e1); b1c=b1*(1-e1) ;c1c=c1*(1-e1)
a2c=a1*(1-e2) ;b2c=b1*(1-e2) ;c2c=c1*(1-e2)
tl=200
def cost(a,b,c,p):
    co=a+b*p+c*p**2
    return co
p11=(b2c-b1c+2*c2c*tl)/(2*(c1c+c2c))
p22=tl-p11
totco=cost(a1c,b1c,c1c,p11)+cost(a2c,b2c,c2c,p22)
print "\npower at station 1 is %dMW \t power at station 2 is %dMW \n total cost on economic critieria method Rs%d per hour"%(p11,p22,totco)
tocoe=cost(a1c,b1c,c1c,tl/2)+cost(a2c,b2c,c2c,tl/2)
eop=tocoe-totco
print "\nextra operating cost due to erroneous scheduling Rs.%d per hour"%(eop)

power at station 1 is 111MW 	 power at station 2 is 89MW 
 total cost on economic critieria method Rs109879 per hour

extra operating cost due to erroneous scheduling Rs.121 per hour

example 10.16 Page 215

In [14]:
#given
ia=32 ;ib=32 ;ic=1.68; f=10**5
wt=18; rt=24-wt
p=30
def inpu(a,b,c,f,t,p):
    In=(a+b*p+c*p**2)*f*t
    return In
hi1=inpu(ia,ib,ic,f,wt,p); hi2=inpu(ia,ib,ic,f,rt,p/2)
print "(a)"
print "for full load condition for %d hours is %ekj \n for half load condition for%d s %ekj \n total load %ekj"%(wt,hi1,rt,hi2,hi1+hi2)
print "(b)"
te=p*wt+(p/2)*rt
ul=te/24
hin=inpu(ia,ib,ic,f,24,ul)
sav=hi1+hi2-hin
savp=sav/(te*1000)
print "\n total energy produced\t%dMW \n unifor load\t%dMW \n heat input under uniform load condition %ekj \n saving in heat energy %ekj \n saving in heat energy per kWh %dkj/kWh"%(te,ul,hin,sav,savp)

(a)
for full load condition for 18 hours is 4.507200e+09kj 
 for half load condition for6 s 5.340000e+08kj 
 total load 5.041200e+09kj
(b)

 total energy produced	630MW 
 unifor load	26MW 
 heat input under uniform load condition 4.799232e+09kj 
 saving in heat energy 2.419680e+08kj 
 saving in heat energy per kWh 384kj/kWh

example 10.17 Page 216

In [16]:
import numpy as np
#given
a1=450 ;b1=6.5 ;c1=0.0013
a2=300 ;b2=7.8 ;c2=0.0019
a3=80 ;b3=8.1 ;c3=0.005
tl=800#total load
ma = [600, 400, 200]
mi = [100, 50, 50]
d=np.mat([[1 ,1, 1] ,[2*c1, -2*c2, 0], [0, -2*c2, 2*c3]])
p=np.mat([[tl], [(b2-b1)], [(b2-b3)]]) 
pp=(d**-1)*p #matrix inversion method
print "\nloads on generaating station by economic creatirian method isp1=%fMW,p2=%fMW,p3=%fMW"%(pp[0],pp[1],pp[2])
for i in range(0,3):
    if pp[i]<mi[i]:
        pp[i]=mi[i]
    
    if pp[i]>mi[i]:
        pp[i]=ma[i]
    

pp[1] = tl-pp[0] - pp[2]
print "\nloads on generating station under critical conditions p1=%dMW p2=%dMW p3=%dMW"%(pp[0],pp[1],pp[2])

loads on generaating station by economic creatirian method isp1=669.734705MW,p2=116.134272MW,p3=14.131023MW

loads on generating station under critical conditions p1=600MW p2=150MW p3=50MW