Ch-6 : Reactive Power Control

exa 6.1 page 201

In [6]:
from math import sqrt, atan, pi, degrees, cos
from scipy.linalg import expm
kV=220 #kV
Z=0.8+1J*0.2 #pu
V=1 #V(Voltage at load terminal)
X=0.2+0.05 #pu(line and transformer reactance)
P=(Z).real #pu
Q=(Z).imag #pu
BaseMVA=100 #MVA
BasekV=220 #kV
I=sqrt((P**2+Q**2)/V**2)*expm([[1J*atan(-(Z).imag/(Z).real)]])[0,0] #pu
Vb=V+I*(X*expm([[1J*pi/2]]))[0,0] #pu(Voltage at 200 kV bus)
fi_p=(atan((Vb).imag/(Vb).real)) #degree(power angle)
Vb=abs(Vb)*kV #kV(Voltage at 200 kV bus)
pf=cos(fi_p+(atan((Z).imag/(Z).real))) #power factor at 220 kV bus
print "Voltage at 220 kV bus = %0.2f kV "%Vb 
print "Power factor at 220 kV bus %0.4f lagging "%pf

Voltage at 220 kV bus = 235.15 kV 
Power factor at 220 kV bus 0.9076 lagging

exa 6.2 page 201

In [7]:
from math import sqrt, atan, pi, degrees, cos
from scipy.linalg import expm
kV=220 #kV
Z=0.8+1J*0.2 #pu
V=1 #V(Voltage at load terminal)
X=0.2+0.05 #pu(line and transformer reactance)
P=(Z).real #pu
Q=(Z).imag #pu
BaseMVA=100 #MVA
BasekV=220 #kV
I=sqrt((P**2+Q**2)/V**2) #pu
Vb=V+I*(X*expm([[1J*pi/2]])) #pu(Voltage at 200 kV bus)
fi_p=degrees(atan((Vb).imag/(Vb).real)) #degree(power angle)
Vb=abs(Vb)*kV #kV(Voltage at 200 kV bus)
pf=cos(fi_p*pi/180) #power factor at 220 kV bus
print "Voltage at 220 kV bus = %0.2f kV "%Vb 
print "Power factor at 220 kV bus = %0.4f lagging "%pf

Voltage at 220 kV bus = 224.63 kV 
Power factor at 220 kV bus = 0.9794 lagging

exa 6.3 page 202

In [8]:
from cmath import sinh, cosh, tan, sin
import cmath as cmt
l=350 #km(length of line)
Z=cmt.rect(180,75*pi/180) #ohm/phase(Total)
Y=cmt.rect(1*10**-3,90*pi/180) #Siemens/phase(Total)
z=Z/l #ohm/km
y=Y/l #Siemens/km
re=l*cmt.sqrt(z*y) #
Zc=cmt.sqrt(z/y) #ohm
print "Part(a) A,B,C,D parameters are : " 
A=cosh(re) #unitless
D=A #unitless
B=Zc*sinh(re) #ohm
C=sinh(re)/Zc #unitless
A_mag=abs(A) #unitless
A_angle=atan((A).imag/(A).real) #radian
B_mag=abs(B) #ohm
B_angle=atan((B).imag/(B).real) #radian
C_mag=abs(C) #unitless
C_angle=atan((C).imag/(C).real) #radian
C_angle=(degrees(C_angle)+180)*180/pi #degree(Converting -ve to +ve angle)
D_mag=abs(D) #unitless
D_angle=atan((D).imag/(D).real) #radian
print "Magnitude of A : %0.2f "% A_mag 
print "Angle of A = %0.2f degree  "%degrees(A_angle) 
print "Magnitude of B = %0.2f ohm"%B_mag 
print "Angle of B = %0.2f degree "%degrees(B_angle) 
print "Magnitude of C : %0.2f"%C_mag 
print "Angle of C = %0.2f degree  "%degrees(C_angle) 
print "Magnitude of D : %0.2f"%D_mag 
print "Angle of D = %0.2f degree" %degrees(D_angle)
#60% series compensation
B=B-1J*60/100*abs(Z)*sin(atan((Z).imag/(Z).real)) #ohm(considering series compensation=60%)
#For Equivalent pi-circuit
print "Part(b) A,B,C,D parameters of compensated line are : " 
Ydash=2/Zc*((cosh(re)-1)/sinh(re)) #S
A=1+B*Ydash/2 #unitless
D=A #unitless
C=2*Ydash/2+B*(Ydash/2)**2 #unitless
A_mag=abs(A) #unitless
A_angle=atan((A).imag/(A).real) #radian
B_mag=abs(B) #ohm
B_angle=atan((B).imag/(B).real) #radian
C_mag=abs(C) #unitless
C_angle=atan((C).imag/(C).real) #radian
C_angle=(degrees(C_angle)+180)*180/pi #radian(Converting -ve to +ve angle)
D_mag=abs(D) #unitless
D_angle=atan((D).imag/(D).real) #degree
print "Magnitude of A : %0.2f "% A_mag 
print "Angle of A = %0.2f degree  "%degrees(A_angle) 
print "Magnitude of B = %0.2f ohm"%B_mag 
print "Angle of B = %0.2f degree "%degrees(B_angle) 
print "Magnitude of C : %0.2f"%C_mag 
print "Angle of C = %0.2f degree  "%degrees(C_angle) 
print "Magnitude of D : %0.2f"%D_mag 
print "Angle of D = %0.2f degree" %degrees(D_angle)
#Answer for some parts are not accurate in the textbook.

Part(a) A,B,C,D parameters are : 
Magnitude of A : 0.91 
Angle of A = 1.42 degree  
Magnitude of B = 174.83 ohm
Angle of B = 75.45 degree 
Magnitude of C : 0.00
Angle of C = 296930.22 degree  
Magnitude of D : 0.91
Angle of D = 1.42 degree
Part(b) A,B,C,D parameters of compensated line are : 
Magnitude of A : 0.97 
Angle of A = 1.33 degree  
Magnitude of B = 78.37 ohm
Angle of B = 55.91 degree 
Magnitude of C : 0.00
Angle of C = 296850.37 degree  
Magnitude of D : 0.97
Angle of D = 1.33 degree

exa 6.4 page 202

In [9]:
from cmath import rect, sqrt,sinh, cosh, tan, sin
from numpy import conj
l=350 #km(length of line)
Z=rect(180,75*pi/180) #ohm/phase(Total)
Y=rect(1*10**-3,90*pi/180) #Siemens/phase(Total)
z=Z/l #ohm/km
y=Y/l #Siemens/km
re=l*sqrt(z*y) #
Zc=sqrt(z/y) #ohm
print "For Uncompensated Line, Constants are :" 
B=Z #ohm#B Parameter
A=1+Z*Y/2 #unitless#A Parameter
D=A #unitless#D Parameter
C=Y*(1+Z*Y/4) #S#C Parameter
A_mag=abs(A) #unitless
A_angle=atan((A).imag/(A).real) #radian
B_mag=abs(B) #ohm
B_angle=atan((B).imag/(B).real) #radian
C_mag=abs(C) #unitless
C_angle=atan((C).imag/(C).real) #radian
C_angle=(degrees(C_angle)+180)*180/pi #radian(Converting -ve to +ve angle)
D_mag=abs(D) #unitless
D_angle=atan((D).imag/(D).real) #degree
print "Magnitude and Angle of B = %0.2f ohm & %0.2f degree "%(B_mag,B_angle) 
print "Magnitude and Angle of A = %0.2f & %0.2f degree "%(A_mag,A_angle) 
print "Magnitude and Angle of D = %0.2f & %0.2f degree "%(D_mag,D_angle) 

print "Magnitude of C is %0.2f"%C_mag 
print "Angle of C = %0.2f degree"%C_angle
print "For Compensated Line, Constants are :" 
B=Z-0.6*1J*406 #ohm#B Parameter
A=1+conj(B)*Y/2 #unitless#A Parameter
D=A #unitless#D Parameter
C=Y*(1+Z*Y/4) #S#C Parameter
A_mag=abs(A) #unitless
A_angle=atan((A).imag/(A).real) #radian
B_mag=abs(B) #ohm
B_angle=atan((B).imag/(B).real) #radian
C_mag=abs(C) #unitless
C_angle=atan((C).imag/(C).real) #radian
C_angle=(degrees(C_angle)+180)*180/pi #radian(Converting -ve to +ve angle)
D_mag=abs(D) #unitless
D_angle=atan((D).imag/(D).real) #degree
print "Magnitude and Angle of B = %0.2f ohm & %0.2f degree "%(B_mag,B_angle) 
print "Magnitude and Angle of A = %0.2f & %0.2f degree "%(A_mag,A_angle) 
print "Magnitude and Angle of D = %0.2f & %0.2f degree "%(D_mag,D_angle) 
print "Magnitude of C is %0.2f"%C_mag 
print "Angle of C = %0.2f degree"%C_angle

For Uncompensated Line, Constants are :
Magnitude and Angle of B = 180.00 ohm & 1.31 degree 
Magnitude and Angle of A = 0.91 & 0.03 degree 
Magnitude and Angle of D = 0.91 & 0.03 degree 
Magnitude of C is 0.00
Angle of C = 5196.59 degree
For Compensated Line, Constants are :
Magnitude and Angle of B = 83.86 ohm & -0.98 degree 
Magnitude and Angle of A = 0.97 & 0.02 degree 
Magnitude and Angle of D = 0.97 & 0.02 degree 
Magnitude of C is 0.00
Angle of C = 5196.59 degree

exa 6.5 page 203

In [10]:
kv1=220 #kv
kv2=132 #kv
baseMVA=200 #MVA
#Base impedence in 132 kv circuit
baseZ2=kv2**2/baseMVA #ohm
z1=1J*75 #ohm
z2=1J*70 #ohm
z3=1J*90 #ohm
z1=z1/baseZ2 #pu
z2=z2/baseZ2 #pu
z3=z3/baseZ2 #pu
X_AD=1J*0.08+z1 #pu#Reactance from A to D
X_BD=1J*0.08+z2 #pu#Reactance from A to D
Zp=z3*X_AD*X_BD/(z3*X_AD+z3*X_BD+X_BD+X_AD) #parallel combination
sc_D=baseMVA/abs(Zp) #MVA#Short Circuit MVA at D
delQBYdelV=sc_D/kv2 #MVA/kv
delQ=delQBYdelV*4 #MVar
print "Var injection at Bus D = %0.2f MVar" %delQ
#Answer in the textbook is not accurate.

Var injection at Bus D = 18.48 MVar

exa 6.6 page 204

In [11]:
from cmath import rect, acos,sinh, cosh, tan, sin
A=rect(0.98,3*pi/180) #Constant
B=rect(110,75*pi/180) #ohm/phase
P=50 #MVA
pf=0.8 #lagging
V=132 #kV
#Formula : Pr=|Vs|*|Vr|/|B|*cosd(Beta-delta)-|A|*|Vr|**2/|B|*cosd(Beta-alfa) :
betaSUBdelta=acos((P*pf+abs(A)*V**2/abs(B)*cos(atan((B).imag/(B).real)-atan((A).imag/(A).real)))/V**2*abs(B)) 
Qr=V**2/abs(B)*sin(betaSUBdelta)-abs(A)*V**2/abs(B)*sin(atan((B).imag/(B).real)-atan((A).imag/(A).real)) #MVar
Qr=P*0.6-Qr #MVar#Since load require lagging component
print "(a) Capacity of shunt compensation equipment = %0.2f MVar   " %Qr.real
#part(b)
#Formula : Pr=|Vs|*|Vr|/|B|*cosd(Beta-delta)-|A|*|Vr|**2/|B|*cosd(Beta-alfa) :
P=0 #MW
betaSUBdelta=acos((P*pf+abs(A)*V**2/abs(B)*cos(atan((B).imag/B.real)-atan((A).imag/(A).real)))/V**2*abs(B)) 
Qr=V**2/abs(B)*sin(betaSUBdelta)-abs(A)*V**2/abs(B)*sin(atan((B).imag/(B).real)-atan((A).imag/(A).real)) #MVar
Qr=P*0.6-Qr #MVar#Since load require lagging component
print "(b) Capacity of shunt compensation equipment = %0.2f MVar "%-Qr.real 

(a) Capacity of shunt compensation equipment = 45.91 MVar   
(b) Capacity of shunt compensation equipment = 3.33 MVar

exa 6.7 page 206

In [12]:
from math import acos, sin, sqrt, acos
V=220 #kV
Z=20+1J*60 #ohm
Pr=100 #MVA
pf=0.8 #lagging pf
P=Pr*10**6*pf/3 #W
theta=acos(pf) #radian
Q=Pr*10**6*sin(theta)/3 #Vars
V1=V/sqrt(3)*1000 #V
V2=V1 #V
#ts**2*[1-(R*P+X*Q)/V1/V2]=V2/V1
ts=sqrt(V2/V1/(1-((Z).real*P+(Z).imag*Q)/V1/V2)) 
tr=1/ts 
print "Tap settings : ts = %0.2f "%ts.real
print "tr = %0.2f "%tr.real

Tap settings : ts = 1.06 
tr = 0.94

exa 6.8 page 207

In [13]:
from __future__ import division
from sympy import symbols, solve
kV1=132 #kV
kV2=33 #kV
kV3=11 #kV
MVA1=75 #MVA
MVA2=50 #MVA
MVA3=25 #MVA
X=0.12 #p.u.
#part(a)
P=60 #MW
V1=125 #kV
V1=V1/kV1 #p.u.
Q=MVA2/MVA1 #p.u.
#V1=Vn+X*Q/Vn
Vn=symbols('Vn') 
eqn=Vn**2-V1*Vn+X*Q
Vn=solve(eqn, Vn) #p.u.
Vn=Vn[0] #p.u.
Vn=Vn*kV1 #kV
k=Vn/kV2 #Transformer ratio
print "Under Load condition, transformer ratio = %0.3f "%k 
#part(b)
V1=140 #kV
V1=V1/kV1 #p.u.
Q=MVA3/MVA1 #p.u.
#V1=Vn+X*Q/Vn
Vn=symbols('Vn') 
eqn=Vn**2-V1*Vn+X*Q
Vn=solve(eqn, Vn) #p.u.
Vn=Vn[0] #p.u.
Vn=Vn*kV1 #kV
k=Vn/kV2 #Transformer ratio
print "Under No Load condition, transformer ratio = %0.3f" % k

Under Load condition, transformer ratio = 0.375 
Under No Load condition, transformer ratio = 0.157

exa 6.9 page 209

In [14]:
from math import tan, acos
V=132 #kV
Z=25+1J*66 #ohm
Pr=100 #MW
pf=0.9 #lagging pf
P=Pr*10**6/3 #W
theta=acos(pf) #radian
Q=Pr*10**6*tan(theta)/3 #vars
V1=V/sqrt(3)*1000 #V
V2=V1 #V
#ts**2*[1-(R*P+X*Q)/V1/V2]=V2/V1
ts=sqrt(V2/V1/(1-((Z).real*P+(Z).imag*Q)/V1/V2)) 
tr=1/ts 
print "Tap settings : ts = %0.2f "%ts.real 
print "tr = %0.2f "%tr.real

Tap settings : ts = 1.22 
tr = 0.82