The load flow problem - Washington State University

1 22 July 20111 LOAD FLOWA. J. Conejo Univ. Castilla La Mancha201122 July 20112 The load flow problem0. caseMatrixGeneral equationsBus classificationVariable typesand limitsBUSY22 July 20113 Theload flow problem3. The gauss -Seidel solution techniqueIntroductionAlgorithm initializationPQ BusesPV BusesStopping criterion22 July 20114 The load flow newton -Raphson solution techniqueIntroductionGeneral fomulationLoad flow caseJacobian matrixSolution outline22 July 20115 The load flow decoupled AC load of load of load flow solution methods22 July 20116 The load flow G mez Exp sito, A. J. Conejo, C. Ca izares. Electric Energy Systems. Analysis and Operation . CRC Press, Boca Raton, Florida, R. Bergen, V. Vittal. Power Systems Analysis . Second Edition. Prentice Hall, Upper Saddle River, New Jersey, July 201171. Introduction22 July 20118 Introduction A snapshot of the system Most used tool in steady State power system analysis Knowning the demand and/or generation of power in each bus, find out: buses voltages load flow in lines and transformers22 July 20119 Introduction The problem is described throught a non-lineal system of equations Need of iterative solution techniques Solution technique: accuracy vs.

2 Computing time22 July 201110 Introduction Applications:1. On-line analyses State estimation Security Economic analyses22 July 201111 Introduction2. Off-line analyses Operation analyses Plannig analyses Network expansion planning Power exchange planning Security and adecuacy analyses-Faults-Stability22 July 2011122. problem formulation22 July 201113 problem formulationTwo-bus caseWe want to find out therelationship betweenand in all buses of the powersystemiiijQPS ijiieVV 2 VGYSY1V2 IGY1 ITransmission line model22 July 201114 problem formulationTwo-bus case BUSBUSBUS212221121121 SGSSSG21S12G22S21G11V YIVVYYYYVVYYYYYYIIY)VV(YVIY)VV(YVI Using Kirchhoff laws:Matrix notation:22 July 201115 problem formulationTwo-bus case Complex power injected in each bus:)VYVY(VIVS)VYVY(VIVSIVjQPSSSIVjQPSSS *2*22*1*212*222*2*12*1*111*111*22222D2G2 *11111D1G1 22 July 201116 problem formulationTwo-bus case continuationNotation:Replacing:)(j21kkk2 222)(j21kkk1111jiijikikk2k2k1k1iikeVYVjQ PeVYVjQPeVVeYY 22 July 201117 problem formulationTwo-bus case continuationTherefore, the no lineal equations for the 2 buses network are.

3 2,1i)(sinVYVQ)(cosVYVP ikki21kkikiiikki21kkikii 22 July 201118 problem formulationMatrixYbus Two bus case SGSSSG22211211 BUSYYYYYYYYYYY22 July 201119 General building rulesMatrix admittance of nodei,,equals the algebraic sum of all the admittances connected to node admittance between nodes i and k, ,equals the negative of the sum of all admittances connecting nodes i and 22 July 201120 problem formulationMatrix Ybus Caracteristics very sparse(>90% for more than 100 buses)BUSYBUSYBUSY22 July 201121 YbusexampleShunt element Series elementline) per (two July 201122 problem formulationGeneral equations 2n equations (static load flow equations) systembusnn,..,1i2)(sinVYVQQQ1)(cosVYVPP P ikkin1kkikiDiGiiikkin1kkikiDiGii 22 July 201123 problem formulationGeneral equationskiikikikikikikikikn1kkiiikikiki kn1kkiiBjGYwhere]2[)cosBsenG(VVQ]1[)senB cosG(VVP Polar representation for voltages and rectangular for admittances 22 July 201124 problem formulationGeneral equations 4n variables If 2n variables are specified, the other 2n are determined by equations [1] and [2]n.

4 ,1i;Q,P,,Viiii 22 July 201125 problem formulationBUS Classification1. PQ busesknown ( known, zero)known ( known, zero)unknownunknowniPiQiVi DiPGiPDiQGiQ22 July 201126 problem formulationBUS Classification2. PV busesknownunknownknown ( specified, known)known (specified)unknown ( unknown, known)unknowniViPiQiPiQDiPGiQGiPDiQi 22 July 201127 problem formulationBUS Classification3. Slack bus, generator with large (specified)known (specified, typically reference)unknown ( known, unknown)unknown ( known, unknown)1V1 1Q1P1DP1DQ1GQ1GP 0122 July 201128 problem formulationVariable types and limits Power balance n1in1iLOSSDiGin1in1iLOSSDiGiQQQPPP22 July 201129 problem formulationVariable types and limits Variable types Control variables(excepting slack bus) Non-control variables State variablesGiPiiV iGiVorQDiDiQP22 July 201130 problem formulationVariable types and limits Variable limits Voltage magnitude Power angle(every existing line) Power limitsmaxiiminiVVV maxkiki max,GiGimin,Gimax,GiGimin,GiQQQPPP gauss -Seidel solution technique22 July 201122 July 201132 gauss -Seidel solution technique )r()1r()r()1r(xx)x(Fx)x(Fx0)x(fNo lineal system.

5 Iteration Stoping rule22 July 201133 gauss -Seidel solution technique ExampleMany iterations! , , , , , , , ,0r2x; )x( )x(f)13()12()11()7()6(2)3(2)2(2)1()0(2)r ()1r(22 22 July 201134 gauss -Seidel solution techniqueAlgorithm beginning1)Known)BusesPV(m,..,2iQ,Q)Buse sPQ(n,..,1miV,VbusslackV)BusesPV(m,..,2i V)BusesPQ(n,..,1miQ)BusesPQ&PV(n,..,2iPm ax,Gimin,Gimaximini1iii 22 July 201135 gauss -Seidel solution techniqueAlgorithm beginning2)Build 3)Initialize voltagesn,..,2in,..,1miVV0ii0ii BUSY22 July 201136 gauss -Seidel solution techniquePQ buses4) PQ busesik;n,..,1k;n,..,1miYYBn,..,1miYjQPA VYVjQPY1 VVYVVYVIVjQPSiiikikiiiiinik1kkik*iiiiiin ik1kkik*iiii*ii*iii*i 22 July 201137 gauss -Seidel solution techniquePQ busesAt iteration (r+1)and bus i, the available values of voltages at previous buses are used:)r(kn1ikik)1r(k1i1kik*)r(ii)1r(iVBV B)V(AV 22 July 201138 gauss -Seidel solution techniquePV buses5)PV buses n1kkik*iin1kkik*ii*iii*iVYVQVYVIVjQPS22 July 201139 gauss -Seidel solution techniquePV busesAt iteration (r+1): 1i1kn1ik)r(kik)1r(kik*)r(i)1r(i)1r(iii)1 r(ii)1r(iiinik)r(kik1i1k*)r(i)1r(kik*)r( i)1r(iVBVB)V(AangleYjQPAV angleVY)V(VY)V(Q22 July 201140 gauss -Seidel solution techniquePV busesBeware of limits!

6 (more on this below)PQbecomesi&QQQQPQ becomesi&QQQQmax,i)1r(imax,i)1r(imin,i)1 r(imin,i)1r(i 22 July 201141 gauss -Seidel solution techniquePV ) Slack bus power (after convergence) n1kkk1*111)r(i)1r(iVYVjQPn,..,2i;VV6) Stop criterion22 July 201142 gauss -Seidel solution techniqueAlgorith final ) Compute line currents (after convergence) iVikVkikSikS0ikYLikY0ikILikILkiI0kiY0kiI 0kik0kiLikikLki0iki0ikLikkiLikYVI;Y)VV(I YVI;Y)VV(I 22 July 201143 gauss -Seidel solution techniqueAlgorith final ) Compute line complex power (after convergence);YVVY)VV(VS;YVVY)VV(VS*0ik*i i*Lik*k*iiik*0ki*kk*Lik*i*kkki iVikVkikSikS0ikYLikY0ikILikILkiI0kiY0kiI 22 July 201144 gauss -Seidel solution techniqueAlgorith final ) Compute losses (after convergence)7) If no convergence, go to step 4. i,kik,losslosskiikik,lossSSi,k;SSS22 July 201145 gauss -Seidel solution techniqueAlgorithm improvement Acceleration factor (in order to decrease the number of iterations): )r(i)1r(i)r(i)1r(iVVVV~ (generally recommended)22 July 201146 gauss -SeidelMatlab codefunction [Vfinal,angfinal,nite,P,Q,errorplot,tiem po]=Gaussgen(m,n,Ybus,Vmodini,Angini,P,Q ,tol,Vmax,Vmin,Qmax,Qmin)%-------------- ---------------------------------------- ---------------------------------------- ---------------------------------------- ------------%-function [Vmod,ang,nite,P,Qerrorplot,tiempo]=Gaus sgen(m,n,Ybus,Vmodini,Angini,P,Q,tol,Vma x,Vmin,Qmax,Qmin)%-Resuelve de forma general problema de carga por el m etodo de gauss -Seidel%-donde:% nudos PV.

7 (m=1 cuando no hay nudos PV)%-n nudos totales%-Ybus matriz de admitancias%-Vmodini tensiones iniciales modulo%-Angini angulos iniciales RADIANES%-P potencia activa inicial%-Q potencia reactiva inicial%-tol tolerancia para error en tension y potencia reactiva%-Vmax Vmin valores limites aceptables para las tensiones%-Qmax Qmin valores limites aceptables para las potencias reactivas%-nite n umero de iteraciones%-Vfinal vector con todas las potencias para cada iteraci on%-angfinal igual pero con los angulos%-tiempo, tiempo invertido en hacer las operaciones%---------------------------- ---------------------------------------- ---------------------------------------- ---------------------------------------% calculo de la matriz BB=zeros(n,n);for t=1:nfor k=1:nB(t,k)=Ybus(t,k)/Ybus(t,t);endend22 July 201147 gauss -SeidelMatlab code 2%calculos los valores en coordenadas cuadrangulares de la tensionV=zeros(n,1);for a=1:nV(a)=Vmodini(a)*exp(i*Angini(a));en dang=Angini;%empieza el bucle:error=1; %valores iniciales para poder entrar en el bucleerrorQ=1;nite=0;tic;while max(abs(error))>tol | max(abs(errorQ))>tolnite=nite+1;Vmod=abs (V);Vini=V;ang=angle(V);Qini=Q;%calculo las reactivas para los nudos PVif m>1for l=2:mAQ=0;for k=1:n AQ=AQ+Ybus(l,k)*V(k);endQ(l)=-imag((V(l) ')*AQ);end22 July 201148 gauss -SeidelMatlab code 3%calculo las A para todos los nudos:for a=1:nA(a)=(P(a)-i*Q(a))/Ybus(a,a);end%ca lculo los angulos nudos PV:for l=2:mAang=0;for k=1:n Aang=Aang+B(l,k)*V(k);endang(l)=angle(A( l)/((V(2))')-Aang+B(l,l)*V(l));endendfor a=1:nA(a)=(P(a)-i*Q(a))/Ybus(a,a);end%ah ora actualizo los voltajes otra vez:for a=1:nV(a)=Vmod(a)*exp(i*ang(a)).

8 End%ahora calculo los nudos PQAV=zeros(n,1);for p=m+1:nfor k=1:nAV(p)=AV(p)+B(p,k)*V(k);endV(p)=A(p )/((V(p))')-AV(p)+B(p,p)*V(p);enderror=V ini-V;errorQ=Qini-Q;errorplot(1,nite)=no rm(abs(error));Vfinal(:,nite)=abs(V);ang final(:,nite)=ang*180/pi; %paso a gradosend 22 July 201149 gauss -SeidelMatlab code 4%calculos los valores de potencia para el nudo slack:S1=0;for t=1:nS1=S1+(V(1)')*(Ybus(1,t)*V(t));end P(1)=real(S1);Q(1)=-imag(S1);tiempo=toc; %alerta por si se sobrepasan valores aceptables:if max(Vmod)>Vmax | min(Vmod)<Vmindisp(' SE SOBREPASAN LIMITES TENSIONES!!')endif max(Q)>Qmax | min(Q)<Qmindisp(' SE SOBREPASAN LIMITES REACTIVA!!')endang=(180/pi)*angle(V);Vmo =abs(V);%represento los errores por iteraci on:plot(1:nite,errorplot,'o');grid on;xlabel('iteraci on');ylabel('error por iteraci on');title('evoluci on error tesi on');%---------------------------------- --------------------------%% %% Realizado por Carlos Ruiz Mora octubre 2006 %% %%-------------------------------------- ----------------------%50 gauss -Seidelexample22 July 201122 July 201151G-S ExampleSolution tolerance is set to July 201152G-S ExampleData below.

9 Base power is SB=100 MVA:BusVoltage ( ) (slack) MW3-PC=100 MWQC=60 MVArLinImpedance ( ) + + + (each)22 July 201153 Solution and unknown:BusTypeDataUnknown1 SlackV1= 1= Q12 PVV2= P2= 2 Q23 PQP3= Q2= 3V322 July 201154 Solution procedure j5525j4020j155j4020j6030j2010j155j2010j3 515 YBUSBUSYONETWOTHREE22 July 201155 Solution procedure3. Voltage magnitude initialization (iteration 0):1VV033 PQ bus00033022 All buses but the reference oneVector form: July 201156 Solution procedurePer iteration:-PV buses (Q, ) i=2,..,m (bus 2) -PQ buses (V, ) i=m+1,..,n (bus 3)-Stopping criterios:a) convergence Sslackand power flows;b) no converge new iteration22 July 201157 Solution procedure4. PV buses: iteration (r+1) : bus 2: 22)1r(22)1r(222)r(323*)r(2)r(222*)r(2)1r (121*)r(2)1r(2 YjQPA]V[ angleVY)V(VY)V(VY)V(Q 22 July 201158 Solution procedure )r(323)1r(121*)1r(2)1r(2)1r(2 VBVB)V(AanglewhereiiikikYYB is a constant22 July 201159 Solution procedure5.

10 Buses PQ iteration (r+1): bus 3:)1r(232)1r(131*)r(33)1r(3 VBVB)V(AV whereiiikik33333 YYBYjQPA Are constants for PQ buses22 July 201160 Solution criterion:310 |QQ|&|VV|)r(j)1r(j)r(i)1r(i 2j ;3 ,2i 22 July 201161 Solution procedureIf )Slack ) Power flows )VYVYVY(VjQPS313212111*111*slack *Likk*iiikY)VV(VS-* 323123132112S,S,S,S,S,S22 July 201162 Solution ) Line losses:23,loss13,loss12,losslosskiikik,l osssSSSS3,2,1i,kSSS 7. If no coveergence, the procedure continues in Step July 201163 ImplementationMATLAB:22 July 201164 Solution11 iterations needed to attain the solutionIteration(pu) ,Q1(slack) July )( 1 )( 2 )( 3 22 July Error 10-7 Max. Error 10-4 Errorsfor|V| & |Q|:22 July 201167 Final SolutionBusP (MW)Q (MVAr)(pu) 0 V 22 July 201168 CheckingChecking using Power-World:ONETWOTHREE22 July 201169 CheckingBus 1 Bus 2 Bus 3 SolutionG-SPWG-SPWG-SPWV (pu) ( ) (MW) (MVAr) July 201170 CheckingVariableV PQError Max.