Buscar

seI7a sergio haffner

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você viu 3, do total de 19 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você viu 6, do total de 19 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você viu 9, do total de 19 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

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

Outros materiais