Baixe o app para aproveitar ainda mais
Prévia do material em texto
Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 1 de 47 VII – Fluxo de carga não linear: algoritmos básicos VII.1 – Formulação do problema básico Para um sistema elétrico, com NB barras, as equações básicas do fluxo de carga para NBk ,,2,1 L= são: ( )∑ ∈ += Km kmkmkmkmmkk BGVVP θθ sencos (VII.1) ( )∑ ∈ −= Km kmkmkmkmmkk BGVVQ θθ cossen (VII.2) constituindo um sistema de NB×2 equações e NB×4 variáveis com as seguintes características: • NB equações do tipo (VII.1); • NB equações do tipo (VII.2); • NB×4 variáveis (para uma dada barra k associam-se: kV , kθ , kP e kQ ). O fluxo de carga básico em redes de energia consiste em resolver este problema definindo (especificando) parte das variáveis ( NB×2 variáveis) e calculando as demais ( NB×2 variáveis). Na formulação clássica, a definição das variáveis que são especificadas e calculadas encontra-se na Tabela VII.1. Tabela VII.1 – Tipos de barra no fluxo de carga convencional. Tipo de barra Notação Dados Incógnitas Barra de carga PQ kP e kQ kV e kθ Tensão controlada PV kP e kV kθ e kQ Referência Vθ kV e kθ kP e kQ Se a barra está representando um ponto do sistema onde é conhecida a injeção líquida de potência (como DGk SSS −= , deve-se conhecer a potência gerada GS e a potência demandada DS ), esta pode ser especificada, restando calcular kV e kθ . Este tipo de barra é denominado barra de carga ou PQ. Se a barra está representando um ponto do sistema onde é possível o controle da magnitude da tensão kV , através do controle da injeção líquida de potência reativa DGk QQQ −= ou por intermédio do ajuste do tap de algum transformador, esta tensão pode ser especificada, restando calcular kθ e kQ . Este tipo de barra é denominado barra de tensão controlada ou PV. Como as perdas não são conhecidas a priori, deve-se deixar pelo menos uma barra onde a injeção líquida de potência DGk SSS −= seja calculada. Este tipo de barra é denominado referência ou Vθ, sendo definidos kV e kθ . Observar que pelo menos uma injeção líquida de potência ativa DGk PPP −= e outra de potência reativa DGk QQQ −= precisam ser calculadas para fechar os balanços de potência ativa e reativa. Embora geralmente tais injeções de balanço sejam associadas à barra de referência angular, não é mandatório que sejam associadas a uma única barra. Uma vez resolvido o problema do fluxo de carga, é conhecido o estado da rede, ou seja, são conhecidos os fasores tensão ( )kkV θ, de todas as barras ( )NBk ,,2,1 L= . Isto torna possível o cálculo de outras variáveis de interesse como as injeções de potência nas barras e os fluxos de potência nos ramos (linhas e transformadores). Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 2 de 47 Sejam NPQ e NPV, respectivamente, o número de barras PQ e PV da rede. O fluxo de carga pode ser decomposto em dois subsistemas de equações algébricas: • Subsistema 1 (dimensão NPVNPQ +×2 ) – Neste subsistema são dados: iP e iQ , { }PQ barras∈i : espii PP = espii QQ = jP e jV , { }PV barras∈j : espjj PP = espjj VV = kV e kθ , { }Vθ barra∈k : espkk VV = espkk θθ = e deseja-se determinar: iV e iθ , { }PQ barras∈i jθ , { }PV barras∈j resultando em um sistema de NPVNPQ +×2 equações algébricas não-lineares (funções quadráticas e trigonométricas) onde parte das incógnitas aparece de forma implícita ( )mkkm θθθ −= . O Subsistema 1 é constituído pelas seguintes equações: ( ) ( ) { } ( ) { } ∈=−− ∈=+− ∑ ∑ ∈ ∈ PQ barras0cossen PV e PQ barras0sencos S1 esp esp kBGVVQ kBGVVP Km kmkmkmkmmkk Km kmkmkmkmmkk θθ θθ • Subsistema 2 (dimensão 2+NPV ) – É resolvido após a solução do Subsistema 1. No Subsistema 2 são dados: kV e kθ , NBk ,,2,1 L= e deseja-se determinar: iP e iQ , { }Vθ barra∈i jQ , { }PV barras∈j resultando em um sistema de 2+NPV equações algébricas não-lineares onde todas as incógnitas ( )kk QP e aparecem isoladas de forma explícita. O Subsistema 2 é constituído pelas seguintes equações: ( ) ( ) { } ( ) { } ∈−= ∈+= ∑ ∑ ∈ ∈ Vθ e PV barrascossen Vθ barrasencos S2 kBGVVQ kBGVVP Km kmkmkmkmmkk Km kmkmkmkmmkk θθ θθ A Tabela VII.2 resume as características dos dois subsistemas. Observar que o Subsistema 2 pode ser facilmente resolvido após a solução do Subsistema 1, ou seja, após a determinação do estado ( )kkV θ, da rede. Por outro lado, a solução do Subsistema 1 exige um processo iterativo pois este é formado por equações não-lineares. Assim, os métodos de resolução do problema do fluxo de carga concentram-se na solução do Subsistema 1. Tabela VII.2 – Características dos subsistemas que constituem o fluxo de carga. Variáveis Subsistema Dimensão Especificadas Calculadas S1 NPVNPQ +×2 iP e iQ , { }PQi barras∈ jP e jV , { }PVj barras∈ kV e kθ , { }Vθ barra∈k iV e iθ , { }PQi barras∈ jθ , { }PVj barras∈ S2 2+NPV kV e kθ , NBk ,,2,1 L= iP e iQ , { }θVi barras∈ jQ , { }PVj barras∈ Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 3 de 47 Exemplo VII.1 – Considerando o sistema elétrico utilizado no Exemplo VI.1, cujos dados encontram-se na Figura VII.1, formular as equações referentes ao Subsistema 1 do fluxo de carga. pu 011 =V 1 2 ( ) pu 4,08,02 jS += 222 θVV = ( ) pu 1,001,0 jZ LT += 1S 12I Figura VII.1 – Sistema elétrico de duas barras. Solução Exemplo VII.1: Conforme as definições já apresentadas para o problema do fluxo de carga, a barra 2 é uma barra de carga pois sua injeção de potência é conhecida; a barra 1 é a referência. Assim, as variáveis especificadas (conhecidas) para o Subsistema 1 são as seguintes: pu 8,08,0022 esp 2 −=−=−= DG PPP pu 4,04,0022 esp 2 −=−=−= DG QQQ pu 0,1esp1 =V rad 0esp1 =θ A matriz admitância da rede é dada por: pu 29,849504,99010,99901,0 1,001,0 11 12 12 o −=−≈ + == jjZY − − = − − =⇒ − − = 9010,99010,9 9010,99010,9 e 9901,09901,0 9901,09901,0 1212 1212 BG YY YYY As equações do Subsistema 1 são as seguintes: ( ) ( )[ ]( )[ ] =−−− =++− 0cossen 0sencosS1 2222121212112 esp 2 2222121212112 esp 2 BVBGVVQ GVBGVVP θθ θθ ( ) ( )( ) =+−−−− =++−−− 09010,9cos9010,9sen9901,04,0 09901,0sen9010,9cos9901,08,0 S1 2222 2222 VV VV θθ θθ Observar que o Subsistema 1 é formado por duas equações não lineares e possui duas incógnitas: 2V e 2θ . As expressões do Subsistema 1 têm origem na equação de balanço de potência da Barra 2 (a potência injetada menos a potência que flui através da linha de transmissão é igual a zero) e diferem das equações obtidas na solução do Exemplo VI.1, oriundas da aplicação de uma equação de malha (soma de tensões), entretanto, a solução de ambos sistemas é mesma1: pu 9460,02 ≈V .e o61,4rad 0804,02 −=−≈θ . As incógnitas do Subsistema 1 podem ser agrupadas no vetor x dado por: } } NPQ NPVNPQ V x + = θ em que θ é o vetor dos ângulos das tensões nodais das barras PQ e PV e V é o vetor das magnitudes das tensões nodais das barras PQ. De forma mais compacta,o Subsistema 1 pode ser rescrito como: ( ) ( ) { }( ) { } ∈=−=∆ ∈=−=∆ PQ barras0, PV e PQ barras0,S1 esp esp kVQQQ kVPPP kkk kkk θ θ 1 Isto pode ser comprovado substituindo este valor nas expressões anteriores. Como os números apresentam alguns arredondamentos, o resultado fica próximo a zero, na ordem de 10-4. Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 4 de 47 colocando as funções kP∆ e kQ∆ na forma vetorial ( ) ( ) 0, 0, esp esp =−=∆ =−=∆ θ θ VQQQ VPPP em que P é o vetor das injeções de potência ativa nas barras PQ e PV e Q é o vetor das injeções de potência reativa nas barras PQ. Definindo a função vetorial ( )xg : ( ) }} NPQ NPVNPQ Q P xg + ∆ ∆ = o Subsistema 1 pode ser rescrito de forma simplificada através da seguinte expressão: ( ) 0=xg Este sistema de equações não-lineares pode ser resolvido por um número muito grande de métodos, sendo que os mais eficientes são os métodos de Newton e o método Desacoplado Rápido. Exemplo VII.2 – Empregando a notação na forma vetorial, determinar as variáveis e equações do Subsistema 1 do problema definido no Exemplo VII.1. Solução Exemplo VII.2: As variáveis e equações do problema são dadas respectivamente por: = = 2 2 VV x θθ ( ) ( )( ) ( ) ( ) = ∆ ∆ = ∆ ∆ = 0 0 , , 222 222 VQ VP xQ xP xg θ θ ( ) ( ) ( ) ( )2222222 2222222 9010,9cos9010,9sen9901,04,0, 9901,0sen9010,9cos9901,08,0, VVVQ VVVP +−−−−=∆ ++−−−=∆ θθθ θθθ Na figura a seguir encontra-se o gráfico das funções ( )222 ,VP θ∆ e ( )222 ,VQ θ∆ para valores em torno do ponto ( ) ( )9460,0;0804,0, 22 −=Vθ que corresponde a solução do problema, pois ( ) 09460,0;0804,02 =−∆P e ( ) 09460,0;0804,02 =−∆Q . A linha pontilhada indica a intersecção entre as duas funções, ou seja, corresponde aos valores para os quais a função ( ) ( )222222 ,, VQVP θθ ∆=∆ . -0.14 -0.12 -0.1 -0.08 -0.06 0.9 0.95 1 -1 -0.5 0 0.5 1 theta2V2 dP2(theta2,V2) dQ2(theta2,V2) dP2=dQ2 A figura anterior foi gerada utilizando o seguinte código MATLAB. Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 5 de 47 Solução Exemplo VII.2 (continuação): % disponivel em: http://slhaffner.phpnet.us/sistemas_de_energia_1/exemplo_VII_2.m clear all; theta2 = -0.15:.002:-0.05; V2 = 0.9:.002:1.0; [T,V] = meshgrid(theta2,V2); g1 = -0.8 - V.*(-0.9901.*cos(T)+9.9010.*sin(T)+0.9901.*V); g2 = -0.4 - V.*(-0.9901.*sin(T)-9.9010.*cos(T)+9.9010.*V); zero = zeros(size(V2,2),size(theta2,2)); figure(1) contour3(theta2,V2,g1,100); %colorbar; xlabel('theta2'); ylabel('V2'); hold on contour3(theta2,V2,g2,100); % determina pontos de interseccao entre as superficies g1 e g2 g11 = []; g22 = []; t = []; v = []; for k=1:size(g1,1) g11 = [g11 g1(k,:)]; g22 = [g22 g2(k,:)]; t = [t T(k,:)]; v = [v V(k,:)]; end i = find(abs(g11-g22)<=0.005); tt = t(i); vv = v(i); gg = g11(i); plot3(tt,vv,gg,'k'); Exemplo VII.3 – Para a rede de quatro barras cujos dados estão na Figura VII.2 e nas Tabelas VII.3 e VII.4, determinar as equações do fluxo de carga. 12Y shjb12 shjb12 23Y shjb23 shjb23 13Y shjb13 shjb13 34:1 a 34Y 14:1 ϕje 14Y 1 2 3 4 1S 3S 4S 2S 1V 2V 3V 4V shjb3 Figura VII.2 – Sistema exemplo de 4 barras. Tabela VII.3 – Dados das barras do sistema de 4 barras. Barra Tipo espV [pu] espθ [rad] espP [pu] espQ [pu] shb [pu] 1 Vθ 1,05 0,0 — — — 2 PQ — — –0,4 –0,2 — 3 PV 0,95 — 0,2 — 0,0675 4 PQ — — –0,8 –0,4 — Tabela VII.4 – Dados dos ramos do sistema de 4 barras. k m kmY [pu] shkmb [pu] kma [pu] kmϕ [rad] 1 2 0,3 – j2,0 0,01 — — 1 3 0,3 – j2,0 0,01 1 4 – j2,0 — — 0,15 2 3 0,2 – j1,0 0,02 — — 3 4 – j1,0 — 0,95 — Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 6 de 47 Solução Exemplo VII.3: Como já determinado na solução do Exemplo VI.2, a expressão da matriz admitância é dada por: +−− −+++++−− −+++− −−−++++ = − 3414343414 34343231334 2 3423132313 232312231212 1413121312141312 0 0 14 14 YYYaYe YajbjbjbYaYYYY YjbjbYYY YeYYjbjbYYY Y j shshsh shsh jshsh ϕ ϕ Substituindo os valores da Tabela VII.4, chega-se a: −+− −+−+− +−−+− ++−+−− = 395,009775,12989,0 95,0805,35,012,023,0 012,097,25,023,0 9775,12989,023,023,098,56,0 jjj jjjj jjj jjjj Y 2 − −− −− −− = 0002989,0 05,02,03,0 02,05,03,0 2989,03,03,06,0 G e − − − − = 395,009775,1 95,0805,312 0197,22 9775,12298,5 B As matrizes anteriores podem ser obtidas com o seguinte código MATLAB. % disponivel em: http://slhaffner.phpnet.us/sistemas_de_energia_1/exemplo_VII_3.m ramos = [ 1 2 0.3-2i 0.01 1.0 0; 1 3 0.3-2i 0.01 1.0 0; 1 4 -2i 0 1.0 0.15; 2 3 0.2-1i 0.02 1.0 0; 3 4 -1i 0 0.95 0]; bksh = [ 0; 0; 0.0675; 0]; nr = size(ramos,1); y = ramos(:,3); bsh = 1i*ramos(:,4); a = ramos(:,5); fi = ramos(:,6); Y = diag(1i*bksh); for k = 1:nr k1 = ramos(k,1); k2 = ramos(k,2); Y(k1,k2) = Y(k1,k2) - a(k)*exp(-1i*fi(k))*y(k); Y(k2,k1) = Y(k2,k1) - a(k)*exp(1i*fi(k))*y(k); Y(k1,k1) = Y(k1,k1) + a(k)*a(k)*y(k) + bsh(k); Y(k2,k2) = Y(k2,k2) + y(k) + bsh(k); end Y G = real(Y) B = imag(Y) Para cada uma das quatro barras do sistema, podem ser escritas duas equações, uma para a injeção de potência ativa e outra para a injeção de potência reativa, da forma como segue: ( ) ( )∑ ∑ ∈ ∈ −= += 1 1 1111 esp 11 1111 esp 11 cos cos Km mmmmm Km mmmmm BsenGVVQ senBGVVP θθ θθ { }4,3,21 =Ω { }4,3,2,11 =K ( ) ( )∑ ∑ ∈ ∈ −= += 2 2 22222 esp 2 22222 esp 2 cos cos Km mmmmm Km mmmmm BsenGVVQ senBGVVP θθ θθ { }3,12 =Ω { }3,2,12 =K ( ) ( )∑ ∑ ∈ ∈ −= += 3 3 3333 esp 33 3333 esp 3 esp 3 cossen sencos Km mmmmm Km mmmmm BGVVQ BGVVP θθ θθ { }4,2,13 =Ω { }4,3,2,13 =K 2 Observar que devido à presença de um transformador defasador a matriz admitância não é simétrica. Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 7 de 47 Solução Exemplo VII.3 (continuação): ( ) ( )∑ ∑ ∈ ∈ −= += 4 4 44444 esp 4 44444 esp 4 cossen sencos Km mmmmm Km mmmmm BGVVQ BGVVP θθ θθ { }3,14 =Ω { }4,3,14 =K Observar que as variáveis que se deseja calcular (incógnitas)encontram-se tanto do lado direito ( 2θ , 3θ , 4θ , 2V e 4V ) quanto do lado esquerdo ( 1P , 1Q e 3Q ) da igualdade. O Subsistema 1, de dimensão 52 =+× NPVNPQ , cujas incógnitas aparecem do lado direito ( 2θ , 3θ , 4θ , 2V e 4V ) é dado por: ( ) ( )∑ ∑ ∈ ∈ −= += 2 2 22222 esp 2 22222 esp 2 cossen sencos Km mmmmm Km mmmmm BGVVQ BGVVP θθ θθ { }3,2,12 =K ( )∑ ∈ += 3 3333 esp 3 esp 3 sencos Km mmmmm BGVVP θθ { }4,3,2,13 =K ( ) ( )∑ ∑ ∈ ∈ −= += 4 4 44444 esp 4 44444 esp 4 cossen sencos Km mmmmm Km mmmmm BGVVQ BGVVP θθ θθ { }4,3,14 =K As equações restantes constituem o Subsistema 2, de dimensão 32 =+NPV . Para este subsistema, as incógnitas aparecem do lado esquerdo ( 1P , 1Q e 3Q ) sendo dado por: ( ) ( )∑ ∑ ∈ ∈ −= += 1 1 1111 esp 11 1111 esp 11 cossen sencos Km mmmmm Km mmmmm BGVVQ BGVVP θθ θθ { }4,3,2,11 =K ( )∑ ∈ −= 3 3333 esp 33 cossen Km mmmmm BGVVQ θθ { }4,3,2,13 =K Observar que após a determinação das incógnitas do Subsistema 1 ( 2θ , 3θ , 4θ , 2V e 4V ), o estado da rede será completamente definido pois o fasor tensão da Barra 1 ( esp11 VV = e esp11 θθ = ) e a magnitude da tensão da Barra 3 ( esp33 VV = ) são conhecidos, o que torna a solução Subsistema 2 trivial (basta substituir os valores das tensões e ângulos e calcular os somatórios). Reescrevendo-se as equações do Subsistema 1, tem-se um sistema de cinco equações não lineares e cinco incógnitas ( 2θ , 3θ , 4θ , 2V e 4V ), das quais três aparecem de forma implícita (os ângulos de fase 2θ , 3θ e 4θ ): ( ) ( ) ( )[ ]23232323esp322222222221212121esp12esp2 sencossencossencos θθθθθθ BGVBGVBGVVP +++++= ( ) ( ) ( )[ ( )]343434344 33333333 esp 332323232231313131 esp 1 esp 3 esp 3 sencos sencossencossencos θθ θθθθθθ BGV BGVBGVBGVVP ++ ++++++= ( ) ( ) ( )[ ]44444444443434343esp341414141esp14esp4 sencossencossencos θθθθθθ BGVBGVBGVVP +++++= ( ) ( ) ( )[ ]23232323esp322222222221212121esp12esp2 cossencossencossen θθθθθθ BGVBGVBGVVQ −+−+−= ( ) ( ) ( )[ ]44444444443434343esp341414141esp14esp4 cossencossencossen θθθθθθ BGVBGVBGVVQ −+−+−= Colocando na forma de diferenças de potência ( ( ) 0,esp =−=∆ θVPPP kkk e ( ) 0,esp =−=∆ θVQQQ kkk ) e considerando que 1cos =kkθ e 0sen =kkθ , tem-se: ( ) ( )[ ] 0sencossencos 23232323esp322221212121esp12esp2 =++++− θθθθ BGVGVBGVVP ( ) ( )[ ( )] 0sencos sencossencos 343434344 esp 33332323232231313131 esp 1 esp 3 esp 3 =++ +++++− θθ θθθθ BGV VGBGVBGVVP ( ) ( )[ ] 0sencossencos 44443434343esp341414141esp14esp4 =++++− VGBGVBGVVP θθθθ Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 8 de 47 Solução Exemplo VII.3 (continuação): ( ) ( )[ ] 0cossencossen 23232323esp322221212121esp12esp2 =−+−−− θθθθ BGVVBBGVVQ ( ) ( )[ ] 0cossencossen 44443434343esp341414141esp14esp4 =−−+−− VBBGVBGVVQ θθθθ Substituindo-se os valores conhecidos mostrados na Tabela VII.3 e pela matriz admitância, têm-se as seguintes equações que constituem o Subsistema 1: ( ) ( )[ ] 0sen1cos2,095,05,0sen2cos3,005,14,0 2323221212 =+−+++−−− θθθθ VV ( ) ( )[ ( )] 0sen95,0cos0 95,05,0sen1cos2,0sen2cos3,005,195,02,0 34344 323223131 =++ +×++−++−− θθ θθθθ V V ( ) ( )[ ] 00sen95,0cos095,0sen9775,1cos2989,005,18,0 4434341414 =++++−−− VV θθθθ ( ) ( )[ ] 0cos1sen2,095,097,2cos2sen3,005,12,0 2323221212 =−−++−−−− θθθθ VV ( ) ( )[ ] 03cos95,0sen095,0cos9775,1sen2989,005,14,0 4434341414 =+−+−−−− VV θθθθ Agrupando as variáveis do Subsistema 1 no vetor x e as equações do Subsistema 1 e utilizando ( )xPk e ( )xQk para representar as injeções de potência ativa e reativa calculadas em função das variáveis x , tem-se: = 4 2 4 3 2 V V x θ θ θ ( ) ( ) ( ) ( ) { } { } ( ) ( ) { } { } ==∆ =−=∆ =−=∆ ==∆ =−=∆ =−=∆ =−=∆ 4,2PQ barras0 0 0 4,3,2PV e PQ barras0 0 0 0 S1 4 esp 44 2 esp 22 4 esp 44 3 esp 33 2 esp 22 Q xQQQ xQQQ P xPPP xPPP xPPP ( ) ( ) }} equações 2 equações 312 0S1 = =+=+ = ∆ ∆ = NPQ NPVNPQ Q P xg Na forma compacta, as expressões das equações do Subsistema 2 são dadas por: ( ) ( ) } { } { } ( ) ( ) { } { } = = = == 3,1Vθ e PV barras 1Vθ barras S2 33 11 11 xQQ xQQ xPP Exemplo VII.4 (alternativo, sem transformador defasador) – Para a rede de quatro barras cujos dados estão na Figura VII.3 e nas Tabelas VII.5 e VII.6, determinar as equações do fluxo de carga. 12Y shjb12 shjb12 23Y shjb23 shjb23 13Y shjb13 shjb13 34:1 a 34Y 1:41a 41Y 1 2 3 4 1S 3S 4S 2S1V 2V 3V 4V shjb3 Figura VII.3 – Sistema exemplo de 4 barras (alternativo). Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 9 de 47 Tabela VII.5 – Dados das barras do sistema de 4 barras (alternativo). Barra Tipo espV [pu] espθ [rad] espP [pu] espQ [pu] shb [pu] 1 Vθ 1,05 0,0 — — — 2 PQ — — –0,4 –0,2 — 3 PV 0,95 — 0,2 — 0,0675 4 PQ — — –0,8 –0,4 — Tabela VII.6 – Dados dos ramos do sistema de 4 barras. k m kmY [pu] shkmb [pu] kma [pu] 1 2 0,3 – j2,0 0,01 — 1 3 0,3 – j2,0 0,01 2 3 0,2 – j1,0 0,02 — 3 4 – j1,0 — 0,95 4 1 – j2,0 — 0,90 Solução Exemplo VII.4: Como já determinado na solução do Exemplo VI.3, a expressão da matriz admitância é dada por: +−− −+++++−− −+++− −−−++++ = 3441 2 4134344141 34343231334 2 3423132313 232312231212 414113121312411312 0 0 YYaYaYa YajbjbjbYaYYYY YjbjbYYY YaYYjbjbYYY Y shshsh shsh shsh Substituindo os valores das Tabelas VII.5 e VII.6, chega-se a: − −+−+− +−−+− +−+−− = 62,295,008,1 95,0805,35,012,023,0 012,097,25,023,0 8,123,023,098,56,0 jjj jjjj jjj jjjj Y −− −− −− = 0000 05,02,03,0 02,05,03,0 03,03,06,0 G e − − − − = 62,295,008,1 95,0805,312 0197,22 8,12298,5 B As matrizes anteriores podem ser obtidas com o seguinte código MATLAB. % disponivel em: http://slhaffner.phpnet.us/sistemas_de_energia_1/exemplo_VII_4.m ramos = [ 1 2 0.3-2i 0.01 1.0 0; 1 3 0.3-2i 0.01 1.0 0; 2 3 0.2-1i 0.02 1.0 0; 3 4 -1i 0 0.95 0; 4 1 -2i 0 0.9 0]; bksh = [ 0; 0; 0.0675; 0]; nr = size(ramos,1); y = ramos(:,3); bsh = 1i*ramos(:,4); a = ramos(:,5); fi = ramos(:,6); Y = diag(1i*bksh); for k = 1:nr k1 = ramos(k,1); k2 = ramos(k,2); Y(k1,k2) = Y(k1,k2) - a(k)*exp(-1i*fi(k))*y(k); Y(k2,k1) = Y(k2,k1) - a(k)*exp(1i*fi(k))*y(k); Y(k1,k1) = Y(k1,k1) + a(k)*a(k)*y(k) + bsh(k); Y(k2,k2) = Y(k2,k2) + y(k) + bsh(k); end Y G = real(Y) B = imag(Y) Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio HaffnerVersão: 10/9/2007 Página 10 de 47 Solução Exemplo VII.4 (continuação): Como os tipos das barras e as conexões existentes são as mesmas, as expressões obtidas no exemplo anterior permanecem válidas, havendo diferença apenas nas expressões cujos valores numéricos foram substituídos. Para os valores das Tabelas VII.5 e VII.6 que definem a matriz admitância, têm-se as seguintes equações que constituem o Subsistema 1: :2P∆ ( ) ( )[ ] 0sen1cos2,095,05,0sen2cos3,005,14,0 2323221212 =+−+++−−− θθθθ VV :3P∆ ( ) ( )[ ( )] 0sen95,0cos0 95,05,0sen1cos2,0sen2cos3,005,195,02,0 34344 323223131 =++ +×++−++−− θθ θθθθ V V :4P∆ ( ) ( )[ ] 00sen95,0cos095,0sen8,1cos005,18,0 4434341414 =++++−− VV θθθθ :2Q∆ ( ) ( )[ ] 0cos1sen2,095,097,2cos2sen3,005,12,0 2323221212 =−−−+−−−− θθθθ VV :4Q∆ ( ) ( )[ ] 062,2cos95,0sen0cos8,1sen005,14,0 44343241414 =+−+−−− VVV θθθθ VII.2 – Resolução de sistemas algébricos não lineares pelo método de Newton- Raphson Considere, inicialmente, um sistema unidimensional (neste caso, os vetores ( )xg e x possuem apenas um componente e podem ser representados por escalares) formado por uma equação do tipo: ( ) 0=xg (VII.3) Resolver este sistema significa determinar o valor de x tal que a função ( )xg seja nula. Sendo ( )xg uma função contínua com suas derivadas contínuas a expansão em série de Taylor em torno de um ponto conhecido 0x é dada por3: ( ) ( ) ( )( ) ( )( ) K+− ∂ ∂ +− ∂ ∂ += 20 2 02 0 0 0 !2 1 !1 1 xx x xg xx x xg xgxg Desprezando-se todos os termos após a derivada primeira, isto é, aproximando-se a função ( )xg por uma reta, conforme mostra a Figura VII.4, chega-se a: ( ) ( ) ( ) ( )000 xx x xg xgxg − ∂ ∂ +≈ (VII.4) A partir da expressão da aproximação linear (VII.4) é possível determinar o ponto 1x no qual o valor desta aproximação é nulo, ou seja, ( ) 01 =xg (vide Figura VII.4) da seguinte forma: ( ) ( ) ( )( ) 001001 =− ∂ ∂ += xx x xg xgxg (VII.5) ( ) ( )01001 xg x xg xx − ∂ ∂ −= x0 x1 g(x) ( )0xg ( ) ( )001 xxgxg ∆+= x Equação da reta tangente por x0: ( ) ( ) ( )( )000 xx x xg xgxg − ∂ ∂ += Ponto no qual a reta tangente por x0 é nula: x1 010 xxx −=∆ Figura VII.4 – Interpretação geométrica do método de Newton. 3 Embora a função seja de apenas uma variável x, adotou-se a notação de derivada parcial para facilitar a extensão para o caso multidimensional. Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 11 de 47 Definindo-se 001010 xxxxxx ∆+=⇒−=∆ e substituindo-se na expressão (VII.5), tem-se: ( )} ( ) ( ) ( ) 0 000 01 0 01 =− ∂ ∂ += ∆∆+ 48476 xxx xx x xg xgxg ( ) ( ) ( ) 000000 =∆ ∂ ∂ +=∆+ x x xg xgxxg (VII.6) ( ) ( )0100 xg x xg x − ∂ ∂ −=∆ (VII.7) onde 0x é a aproximação inicial e 001 xxx ∆+= uma primeira aproximação. Observar que 1x não é a solução da equação inicial (VII.3) e sim a solução de uma aproximação linear dada por (VII.4). A solução da equação (VII.3) é obtida repetindo-se este processo até que o módulo da função ( )xg esteja suficientemente próximo de zero (dentro de uma tolerância definida). Isto sugere o seguinte algoritmo denominado método de Newton-Raphson. Algoritmo do método de Newton-Raphson unidimensional i. Fazer 0=ν e escolher uma aproximação inicial 0xx =ν . ii. Calcular o valor da função ( )xg , no ponto νxx = : ( )νxg iii. Comparar o valor calculado ( )νxg com a tolerância especificada ε: se ( ) εν ≤xg , então νxx = será a solução procurada (dentro da faixa de tolerância ±ε); se ( ) εν >xg , prosseguir. iv. Linearizar a função ( )xg em torno do ponto ( )( )νν xgx , . Isto se resume na determinação da seguinte derivada – vide equação (VII.6): ( ) x xg ∂ ∂ ν v. Calcular a correção νx∆ que resolve o problema linearizado (VII.6) – vide equação (VII.7): ( ) ( )ννν xg x xg x 1− ∂ ∂ −=∆ vi. Determinar a nova estimativa de x passa a ser: ννν xxx ∆+=+1 vii. Fazer 1+= νν e voltar para o Passo (ii). Exemplo VII.5 – Utilizando o método de Newton-Raphson, determinar a solução para a equação4 xx sen2 −= , considerando uma tolerância 001,0=ε . Solução Exemplo VII.5: Inicialmente, faz-se: ( ) 0sen2 =+−= xxxg ( ) ( ) xxx xx xg cos1sen2 +=+− ∂ ∂ = ∂ ∂ Considerando uma solução inicial 00 =x obtêm-se os resultados mostrados na Tabela VII.7. 4 x em radianos. Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 12 de 47 Solução Exemplo VII.5 (continuação): Tabela VII.7 – Resultados parciais do processo iterativo – método de Newton-Raphson ( 00 =x ). ν x ν ( )νxg ( ) x xg ∂ ∂ ν νx∆ 0 0 –2 2 1 1 1 –0,1585 1,5403 0,1029 2 1,1029 –0,0046 1,4510 0,0031 3 1,1061 –4,4×10-6 — — Portanto, para uma tolerância 001,0=ε , a solução 1061,1=x foi obtida após 3 iterações, seguindo o processo ilustrado na Figura VII.5. -1 -0.5 0 0.5 1 1.5 2 -4 -3 -2 -1 0 1 x g(x) Figura VII.5 – Processo de convergência para 00 =x . Observar que a escolha da solução inicial afeta o processo de convergência conforme pode ser comprovado pela comparação entre os resultados obtidos com 00 =x e com 20 =x , quando foram necessárias 4 iterações (vide Tabela VII.8 e Figura VII.6). Tabela VII.8 – Resultados parciais do processo iterativo – método de Newton-Raphson ( 20 =x ). ν xν ( )νxg ( ) x xg ∂ ∂ ν νx∆ 0 2 0,9093 0,5838 –1,5574 1 0,4426 –1,1291 1,9036 0,5931 2 1,0357 –0,1040 1,5099 0,0689 3 1,1046 –0,0021 1,4495 0,0014 4 1,1061 –9,10×10-7 — — -1 -0.5 0 0.5 1 1.5 2 -4 -3 -2 -1 0 1 x g(x) Figura VII.6 – Processo de convergência para 20 =x . Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 13 de 47 Uma variante do método de Newton-Raphson é obtida considerando-se a derivada constante, isto é, o Passo (iv) do algoritmo é realizado uma única vez, na primeira iteração quando 0=ν : ( ) ( ) x xg x xg ∂ ∂ = ∂ ∂ 0ν Utilizando-se derivada constante, em geral, o número de iterações (para uma tolerância definida) é maior que no método de Newton original, mas cada uma das iterações se torna mais rápida pois não é necessário recalcular a derivada. Exemplo VII.6 – Utilizando o método de Newton-Raphson com derivada constante (Von Mises), determinar a solução para a equação xx sen2 −= , considerando uma tolerância 001,0=ε . Solução Exemplo VII.6: Considerando uma solução inicial 00 =x obtêm-se os resultados mostrados na Tabela VII.9. Tabela VII.9 – Resultados parciais do processo iterativo – método de Von Mises. ν xν ( )νxg ( ) x xg ∂ ∂ ν ν x∆ 0 0 –2 2 1 1 1 –0,1585 2 0,0793 2 1,0793 –0,0391 2 0,0196 3 1,0988 –0,0105 2 0,0052 4 1,1041 –0,0029 2 0,0014 5 1,1055 –0,0008 — — Portanto, para uma tolerância 001,0=ε , 1055,1=x . Como esperado, com a utilização de derivada constante, foi necessário realizar um número maior de iterações (5 ao invés de 3). Exercício VII.1 – Utilizandoos métodos de Newton-Raphson e de Von Mises, determinar, com uma tolerância 001,0=ε , o valor de x tal que 53sen 2 +−= xxe x . Considere-se, agora, a resolução do seguinte sistema n-dimensional: ( ) 0=xg (VII.8) onde ( )xg é uma função vetorial e x é o vetor das incógnitas, isto é: ( ) ( ) ( ) ( ) = xg xg xg xg n M 2 1 = nx x x x M 2 1 A solução do sistema de equações (VII.8) é obtida a partir de uma extensão do algoritmo de Newton- Raphson anterior, descrito a seguir. Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 14 de 47 Algoritmo do método de Newton-Raphson n-dimensional i. Fazer 0=ν e escolher uma aproximação inicial 0xx =ν . ii. Calcular ( )xg , no ponto νxx = : ( )νxg iii. Testar convergência: se ( ) [ ]nixg i ,1 para ∈≤ εν , então o processo convergiu para a solução νxx = ; caso contrário, prosseguir. iv. Linearizar a função vetorial ( )xg em torno do ponto ( )( )νν xgx , . Isto se resume na determinação da seguinte matriz de derivadas, denominada matriz Jacobiana: ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ = ∂ ∂ = n nnn n n x xg x xg x xg x xg x xg x xg x xg x xg x xg x xg xJ ννν ννν ννν ν ν L MOMM L L 21 2 2 2 1 2 1 2 1 1 1 v. Calcular a correção νx∆ que resolve o problema linearizado: ( )[ ] ( )ννν xgxJx 1−−=∆ vi. Determinar a nova estimativa de x que passa a ser: ννν xxx ∆+=+1 vii. Fazer 1+= νν e voltar para o Passo (ii). Exemplo VII.7 – Utilizando o método de Newton-Raphson, determinar a solução considerando uma tolerância 001,0== yx εε , para o seguinte sistema de equações: 62 42 2 =+ =+ yx yx Solução Exemplo VII.7: Inicialmente, faz-se: ( ) ( ) 062, 042, 2 21212 21211 =−+= =−+= xxxxg xxxxg ( ) = ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ = ∂ ∂ = 2 2 2 1 2 2 1 1 1 22 12 x x g x g x g x g x xg J Considerando uma solução inicial 001 =x e 3 0 2 =x obtêm-se os resultados mostrados na Tabela VII.10. Tabela VII.10 – Resultados parciais do processo iterativo – método de Newton-Raphson n-dimensional. ν ν ν 2 1 x x ( )( )ν ν xg xg 2 1 ( )[ ] 1−− νxJ νν 2 1 x x ∆ ∆ 0 0 3 –1 3 2,02,0 1,06,0 − − 0,9 –0,8 1 0,9 2,2 0 0,64 294,0294,0 147,0647,0 − − 0,094 –0,188 2 0,994 2,012 0 0,0354 331,0331,0 165,0665,0 − − 0,00586 –0,0117 3 0,999977 2,000046 0 0,000137 — — Portanto, para uma tolerância 001,0== yx εε , 99977,01 == xx e 000046,22 == xy . Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 15 de 47 Solução Exemplo VII.7 (continuação): Observar que a partir da primeira iteração a função ( )νxg1 passa a apresentar valor nulo. Isto ocorre porque a aproximação linear empregada para representar esta função corresponde à própria função pois esta é de primeira ordem nas variáveis 1x e 2x . Exemplo VII.8 – Utilizando o método de Von Mises, determinar a solução do sistema de equações do exemplo anterior. Solução Exemplo VII.8: Considerando uma solução inicial 001 =x e 302 =x obtêm-se os resultados mostrados na Tabela VII.11. Tabela VII.11 – Resultados parciais do processo iterativo – método de Von Mises n-dimensional. ν ν ν 2 1 x x ( )( )ν ν xg xg 2 1 ( )[ ] 1−− νxJ νν 2 1 x x ∆ ∆ 0 0 3 –1 3 2,02,0 1,06,0 − − 0,9 –0,8 1 0,9 2,2 0 0,64 2,02,0 1,06,0 − − 0,064 –0,128 2 0,964 2,072 0 0,2212 2,02,0 1,06,0 − − 0,00221 –0,0442 3 0,9861 2,0278 0 0,08406 2,02,0 1,06,0 − − 0,0084 –0,0168 4 0,99452 2,01095 0 0,03297 2,02,0 1,06,0 − − 0,0033 –0,0066 5 0,99782 2,00436 0 0,01309 2,02,0 1,06,0 − − 0,00131 –0,00262 6 0,99913 2,00174 0 0,00522 2,02,0 1,06,0 − − 0,00052 –0,00104 7 0,99965 2,00069 0 0,00208 2,02,0 1,06,0 − − 0,00021 –0,00042 8 0,99986 2,00028 0 0,000834 — — Portanto, para uma tolerância 001,0== yx εε , 99986,01 == xx e 00028,22 == xy . Exercício VII.2 – Utilizando os métodos de Newton-Raphson e de Von Mises, determinar, com uma tolerância 001,0=ε , a solução do seguinte sistema de equações: 52 43 2 2 −=− =+ yxy xyx VII.3 – Fluxo de carga pelo método de Newton-Raphson O método de Newton-Raphson é aplicado para a resolução do Subsistema 1 (S1) sendo dado por: ( ) ( )( ) ( ) { } ( ) { } ∈=−=∆ ∈=−=∆ ⇒ =−=∆ =−=∆ = PQ barras0, PV e PQ barras0, 0, 0, S1 esp esp esp esp kVQQQ kVPPP VQQQ VPPP kkk kkk θ θ θ θ (VII.9) De acordo com o algoritmo de Newton-Raphson deseja-se determinar o vetor das correções x∆ , o que exige a resolução do sistema linear dado por: Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 16 de 47 ( ) ( ) υυυ xxJxg ∆−= ⇒ ( )[ ] ( )ννν xgxJx 1−−=∆ onde ( ) PQ PV PQ ← +← ∆ ∆ = υ υ υ Q P xg PQ PV PQ ← +← = υ υ υ θ V x PQ PV PQ ← +← ∆ ∆ =∆ υ υ υ θ V x ( ) ( ) ( ) ( ) ( ) ( ) PQ PV PQ PQPVPQ ← +← ↑ + ↑ ∂ ∆∂ ∂ ∆∂ ∂ ∆∂ ∂ ∆∂ = ∂ ∂ = υ υ υ θ θ V QQ V PP x xg xJ Considerando as expressões dos vetores P∆ e Q∆ e que espP e espQ são constantes, a matriz Jacobiana pode ser rescrita da seguinte maneira: ( ) ( ) ( ) ( ) ( ) PQ PV PQ PQPVPQ ,, ,, ← +← ↑ + ↑ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ −= υ υ θ θ θ θ θ θ V VQVQ V VPVP xJ sendo as submatrizes representadas por: ( ) θ θ ∂ ∂ = ,VPH ( ) V VPN ∂ ∂ = θ, ( ) θ θ ∂ ∂ = ,VQ M ( ) V VQ L ∂ ∂ = θ, Assim, a equação que define a aplicação do método de Newton ao fluxo de carga fica sendo: ∆ ∆ = ∆ ∆ υ υυ υ υ θ VLM NH Q P (VII.10) Considerando que: ( ) ( )∑ ∑ Ω∈ ∈ ++= =+= km kmkmkmkmmkkkk Km kmkmkmkmmkk BGVVGV BGVVP θθ θθ sencos sencos 2 ( ) ( )∑ ∑ Ω∈ ∈ −+−= =−= km kmkmkmkmmkkkk Km kmkmkmkmmkk BGVVBV BGVVQ θθ θθ cossen cossen 2 as submatrizes que compõem a matriz Jacobiana são dadas por: ( ) ( ) ( ) Ω∉= Ω∈−= ∂ ∂ = +−= ∂ ∂ = ∂ ∂ = ∑ Ω∈ kkl kklklklkllkl k kl m kmkmkmkmmk k k kk lH lBGVVPH BGVVPH VPH k 0 cossen cossen , θθ θ θθ θ θ θ Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 17 de 47 ( ) ( ) ( ) Ω∉= Ω∈+= ∂ ∂ = ++= ∂ ∂ = ∂ ∂ = ∑ Ω∈ kkl kklklklklk l k kl m kmkmkmkmmkkk k k kk lN lBGV V PN BGVGV V PN V VPN k 0 sencos sencos2 , θθ θθ θ ( ) ( ) ( ) Ω∉= Ω∈+−= ∂ ∂ = += ∂ ∂ = ∂ ∂ = ∑ Ω∈ kkl kklklklkllk l k kl m kmkmkmkmmk k k kk lM lBGVVQM BGVVQM VQ M k 0 sencos sencos , θθ θ θθ θ θ θ ( ) ( ) ( ) Ω∉= Ω∈−= ∂ ∂ = −+−= ∂ ∂ = ∂ ∂ = ∑ Ω∈ kkl kklklklklk l k kl m kmkmkmkmmkkk k k kk lL lBGV V Q L BGVBV V Q L V VQ L k 0 cossen cossen2 , θθ θθ θ Assim, os passos a serem executados para solução do fluxo de carga pelo método de Newton são os seguintes: Fluxo de carga pelo método de Newton-Raphson – Algoritmo i. Fazer 0=υ e escolher os valores iniciais dos ângulos das tensões das barras PQ e PV ( )0θθθ υ == e as magnitudes das tensões das barras PQ ( )0VVV == υ . ii. Calcular: ( )θ,VPk para as barras PQ e PV ( )θ,VQk para as barras PQ e determinar o vetor dos resíduos (“mismatches”) υP∆ e υQ∆ . iii. Testar a convergência: se { }{ } Pkk P ευ ≤∆+∈ PVPQmax e { }{ } Qkk Q ευ ≤∆∈ PVmax , o processo convergiu para a solução ( )υυ θ,V ; caso contrário, continuar. iv. Calcular a matriz Jacobiana: ( ) ( ) ( )( ) ( ) −= υυυυ υυυυ υυ θθ θθθ ,, ,, , VLVM VNVHVJ v. Determinar a nova solução ( )11 , ++ υυ θV , onde: υυυ υυυ θθθ VVV ∆+= ∆+= + + 1 1 sendo υV∆ e υθ∆ obtidos com a solução do seguinte sistema linear: ( ) ( )( ) ( ) ∆ ∆ = ∆ ∆ υ υυ υυυυ υυυυ υ υ θ θθ θθ VVLVM VNVH Q P ,, ,, vi. Fazer 1+= υυ e voltar para o Passo (ii). Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 18 de 47 Após a determinação do fasor tensão de todas as barras, a solução do Subsistema 2 (S2) é trivial, sendo obtida através das expressões: ( ) ( ) { } ( ) { } ∈−= =+= = ∑ ∑ ∈ ∈ referência e PV barrascossen referência de barrasencos 2S kBGVVQ kBGVVP Km kmkmkmkmmkk Km kmkmkmkmmkk θθ θθ (VII.11) Exemplo VII.9 – Utilizando o método Newton, determinar a solução do problema do fluxo de carga correspondente ao sistema elétrico de duas barras utilizado no Exemplo VII.1, considerando uma tolerância 001,0== QP εε . Solução Exemplo VII.9: Conforme já determinado, as incógnitas e equações do Subsistema 1 são as seguintes: = 2 2 V x θ ( ) ( )( ) =+−−−−=∆ =++−−−=∆ 09010,9cos9010,9sen9901,04,0 09901,0sen9010,9cos9901,08,0 S1 22222 22222 VVQ VVP θθ θθ Para este problema a matriz Jacobiana apresenta os seguintes elementos, para os quais { }12 =Ω : ( ) ( ) ( )212121211222222 2 2 22 cossencossen , 2 θθθθ θθ θ BGVVBGVVPHVPH m mmmmm +−=+−=∂ ∂ == ∂ ∂ = ∑ Ω∈ ( ) ( ) ( )2121212112222222222 2 2 22 sencos2sencos2 , 2 θθθθθ BGVGVBGVGV V PN V VPN m mmmmm ++=++=∂ ∂ == ∂ ∂ = ∑ Ω∈ ( ) ( ) ( )212121211222222 2 2 22 sencossencos , 2 θθθθ θθ θ BGVVBGVVQMVQM m mmmmm +=+=∂ ∂ == ∂ ∂ = ∑ Ω∈ ( ) ( ) ( )2121212112222222222 2 2 22 cossen2cossen2 , 2 θθθθ θ BGVBVBGVBV V QL V VQ L m mmmmm −+−=−+−=∂ ∂ == ∂ ∂ = ∑ Ω∈ Levando em conta a matriz admitância da rede (já calculada na solução do Exercício VII.1), e os valores especificados para pu 1esp11 == VV e rad 0 esp 11 == θθ , têm-se os seguintes elementos: ( )22222 cos9010,9sen9901,0 θθ += VH ( )22222 sen9010,9cos9901,09802,1 θθ +−+= VN ( )22222 sen9010,9cos9901,0 θθ +−= VM ( )22222 cos9010,9sen9901,0802,19 θθ −−+= VL Neste caso, como a matriz Jacobiana é dada por: −= 2222 2222 LM NH J então sua inversa é dada por: − − − − = −= − − 2222 2222 22222222 1 2222 22221 1 HM NL MNLHLM NH J Considerando uma solução inicial rad 002 =θ e pu 102 =V , obtêm-se os resultados mostrados na Tabela VII.12. Portanto, para uma tolerância 001,0== QP εε , a solução do Subsistema 1 é dada por: pu 9461,02 =V e o61,4rad 0804,02 −=−=θ , mesmo valor obtido na solução do Exemplo VI.1. Na Tabela VII.12 é importante observar os valores obtidos nas submatrizes que constituem o Jacobiano. Enquanto os elementos das submatrizes H e L possuem maior valor absoluto e apresentam variações pequenas ao longo do processo iterativo, os elementos das submatrizes N e M possuem menor valor absoluto e variam de forma significativa. Modelagem e Análise de Sistemas Elétricos em Regime Permanente Fluxo de carga não linear: algoritmos básicos – Sérgio Haffner Versão: 10/9/2007 Página 19 de 47 Solução Exemplo VII.9 (continuação): Tabela VII.12 – Resultados parciais do processo iterativo – fluxo de carga Newton. ν ν νθ 2 2 V ( )( )ν ν xQ xP 2 2 ∆ ∆ ( )[ ]νxJ ( )[ ] 1−− νxJ ννθ 2 2 V∆ ∆ 0 0 1 –0,8 –0,4 9010,99910,0 9901,09010,9 − −− 1000,00100,0 0100,01000,0 − –0,0760 –0,0480 1 –0,0760 0,9520 –0,0418 –0,0463 0543,96555,1 1462,03270,9 − −− 1101,00195,0 0017,01069,0 − –0,0044 –0,0059 2 –0,0804 0,9461 –0,0003 –0,0004 — — — Os resultados mostrados na Tabela VII.12, foram obtidos executando-se a seguinte rotina em MATLAB. % disponivel em: http://slhaffner.phpnet.us/sistemas_de_energia_1/exemplo_VII_9.m clear all; saida=fopen('saida.txt','w'); p2=-0.8; q2=-0.4; v1=1; t1=0; v2=1; t2=0; x=[t2; v2]; G11=0.9901; G12=-0.9901; G21=-0.9901; G22=0.9901; B11=-9.901; B12=9.901; B21=9.901; B22=-9.901; for k=0:5, dp2=p2-v2*(v1*(G21*cos(t2)+B21*sin(t2))+G22*v2); dq2=q2-v2*(v1*(G21*sin(t2)-B21*cos(t2))-B22*v2); gx=[dp2; dq2]; h22=v2*v1*(-G21*sin(t2)+B21*cos(t2)); n22=2*v2*G22+v1*(G21*cos(t2)+B21*sin(t2)); m22=v2*v1*(G21*cos(t2)+B21*sin(t2)); l22=-2*v2*B22+v1*(G21*sin(t2)-B21*cos(t2)); Jac=-[h22 n22 m22 l22]; Jac1=-inv(Jac); dx=Jac1*gx; y=[k x(1) gx(1) Jac(1,1) Jac(1,2) Jac1(1,1) Jac1(1,2) dx(1)]; fprintf(saida,'%2.0f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f\n',y); y=[x(2) gx(2) Jac(2,1) Jac(2,2) Jac1(2,1) Jac1(2,2) dx(2)]; fprintf(saida,' %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f\n\n',y); x=x+dx; t2=x(1); v2=x(2); end p1=v1*(v2*(G12*cos(-t2)+B12*sin(-t2))+G11*v1); q1=v1*(v2*(G12*sin(-t2)-B12*cos(-t2))-B11*v1); y=[p1 q1]; fprintf(saida,'%8.4f %8.4f',y); fclose(saida); Por outro lado, o Subsistema 2 corresponde ao cálculo da injeção de potência na barra de referência: ( ) ( ) ( ) −= += = ∑ ∑ ∈ ∈ 1 1 111111 111111 cossen sencos 2S Km mmmmm Kmmmmmm BGVVQ BGVVP θθ θθ ( ) ( )[ ]( )[ ] −−= ++= 11112121212211 11112121212211 cossen sencos S2 BVBGVVQ GVBGVVP θθ θθ Substituindo os valores conhecidos, chega-se a: ( ) = = pu 4894,0 pu 8089,0 S2 1 1 Q P
Compartilhar