UNIVERSIDADE DE SÃO PAULO

ESCOLA DE ENGENHARIA DE SÃO CARLOS

PROGRAMA DE PÓS-GRADUAÇÃO EM ENGENHARIA ELÉTRICA

JÚLIA OLIVEIRA FERNANDES

ATIVIDADE FINAL – Fluxo de Potência Ótimo

SEL 5747 – Análise Estática de Sistemas de Energia Elétrica

São Carlos
2023
01- O Matpower resolve o FPO considerando as seguintes variáveis: P_gerado,
Q_gerado, V e θ. Sendo assim. Ainda ficam faltando as variáveis de tap e shunt.
Vamos trabalhar com o sistema de 14 barras do IEEE.

Para verificar a influência de todas as variáveis no problema de FPO,

vamos considerar que o shunt da barra tenha as seguintes opções:
17, 19, e 21 MVAr.
E que os taps dos transformadores das linhas 4 – 7, 4 – 9, e 5 – 6 tenham
as seguintes opções:
4 – 7 → 0.968, 0.978 e 0.988
4 – 9 → 0.949, 0.969 e 0.989
5 – 6 → 0.932, 0.962 e 0.992
É possível verificar que o problema é de ordem combinatória!
a) Encontrar qual o melhor valor da função custo de geração, considerando os
dados apresentados.
b) Encontrar qual o melhor valor para a função perdas.

Conforme dito no problema sugerido, o cálculo do FPO deve ser

implementado utilizando o pacote computacional Matpower, e o próprio já
apresenta o sistema de 14 barras em seus dados padrões.
Com isso, o primeiro passo é carregar o alimentador utilizado e apenas
editar os parâmetros desejados de acordo exposto no problema.
Inicialmente, como o problema é combinatorial, é necessário que a
quantidade total de possibilidades a serem analisadas, seja conhecida. Com
isso, o seguinte cálculo deve ser feito:

𝑃𝑡𝑜𝑡𝑎𝑙 = 𝑁𝑠ℎ𝑢𝑛𝑡 × 𝑁𝑡𝑎𝑝 4−7 × 𝑁𝑡𝑎𝑝 4−9 × 𝑁𝑡𝑎𝑝 5−6 (1)

𝑃𝑡𝑜𝑡𝑎𝑙 = 3 × 3 × 3 × 3 = 81 (2)

Sendo 𝑃𝑡𝑜𝑡𝑎𝑙 a quantidade total de possibilidades, 𝑁𝑠ℎ𝑢𝑛𝑡 o número de

possibilidades de valores para a barra shunt, e os demais 𝑁𝑡𝑎𝑝 𝑖−𝑗 o número de
possibilidades de valores dos taps entre as barras selecionadas (i – j),
resultando em um total de 81 possibilidades.
Devido ao número considerável de possibilidades e a não limitação do
tempo para uma solução ótima do problema, o método escolhido para o
cálculo do FPO, foi a busca exaustiva, que irá testar todas as possíveis
possibilidades em busca da melhor combinação.
 Partindo então para a solução do item “a”, e conhecendo a função objetivo
padrão do Matpower, que apresenta uma função de custo dependente das
potências ativa e reativa totais geradas, e estima o custo financeiro em dólar,
temos a seguinte expressão:

𝑛𝑔

𝑓(𝑃𝑔 , 𝑄𝑔 ) = ∑ 𝑓𝑃𝑖 (𝑃𝑔𝑖 ) + 𝑓𝑄𝑖 (𝑄𝑔𝑖 ) (3)

𝑖=𝑖

Com, 𝑓𝑃𝑖 e 𝑓𝑄𝑖 as funções custos de geração individual de cada gerador i,

relacionando as potências em MW ou MVAr a um valor financeiro em US$/h.

Após a execução do algoritmo, temos a seguinte solução ótima para o

problema sugerido no item “a”:

𝑸𝒔𝒉𝒖𝒏𝒕 𝜶𝒕𝒂𝒑 𝟒−𝟕 𝜶𝒕𝒂𝒑 𝟒−𝟗 𝜶𝒕𝒂𝒑 𝟓−𝟔

21MVAr 0.978 0.949 0.962

A partir desta combinação encontramos o menor custo financeiro 8078.97

US$/h. É possível observar através da Figura 1 apresentada abaixo, que tal
combinação resulta no menor valor obtido dentre as possíveis combinações.

8083

8082.5

8082

8081.5

8081

8080.5

8080

8079.5

8079

8078.5
10 20 30 40 50 60 70 80
Figura 1 – Custo financeiro.

O código utilizado para os devidos cálculos, construção do gráfico, e

detalhes da rede estão em anexo ao final da atividade.
 Para o item “b”, o manual oficial do Matpower, sugere que alterações nos
coeficientes da função custo sejam realizadas afim de torna-la uma função de
primeira ordem sem termos independentes, além de que todos os geradores
possuam custos igualitários. Portanto, temos a seguinte função:

𝑛𝑔 𝑛𝑔

𝑓(𝑃𝑔 ) = ∑ 𝑓𝑃𝑖 (𝑃𝑔𝑖 ) = ∑ 𝐶𝑖 × 𝑃𝑔𝑖 (4)

𝑖=𝑖 𝑖=𝑖

Como dito anteriormente, o custo dos geradores foram considerados iguais

a 1 para todos, logo a expressão a cima fica da seguinte forma:

𝑛𝑔

𝑓(𝑃𝑔 ) = ∑(𝑃𝑔𝑖 ) = 𝑃𝑔 (5)

𝑖=𝑖

Portanto a função objetivo é exatamente a potência ativa total gerada.

Sabendo que a perda ativa total é dada por:

𝑃𝐿 = 𝑃𝑔 − 𝑃𝑐 (6)

Considerando um cenário de carga fixa, quando minimizamos a potência

gerada, minimizamos também as perdas. É importante salientar que a função
objetivo não é exatamente a perda do sistema, necessitando o seu cálculo seja
feito posteriormente.
Após a execução do algoritmo, temos a seguinte solução ótima para o
problema sugerido no item “b”:

𝑸𝒔𝒉𝒖𝒏𝒕 𝜶𝒕𝒂𝒑 𝟒−𝟕 𝜶𝒕𝒂𝒑 𝟒−𝟗 𝜶𝒕𝒂𝒑 𝟓−𝟔

21MVAr 0.988 0.949 0.992

A partir desta combinação encontramos então o valor da função objetivo

igual a 259.49, sendo este o valor de geração ativa total que irá minimizar as
perdas. Portando para o cálculo do valor as perdas, utilizamos a Equação 6
apresentada anteriormente.

𝑃𝐿 = 259.49 − 259.00 = 0.49 𝑀𝑊 (7)

É possível ainda, analisar na Figura 2 abaixo que o menor valor de geração

ativa e também das perdas, é resultado da melhor combinação apresentada
anteriormente.
 259.58

259.57

259.56

259.55

259.54

259.53

259.52

259.51

259.5

259.49
10 20 30 40 50 60 70 80

Figura 2 – Geração de potência ativa.

Por fim, através das análises feitas para os itens “a” e “b”, é possível
observar que são conjuntos distintos, uma vez que o 𝛼𝑡𝑎𝑝 4−9 permaneceu o
mesmo para ambas as soluções. Conclui-se também que a utilização do fluxo
de potência ótimo é capaz de otimizar o sistema afim de proporcionar uma
maior economia ao mesmo.

O código utilizado para os devidos cálculos, construção do gráfico, e

detalhes da rede estão em anexo ao final da atividade.
 ANEXO ITEM “A”
%|-------------------------------------------|
%| FLUXO DE POTÊNCIA ÓTIMO |
%| SISTEMA 14 BARRAS IEEE |
%| Item A |
%| Júlia O. Fernandes |
%|-------------------------------------------|
% FLUXO DE POTÊNCIA ÓTIMO
clc % limpa a tela do workspace
close all % fecha qualquer figura previamente apresentada
clear all % limpa todas as variáveis do problema

%format shortG;
%Definindo constantes
define_constants;

%"mpopt" = compatibilidade com outras versões do matpower

mpopt = mpoption('verbose', 0, 'out.all', 0);

%Carregando o arquivo do sistema de 14 barras

mpc = loadcase('case14');

%Possíveis valores de shunt e tap

shunt = [17.0, 19.0, 21.0];
tap4_7 = [0.968, 0.978, 0.988];
tap4_9 = [0.949, 0.969, 0.989];
tap5_6 = [0.932, 0.962, 0.992];

%Cálculo da quantidade de combinações totais

P_total = length(shunt) * length(tap4_7) * length(tap4_9) *
length(tap5_6);

%Matriz custos ref. as combinações

custo = zeros(1, P_total) + Inf;

%Inicializando o contador
i = 1;

% mpc.bus = dados das barras

% mpc.branch = dados dos ramos
for j0 = 1:1:length(shunt)
mpc.bus(9, BS) = shunt(1, j0); %BS susceptância shunt do arquivo
carregado
for j1 = 1:1:length(tap4_7)
mpc.branch(8, TAP) = tap4_7(1, j1); %TAP do transformador do
arquivo carregado
for j2 = 1:1:length(tap4_9)
mpc.branch(9, TAP) = tap4_9(1, j2);
for j3 = 1:1:length(tap5_6)
mpc.branch(10, TAP) = tap5_6(1, j3);
results = runopf(mpc, mpopt);
if results.f < min(custo)
best_comb = [shunt(1, j0), tap4_7(1, j1),
tap4_9(1, j2), tap5_6(1, j3)];
end
custo(1,i) = results.f;
i = i + 1;
end
 end
end
end

%Vetor menor custo

min_custo = min(custo)

%Melhor combinação
best_comb

%Plotando o gráfico custos

plot(sort(custo))
grid; %Plano quadriculado para o gráfico
xlim([1 P_total])

mpc.bus(9, BS) = best_comb(1);

mpc.branch(8, TAP) = best_comb(2);
mpc.branch(9, TAP) = best_comb(3);
mpc.branch(10, TAP) = best_comb(4);

%Calculo do fluxo de potência ótimo runopf

results = runopf(mpc);
______________________________________________________________________
min_custo =

8079

best_comb =

21 0.978 0.949 0.962

MATPOWER Version 7.1, 08-Oct-2020 -- AC Optimal Power Flow
AC OPF formulation: polar voltages, power balance eqns
MATPOWER Interior Point Solver -- MIPS, Version 1.4, 08-Oct-2020
(using built-in linear solver)
Converged!

Converged in 0.02 seconds

Objective Function Value = 8078.97 $/hr
================================================================================
| System Summary |
================================================================================

How many? How much? P (MW) Q (MVAr)

--------------------- ------------------- ------------- -----------------
Buses 14 Total Gen Capacity 772.4 -52.0 to 148.0
Generators 5 On-line Capacity 772.4 -52.0 to 148.0
Committed Gens 5 Generation (actual) 268.2 64.6
Loads 11 Load 259.0 73.5
Fixed 11 Fixed 259.0 73.5
Dispatchable 0 Dispatchable -0.0 of -0.0 -0.0
Shunts 1 Shunt (inj) -0.0 23.2
Branches 20 Perdas (I^2 * Z) 9.23 38.72
Transformers 3 Branch Charging (inj) - 24.4
Inter-ties 0 Total Inter-tie Flow 0.0 0.0
Areas 1

Minimum Maximum
------------------------- --------------------------------
Voltage Magnitude 1.015 p.u. @ bus 3 1.060 p.u. @ bus 1
Voltage Angle -14.36 deg @ bus 14 0.00 deg @ bus 1
P Perdas (I^2*R) - 2.89 MW @ line 1-2
Q Perdas (I^2*X) - 8.83 MVAr @ line 1-2
Lambda P 36.73 $/MWh @ bus 1 41.17 $/MWh @ bus 14
Lambda Q -0.15 $/MWh @ bus 1 0.49 $/MWh @ bus 14
===============================================================================
| Bus Data |
================================================================================
Bus Voltage Generation Load Lambda($/MVA-hr)
# Mag(pu) Ang(deg) P (MW) Q (MVAr) P (MW) Q (MVAr) P Q
----- ------- -------- -------- -------- -------- -------- ------- -------
1 1.060 0.000* 194.43 0.00 - - 36.733 -0.149
2 1.039 -3.996 36.74 14.81 21.70 12.70 38.368 -
3 1.015 -9.905 28.91 22.90 94.20 19.00 40.578 -
4 1.017 -8.709 - - 47.80 -3.90 40.195 0.048
5 1.021 -7.499 - - 7.60 1.60 39.664 0.081
6 1.060 -12.771 0.00 21.15 11.20 7.50 39.740 -
7 1.050 -11.289 - - - - 40.164 0.049
8 1.060 -10.550 8.16 5.78 - - 40.163 -
9 1.051 -13.097 - - 29.50 16.60 40.149 0.081
10 1.045 -13.330 - - 9.00 5.80 40.300 0.206
11 1.049 -13.184 - - 3.50 1.80 40.149 0.171
12 1.045 -13.610 - - 6.10 1.60 40.383 0.210
13 1.041 -13.667 - - 13.50 5.80 40.574 0.332
14 1.029 -14.359 - - 14.90 5.00 41.174 0.490
-------- -------- -------- --------
Total: 268.23 64.64 259.00 73.50

================================================================================
| Branch Data |
================================================================================
Brnch From To From Bus Injection To Bus Injection Loss (I^2 * Z)
# Bus Bus P (MW) Q (MVAr) P (MW) Q (MVAr) P (MW) Q (MVAr)
----- ----- ----- -------- -------- -------- -------- -------- --------
1 1 2 129.45 -3.94 -126.56 6.95 2.891 8.83
2 1 5 64.98 3.94 -62.93 -0.80 2.052 8.47
3 2 3 55.43 0.01 -54.09 1.00 1.339 5.64
4 2 4 48.98 -2.65 -47.69 2.97 1.291 3.92
5 2 5 37.19 -2.21 -36.46 0.76 0.729 2.23
6 3 4 -11.20 2.90 11.29 -3.99 0.090 0.23
7 4 5 -50.42 5.27 50.75 -4.22 0.332 1.05
8 4 7 23.51 -4.83 -23.51 5.94 -0.000 1.11
9 4 9 15.50 4.47 -15.50 -3.21 0.000 1.26
10 5 6 41.04 2.67 -41.04 1.11 -0.000 3.78
11 6 11 5.51 3.10 -5.48 -3.03 0.034 0.07
12 6 12 7.53 2.48 -7.46 -2.34 0.069 0.14
13 6 13 16.79 6.96 -16.60 -6.58 0.195 0.38
14 7 8 -8.16 -5.62 8.16 5.78 0.000 0.16
15 7 9 31.67 -0.32 -31.67 1.32 0.000 1.00
16 9 10 7.04 4.64 -7.02 -4.58 0.020 0.05
17 9 14 10.63 3.87 -10.48 -3.56 0.147 0.31
18 10 11 -1.98 -1.22 1.98 1.23 0.004 0.01
19 12 13 1.36 0.74 -1.36 -0.73 0.005 0.00
20 13 14 4.46 1.51 -4.42 -1.44 0.035 0.07
-------- --------
Total: 9.231 38.72

================================================================================
| Voltage Constraints |
================================================================================
Bus # Vmin mu Vmin |V| Vmax Vmax mu
----- -------- ----- ----- ----- --------
1 - 0.940 1.060 1.060 633.471
6 - 0.940 1.060 1.060 42.264
8 - 0.940 1.060 1.060 29.064

================================================================================
| Generation Constraints |
================================================================================
Gen Bus Active Power Limits
# # Pmin mu Pmin Pg Pmax Pmax mu
---- ----- ------- -------- -------- -------- -------
4 6 0.260 0.00 0.00 100.00 -

Gen Bus Reactive Power Limits

# # Qmin mu Qmin Qg Qmax Qmax mu
---- ----- ------- -------- -------- -------- -------
1 1 0.149 0.00 0.00 10.00 -
 ANEXO ITEM “B”
%|-------------------------------------------|
%| FLUXO DE POTÊNCIA ÓTIMO |
%| SISTEMA 14 BARRAS IEEE |
%| Item B |
%| Júlia O. Fernandes |
%|-------------------------------------------|
% FLUXO DE POTÊNCIA ÓTIMO
clc % limpa a tela do workspace
close all % fecha qualquer figura previamente apresentada
clear all % limpa todas as variáveis do problema

%format shortG;
%Definindo constantes
define_constants;

%"mpopt" = compatibilidade com outras versões do matpower

mpopt = mpoption('verbose', 0, 'out.all', 0);

%Carregando o arquivo do sistema de 14 barras

mpc = loadcase('case14');

%funções de geração
mpc.gencusto(:, NCUSTO) = 2;
mpc.gencusto(:, CUSTO:CUSTO+1) = [1 0;1 0;1 0;1 0;1 0];

%Possíveis valores de shunt e tap

shunt = [17.0, 19.0, 21.0];
tap4_7 = [0.968, 0.978, 0.988];
tap4_9 = [0.949, 0.969, 0.989];
tap5_6 = [0.932, 0.962, 0.992];

%Cálculo da quantidade de combinações totais

P_total = length(shunt) * length(tap4_7) * length(tap4_9) *
length(tap5_6);

%Matriz perdas ref. as combinações

perdas = zeros(1, P_total) + Inf;

%Inicializando o contador
i = 1;

% mpc.bus = dados das barras

% mpc.branch = dados dos ramos
for j0 = 1:1:length(shunt)
mpc.bus(9, BS) = shunt(1, j0);
for j1 = 1:1:length(tap4_7)
mpc.branch(8, TAP) = tap4_7(1, j1);
for j2 = 1:1:length(tap4_9)
mpc.branch(9, TAP) = tap4_9(1, j2);
for j3 = 1:1:length(tap5_6)
mpc.branch(10, TAP) = tap5_6(1, j3);
results = runopf(mpc, mpopt);
if results.f < min(perdas)
best_comb = [shunt(1, j0), tap4_7(1, j1),
tap4_9(1, j2), tap5_6(1, j3)];
end
perdas(1,i) = results.f;
 i = i + 1;
end
end
end
end

%Vetor menor perda

min_perdas = min(perdas)

%Melhor combinação
best_comb

%Plotando o gráfico perdas

plot(sort(perdas))
grid; %Plano quadriculado para o gráfico
xlim([1 P_total])

mpc.bus(9, BS) = best_comb(1);

mpc.branch(8, TAP) = best_comb(2);
mpc.branch(9, TAP) = best_comb(3);
mpc.branch(10, TAP) = best_comb(4);

%Calculo do fluxo de potência ótimo runopf

results = runopf(mpc);
min_perdas =

259.49

best_comb =

21 0.988 0.949 0.992

MATPOWER Version 7.1, 08-Oct-2020 -- AC Optimal Power Flow

AC OPF formulation: polar voltages, power balance eqns
MATPOWER Interior Point Solver -- MIPS, Version 1.4, 08-Oct-2020
(using built-in linear solver)
Converged!

Converged in 0.02 seconds

Objective Function Value = 259.49 $/hr
================================================================================
| System Summary |
================================================================================

How many? How much? P (MW) Q (MVAr)

--------------------- ------------------- ------------- -----------------
Buses 14 Total Gen Capacity 772.4 -52.0 to 148.0
Generators 5 On-line Capacity 772.4 -52.0 to 148.0
Committed Gens 5 Generation (actual) 259.5 47.6
Loads 11 Load 259.0 73.5
Fixed 11 Fixed 259.0 73.5
Dispatchable 0 Dispatchable -0.0 of -0.0 -0.0
Shunts 1 Shunt (inj) -0.0 23.5
Branches 20 Perdas (I^2 * Z) 0.49 22.64
Transformers 3 Branch Charging (inj) - 25.0
Inter-ties 0 Total Inter-tie Flow 0.0 0.0
Areas 1

Minimum Maximum
------------------------- --------------------------------
Voltage Magnitude 1.034 p.u. @ bus 14 1.060 p.u. @ bus 8
Voltage Angle -0.14 deg @ bus 14 13.15 deg @ bus 8
P Perdas (I^2*R) - 0.19 MW @ line 9-14
Q Perdas (I^2*X) - 14.75 MVAr @ line 7-8
Lambda P 1.00 $/MWh @ bus 1 1.03 $/MWh @ bus 14
Lambda Q -0.00 $/MWh @ bus 1 0.01 $/MWh @ bus 14

================================================================================
| Bus Data |
================================================================================
Bus Voltage Generation Load Lambda($/MVA-hr)
# Mag(pu) Ang(deg) P (MW) Q (MVAr) P (MW) Q (MVAr) P Q
----- ------- -------- -------- -------- -------- -------- ------- -------
1 1.052 0.000* 1.23 0.01 - - 1.000 -0.002
2 1.049 0.013 21.10 0.30 21.70 12.70 1.000 -
3 1.049 0.009 94.47 16.75 94.20 19.00 1.000 -
4 1.047 0.013 - - 47.80 -3.90 1.000 0.001
5 1.049 0.011 - - 7.60 1.60 1.000 -
6 1.060 1.156 46.31 19.74 11.20 7.50 1.000 -
7 1.054 4.410 - - - - 1.000 0.001
8 1.060 13.149 96.39 10.76 - - 1.000 -
9 1.059 1.280 - - 29.50 16.60 1.000 0.001
10 1.052 0.968 - - 9.00 5.80 1.006 0.005
11 1.052 0.929 - - 3.50 1.80 1.006 0.004
12 1.046 0.354 - - 6.10 1.60 1.016 0.005
13 1.042 0.323 - - 13.50 5.80 1.019 0.008
14 1.034 -0.140 - - 14.90 5.00 1.029 0.012
-------- -------- -------- --------
Total: 259.49 47.56 259.00 73.50

================================================================================
| Branch Data |
================================================================================
Brnch From To From Bus Injection To Bus Injection Loss (I^2 * Z)
# Bus Bus P (MW) Q (MVAr) P (MW) Q (MVAr) P (MW) Q (MVAr)
----- ----- ----- -------- -------- -------- -------- -------- --------
1 1 2 1.04 1.55 -1.04 -7.37 0.004 0.01
2 1 5 0.19 -1.54 -0.19 -3.89 0.001 0.00
3 2 3 0.07 -2.30 -0.07 -2.53 0.000 0.00
4 2 4 0.37 -0.77 -0.36 -2.96 0.001 0.00
5 2 5 0.00 -1.96 -0.00 -1.85 0.000 0.00
6 3 4 0.34 0.27 -0.34 -1.68 0.001 0.00
7 4 5 -1.47 -4.87 1.48 4.88 0.003 0.01
8 4 7 -40.98 4.46 40.98 -1.30 0.000 3.16
9 4 9 -4.65 8.95 4.65 -8.49 0.000 0.46
10 5 6 -8.89 -0.74 8.89 0.92 0.000 0.18
11 6 11 3.34 2.35 -3.33 -2.32 0.014 0.03
12 6 12 7.23 2.43 -7.16 -2.29 0.064 0.13
13 6 13 15.65 6.55 -15.48 -6.22 0.170 0.33
14 7 8 -96.39 3.99 96.39 10.76 0.000 14.75
15 7 9 55.41 -2.69 -55.41 5.74 -0.000 3.05
16 9 10 9.20 5.37 -9.17 -5.28 0.032 0.09
17 9 14 12.06 4.32 -11.87 -3.93 0.186 0.40
18 10 11 0.17 -0.52 -0.17 0.52 0.000 0.00
19 12 13 1.06 0.69 -1.06 -0.69 0.003 0.00
20 13 14 3.04 1.11 -3.03 -1.07 0.016 0.03
-------- --------
Total: 0.494 22.64

================================================================================
| Voltage Constraints |
================================================================================
Bus # Vmin mu Vmin |V| Vmax Vmax mu
----- -------- ----- ----- ----- --------
6 - 0.940 1.060 1.060 0.421
8 - 0.940 1.060 1.060 0.582

================================================================================
| Generation Constraints |
================================================================================
Gen Bus Reactive Power Limits
# # Qmin mu Qmin Qg Qmax Qmax mu
---- ----- ------- -------- -------- -------- -------
1 1 0.002 0.00 0.01 10.00 -