from __future__ import division
from numpy import sqrt
#To Determine the most economical current density
OFC=15 #Cost of single phase overhead feeder per km per unit area
AIDC= 10*OFC/100 #Annual Interest and depriciation charges
R=1/58 #Resistance of 1m length and 1 sq.mm
CLPUL= 2*R*1000 # Copper Loss per unit length per unit area per unit square current
ACL= CLPUL*365*24*0.5 # Annual Copper Loss
GC= 5/100 # Genereating Cost per unit.
ACEL= ACL*GC/1000 # Annual cost of energy loss
CD= sqrt(AIDC/ACEL) #Current Density
print 'The Most Economical Current Density for this Case is %0.2f A/sq.mm'%CD
#To determine the most economical cross sectional area
OFC=180 #Cost of single phase overhead feeder per km per unit area + 1200
AIDC= 10*OFC/100 #Annual Interest and depriciation charges + 120
R=1/58 #Resistance of 1m length and 1 sq.mm
I=200 # Maximum Current Flowing
CLPUL= 2*(I**2)*R*1000 # Copper Loss per unit length per unit area
ACL= CLPUL*365*24*(8/12) # Annual Copper Loss
GC= 5/100 # Genereating Cost per unit.
ACEL= ACL*GC/1000 # Annual cost of energy loss
A= sqrt(ACEL/AIDC) #Cross Sectional Area of the cable
print 'The Most Economical cross sectional area of the cable for this Case is %0.2f sq.mm'%A
# To determine the Most Economical Cross Sectional Area to supply a 3 Phase Load
LLF = lambda b:(0.25*b)+(0.75*(b**2)) # Function to determine the Loss Load Factor
OFC=0.20 #Cost of single phase overhead feeder per m per unit area + 10
AIDC= 10*OFC/100 #Annual Interest and depriciation charges + 1
TE= 2.5*(10**6) # Total energy to be supplied per annum
CEW=10/100 # Cost of energy wasted per unit
LFS= TE/(1000*365*24) # Load factor of supply
Llf=LLF(LFS) # Load Loss factor
R=1/58 # Resistance of the cable per unit length
PF=1 # Unity power factor
MD= 1*(10**6) # Maximum Demand
V=11*(10**3) # Voltage of the feeder
I=MD/(sqrt(3)*V*PF) # Full Load Current
FLCL= 3*(I**2)*R # Full Load Copper Loss per Metre
ACL= Llf*FLCL # Actual Copper Loss
CCL=ACL*(365*24*CEW/1000) # Cost of Copper Loss
A=sqrt(CCL/AIDC)
print 'The Most Economical Cross sectional area for this Case is %g A/sq.mm'%A
# Calculation Mistake in the Book. Hence according to the concepts in the book the answer is as calculated. Please Note.
#To determine the most economical cross section for a 3 Phase line 8 km long
OFC=5000 #Cost of feeder per km per unit area + 6250
AIDC= 8*OFC/100 #Annual Interest and depriciation charges + 625
AIDCPD= AIDC/365 # Annual Interest and depriciation charges per day
R=1/58 #Resistance of 1m length and 1 sq.mm
CEU= 5/100 # Cost of energy per unit
V=33*(10**3) # Voltage of the feeder
L10=3000*(10**3) # Load for 10 hrs at unity power factor.
L8=1000*(10**3) # Load for 8 hrs at unity power factor.
L6=2000*(10**3) # Load for 6 hrs at unity power factor.
LC = lambda b:b/(sqrt(3)*V) # Funtion to determine the Load Current
I10=LC(L10) # Load Current for 10 hrs at unity power factor.
I8=LC(L8) # Load Current for 8 hrs at unity power factor.
I6=LC(L6) # Load Current for 6 hrs at unity power factor.
ELPD=3*(((I10**2)*10)+((I8**2)*8)+((I6**2)*6))*1000/(100*58) # Energy Loss Over the DAy per km (Area in sq.cm)
CEL=ELPD*CEU/1000 # Cost of energy loass per km
A=sqrt(CEL/AIDCPD) # Cross Sectional Area
print 'The Most Economical Cross Sectional Area for a 3 Phase line 8 km Long is %0.4f sq..cm'%A
# To Calculate Voltage between middle wire and outer at each load point
Vs=220 #Supply Voltage
#Resistances of the respective Sections
Rab=0.2
Rbc=0.2
Rdh=0.2
Rhe=0.2
Ref=0.3
Rfl=0.1
Rlg=0.3
Rij=0.25
Rjk=0.2
Ra=0.2
Rd=0.4
Ri=0.3
#Currents following through the respective Section, Found using KCL
Iad=5
Ibe=10
Icl=12
Ihi=15
Ifj=5
Igk=15
Ia=Iad+Ibe+Icl # Current through the positive wire
Iab=Ia-Iad
Ibc=Ia-Iad-Ibe
Ii=Ihi+Ifj+Igk # Current through the negative wire
Id=Ii-Ia #Current through the Middle wire
Idh=Iad+Id
Ihe=Ihi-Idh
Ief=Ibe-Ihe
Ifl=Ief-Ifj
Ilg=Igk
Iij=Ii-Ihi
Ijk=Ii-Ihi-Ifj
#Voltage drops across each section
Va= Ra*Ia
Vab=Iab*Rab
Vbc=Ibc*Rbc
Vi=Ii*Ri
Vd=Id*Rd
Vdh=Idh*Rdh
Vhe=Ihe*Rhe
Vef=Ief*Ref
Vfl=Ifl*Rfl
Vlg=Ilg*Rlg
Vij=Iij*Rij
Vjk=Ijk*Rjk
#Voltage across the middle wire and the outer load points Using KVL
Vad=Vs-Va+Vd
Vbe=Vad-Vab-Vhe+Vdh
Vcl=Vbe-Vbc+Vfl+Vef
Vhi=Vs-Vd-Vdh-Vi
Vfj=Vhi-Vef+Vhe-Vij
Vgk=Vfj-Vfl-Vlg-Vjk
print 'The Voltages between middle wire and outer wire at each load point are:'
print ' 1. Vad = %g V'%(Vad)
print ' 2. Vbe = %g V'%(Vbe)
print ' 3. Vcl = %g V'%(Vcl)
print ' 4. Vhi = %g V'%(Vhi)
print ' 5. Vfj = %g V'%(Vfj)
print ' 6. Vgk = %g V'%(Vgk)
#To obtain the voltages at the far end of the positive and negative wires
Vs=220 # Supply Voltage
#Resistances of the given sections
r1=0.015
r2=0.035
r3=0.02
r4=0.01
r5=0.025
r6=0.015
r7=0.03
r8=0.01
r9=0.02
r10=0.03
r11=0.04
#Currents follwing in between the two wires
I1=10
I2=20
I3=25
I4=5
I5=15
I6=15
I7=15
I8=15
I9=18
I10=30
I11=15
#Effective resistances for the above currents
R1=r1
R2=r1+r2
R3=r1+r2+r3
R4=r1+r2+r3+r4
R5=r1+r2+r3+r4+r5
R6=r1+r2+r3+r4+r5+r6
R7=r7
R8=r7+r8
R9=r7+r8+r9
R10=r7+r8+r9+r10
R11=r7+r8+r9+r10+r11
#Voltage drop in the respective wires
#Outer positive wire
Vop=(I1*R1)+(I2*R2)+(I3*R3)+(I4*R4)+(I5*R5)+(I6*R6)
#Outer Negative wire
Vnp=(I7*R7)+(I8*R8)+(I9*R9)+(I10*R10)+(I11*R11)
# Net drop in neutral wire
Vn=Vop-Vnp
#Effective drop in the respective wires
#Outer positive wire
Veop=Vop+Vn
#Outer negative wire
Venp=Vnp-Vn
#Voltage at far end of the wires
#Positive Wire
Vpf=Vs-Veop
#Negaitve Wire
Vnf=Vs-Venp
print "The voltage drop at the far end of the wires are:"
print "Positive Wire: %g V"%(Vpf)
print "Negative Wire: %g V"%(Vnf)
#To determine the cross section of the conductor for a minimum consumer voltage
from sympy import symbols, solve
#Unknown Variable obtained in the equation
Ib, r=symbols("Ib r")
#Voltages at the respective ends
Va=230
Vb=230
#Minimum Consumers's Voltage
Vc=220
#Lenghths of the segments
r1=200
r2=200
r3=100
r4=300
r5=300
r6=100
#Effective length of segments
R1=r1
R2=r1+r2
R3=r1+r2+r3
R4=r1+r2+r3+r4
R5=r1+r2+r3+r4+r5
R6=r1+r2+r3+r4+r5+r6
#Current drawn by different loads
I1=25
I2=20
I3=25
I4=25
I5=10
#By the law of momemts of currents
Vnet=2*((I1*R1)+(I2*R2)+(I3*R3)+(I4*R4)+(I5*R5)-(Ib*R6)) #Since the equation is equated to zero, r vanishes
X=Vnet-(Va-Vb) #Polynomial Equation to Find Ib
Ib=float(solve(X, Ib)[0]) # Numerical Value of Ib
Ia=(I1+I2+I3+I4+I5)-Ib
#From the given figure it is clear that point M is the point of Minimum Potential
Vd=Va-Vc
X=(2*r*((Ia*r1)+((Ia-I1)*r2)+((Ia-I1-I2)*r3)))-Vd #Polynomial Equation to Find r
r=float(solve(X, r)[0]) # Resistance per unit length (Numerical Value)
Rstd=(1/58) # Resistance for 1m and 1 sq.mm
A=Rstd/r # Cross Section of Conductor required.
print 'The Cross Section of the Conductor to provide minimal consumer''s voltage is %0.f sq.mm '%A
#To determine the cross section of the conductor for a minimum consumer voltage
from sympy import symbols, solve
#Unknown Variable obtained in the equation
x,r = symbols("x r")
#Voltages at the respective ends
Va=235
Vb=230
#Minimum Consumers's Voltage
Vc=220
#Lenghths of the segments
r1=200
r2=200
r3=100
r4=300
r5=300
r6=100
#Effective length of segments
R1=r1
R2=r1+r2
R3=r1+r2+r3
R4=r1+r2+r3+r4
R5=r1+r2+r3+r4+r5
R6=r1+r2+r3+r4+r5+r6
#Current drawn by different loads
I1=25
I2=20
I3=25
I4=25
I5=10
#The Minimum Point assumed is N, Hence the current following to Point N is given by 'x'
#The Effective Drops
Van=Va-Vc
Vbn=Vb-Vc
V=Van/Vbn
A=2*(((x+I1)*r1)+(r2*x)) #Wrt to Van
B=2*(((I2-x)*r3)+((I2+I3-x)*r4)+((I2+I3+I4-x)*r5)+((I2+I3+I4+I5-x)*r6)) #Wrt to Vbn
C=A/B # The 'r' term gets eliminated
X=C-V # Polynomial Eqaution to find x
x=solve(X,x)[0]
Ia=I1+x #Current Supplied at end A
#It is clear the above assumed Minimal Point is wrong and it has to Point M
#Therefore finding drop at Point M
X=(2*r*(((x+I1)*r1)+(r2*x)+((x-20)*r3)))-Van #Polynomial Equation to find r
r=solve(X,r)[0] # Numerical Value of r
Rstd=(1/58) # Resistance for 1m and 1 sq.mm
A=Rstd/r # Cross Section of Conductor required.
print 'The Cross Section of the Conductor to provide minimal consumer''s voltage is %0.2f sq.mm'%A
# Please note the calculation mistake in the book. The value of r found out is wrong in the text book.
#To Determine the location and magnitude of minimum voltage
Vs=220 #Supply Voltage at End A and B
#Different Conductor Lenghts
#From End A
L1=100
L2=50
L3=50
L4=400 #Length of uniform loading
A=0.5 #Uniforming loading spread over 400m
r=0.05 # Resistance of Conductor per Km
#Different Currents drawn by various loads
I1=50
I2=75
I3=A*L4
#Taking moments of all currents at A
Ib=((I1*L1)+((L1+L2)*I2)+((L1+L2+L3+(L4/2))*I3))/(L1+L2+L3+L4)
Ia=I1+I2+I3-Ib
#Minimum Potential Point in this case is the point where All the current from B is drawn
X=Ib/A # Distance from B
Y=(L1+L2+L3+L4)-X #Distance from A
#Minimum Potential Drop
Vmind=(2*r/1000)*((Ia*L1)+((Ia-I1)*L2)+((Ia-I1-I2)*L3)+(((Y-L1-L2-L3)*A)*((Y-L1-L2-L3)/2)))
Vmin=Vs-Vmind #Minimum Potential
print 'The Location of The Minimum Voltage is %g m from side A and its magnitude is %0.1f V'%(Y,Vmin)
from __future__ import division
#To determine the currents supplied to the ring main from A and B
from numpy import array, linalg
#Currents in the ring scheme going clockwise
I1=50
I2=15
I3=25
I4=10
It = I1+I2+I3+I4 # Total Current
#Resistances of the respective segments going clockwise
R1=0.1
R2=0.12
R3=0.04
R4=0.08
R5=0.06
R6=0.02
Ra=0.04 #Resistance at A
Rb=0.06 #Resistance at B
# Va=Vb Net Voltage is Zero
#Dividing the total current from A as x and y
#Taking voltage drops in clockwise direction and anticlockwise directions
# We get two simultaneous equations
# 3.2x + y = 120
# x + 3y = 114
R=array([[3.2,1],[1,3] ])
V=array([[120],[114] ])
I=linalg.inv(R)*V
#To calculate and seperate the respective currents from the above matrix eqution
x=I[0]
y=I[1]
Ia=x+y
Ib=It-Ia
print 'a) When Va = Vb Ia = %g A and Ib =%g'%(Ia[0],Ib[0])
# Va+5=Vb Net Voltage is 5V
#Dividing the total current from A as x and y
#Taking voltage drops in clockwise direction and anticlockwise directions
# We get two simultaneous equations
# 3.2x + y = 120 -(5*NetVoltage) = 70
# x + 3y = 114 -(5*NetVoltage) = 64
R=array([[3.2,1],[1,3]] )
V=array([[70],[64] ])
I=linalg.inv(R)*V
#To calculate and seperate the respective currents from the above matrix eqution
x=I[0]
y=I[1]
Ia=x+y
Ib=It-Ia
print 'b) When Va + 5 = Vb Ia = %g A and Ib =%g'%(Ia[0],Ib[0])
#Highly Accurate Answers, Text Book answers are rounded off.
#To Compare the volume of copper required
from sympy import symbols, solve
#Unknown Resistances
r, r1, x=symbols('r r1 x')
#Lengths of the segements of the ring scheme
L1=100
L2=200
L3=200
L4=150
L5=150
#Currents taken by respective loads
I1=40
I2=20
I3=100
I4=40
It=I1+I2+I3+I4 #Total Current
#Without Interconnector
#Let x be the current flowing through the entire ring
Eq=(L1*x)+(L2*(x-I1))+(L3*(x-I1-I2))+(L4*(x-I1-I2-I3))+(L4*(x-I1-I2-I3-I4)) #Polynomial Equation to find x
x=float(solve(Eq, x)[0])
x1=It-x #Current flowing in the other direction
Vac1=((x1*L5)+((x1-I4)*L4))*r # Voltage across AC without the connector
MVac1=((x1*L5)+((x1-I4)*L4)) # Magnitude of Vac1
#With Interconnector
#Considering x amount of current to flow clockwise through segment AE
#Considering y amount of current to flow anticlockwise through segment AB
#Considering 200-(x+y) amount if current to flow through the segment AC
# Mesh Analysis of ABCDE gives 5x - 3y = 140
# Mesh Analysis of ABC gives 5x + 11y = 1120
R=array([[5,-3],[5,11]])
V=array([[140],[1120] ])
I=linalg.inv(R)*V
x=I[0]
y=I[1]
Vac2=(It-(x+y))*250*r1 # Voltage across AC with connector
MVac2=(It-(x+y))*250 # Magnitude of Vac2
print 'The Voltage drop across AC in both case is the same.'
print Vac2[0],'Is Equal to',Vac1
#To Compute the Numerical Values of the Ratio of resistances
RatioA = MVac1/MVac2
print RatioA[0],'is',r1,'divided by',r
#Effective Length of both the cases
Leff=L1+L2+L3+L4+L5
LeffC=Leff+250
#Volume is Length * Area
RatioV=Leff*RatioA/LeffC
print 'The Volume of copper without the connector is %g times the volume required with connector'%RatioV[0]
#To Determine the Voltage at the far end
from math import sin, acos
r=0.3 # Loop Resistance per km
x=0.15 # Loop Reactance per km
reactive = lambda y,z,v:(sin(acos(y)))*z*v*x/1000 #To find the reactive drop of the current
active = lambda y,z,v:y*z*v*r/1000 #To find the reactive drop of the current
#Power Factors of the loads from left to right
pf1=0.707
pf2=1
pf3=0.8
#Currents delivering the respective loads from left to right
I1=50
I2=60
I3=40
#Length of the conductors in metres
l1=200
l2=300
l3=300
#Effective length of the conductors in metres
L1=l1
L2=l1+l2
L3=l1+l2+l3
#Active component drops of the respective currents
Va1=active(pf1,I1,L1)
Va2=active(pf2,I2,L2)
Va3=active(pf3,I3,L3)
Vat=Va1+Va2+Va3 #Total Active Component Drop
#Reactive component drops of the respective currents
Vr1=reactive(pf1,I1,L1)
Vr2=reactive(pf2,I2,L2)
Vr3=reactive(pf3,I3,L3)
Vrt=Vr1+Vr2+Vr3 #Total Reactive Component Drop
Vt=Vrt+Vat # Total voltage drop
print 'The voltage drop at the far end is %g V '%(Vt)
#To determine line currents and neutral currents
from math import sqrt, acos, exp, atan, degrees, pi
from cmath import exp
V=400 # Supply Voltage
Vph=400/sqrt(3) #Phase Voltage
L=480*(10**3) #Balanced Load
#Loads at unit power factor
Lr=50*(10**3) #Load in R phase
Ly=150*(10**3) #Load in Y phase
Lb=200*(10**3) #Load in B phase
pf=0.8 #Power Factor lagging
theta=-1 * acos(pf) #Angle in radians of the balanced LAGGING current
I=L/(sqrt(3)*V*pf)*exp(1J*theta) #Balanced Current Magnitude
Ir=Lr/Vph #Magnitude of R phase voltage
Iy=Ly/Vph #Magnitude of Y phase voltage
Ib=Lb/Vph #Magnitude of B phase voltage
#Vr is taken as reference
#Angles of the phase voltages wrt to the reference
Thetay=-120
Thetab=120
Thetar=0
#Net Currents in the respective phases (RYB)
Irf=Ir+I
Iyf=Iy+I
Ibf=Ib+I
#Angle of the above currents in degrees
r=degrees(atan((Irf.imag)/(Irf.real)) )
y=degrees(atan((Iyf.imag)/(Iyf.real)) )
b=degrees(atan((Ibf.imag)/(Ibf.real)) )
#Angles of the above currents with respect to the reference
rf=Thetar+r
yf=Thetay+y
bf=Thetab+b
#Effective Currents wrt to the reference voltage
Irn=abs(Irf)*exp(1J*pi*rf/180)
Iyn=abs(Iyf)*exp(1J*pi*yf/180)
Ibn=abs(Ibf)*exp(1J*pi*bf/180)
In=Irn+Iyn+Ibn # Neutral Current
#Note Take Vr as reference
print 'The Net Current in phase R is %g/_%g A'%(abs(Irf),rf)
print 'The Net Current in phase Y is %g/_%g A'%(abs(Iyf),yf)
print 'The Net Current in phase B is %g/_%g A '%(abs(Ibf),bf)
print 'The Net Neutral Current is %g/_%g A\n'%(abs(In),degrees(atan((In.imag)/(In.real))))
#To determine the best location of the substation for a given set of loads
#Various loads and how they are positioned on the corners of a square of length 25km
L1=5000 #(0,25)
L2=8000 #(25,25)
L3=3000 #(25,0)
L4=6000 #(0,0)
L=25 #Length of the square
TL=L1+L2+L3+L4 # Total load
X=((L1*0)+(L2*L)+(L3*L)+(L4*0))/TL # X coordinate
Y=((L1*25)+(L2*L)+(L3*0)+(L4*0))/TL# Y coordinate
print 'The Susbstation must be located at (%g km,%0.2f km) '%(X,Y)