UMA APROXIMAÇÃO LINEAR PARA MODELAGEM MATEMÁTICA DO
BALANCEAMENTO DE CARGAS EM REDES DE DISTRIBUIÇÃO
Bettoni, Luiz M. M. a
Lara, Lucas El G. de a
[email protected]
[email protected]
Passarin, Thiago A. R. a
Oliveira, Rogério P. de b
[email protected]
[email protected]
Arruda, Lúcia V. R. de a
Magatão, Leandro a
[email protected]
[email protected]
Neves Jr, Flávio a
Stebel, Sérgio L. a
[email protected]
[email protected]
RESUMO
Em sistemas de distribuição de energia elétrica é comum a prática do balanceamento de cargas
entre as fases de um circuito a fim de melhorar suas condições de operação, aumentar a vida útil
dos equipamentos ou mesmo elevar os níveis de tensão ao longo da rede. O presente artigo
apresenta um modelo matemático de otimização para o problema, desenvolvido em Programação
Linear Inteira Mista (PLIM), caracterizado por duas aproximações lineares: o balanceamento,
dado pela limitação do vetor de desequilíbrio entre fases; e o fluxo de potência, baseado em uma
iteração do método backward-forward sweep. O modelo resultante pretende reduzir custos
operacionais ao sugerir o menor conjunto de trocas de cargas capaz de atender aos índices
solicitados de balanceamento mínimo e queda máxima de tensão, respeitando as demais
considerações técnicas estipuladas. A proposta foi implementada computacionalmente e validada
em estudos de caso com dados de circuitos reais.
PALAVRAS CHAVE. Balanceamento de Cargas. Redes de Distribuição de Energia Elétrica.
Programação Linear Inteira Mista. AE - Aplicações a Energia.
ABSTRACT
In electrical energy distribution systems a common practice is the balance of the loads between
circuit phases aiming at the improvement of its operational condition, the increase of the useful
life of equipment or even the elevation of voltage levels over the network. This paper presents an
optimization mathematical model for the problem, developed in Mixed Integer Linear
Programming (MILP) and characterized by two linear approaches: the load balancing, given by a
limitation of the vector of unbalance between phases, and the power flow, based on backwardforward sweep method with one iteration. The resulting model aims to reduce operational costs
suggesting the lowest number of load changes able to comply with required parameters of
minimum balance and maximum voltage drop, following other specified technical considerations.
This proposal was computationally implemented and was validated in case studies with data from
real circuits.
KEYWORDS. Load Balancing. Electrical Energy Distribution Networks. Mixed Integer Linear
Programming. AE - Applications to Energy.
a UTFPR - Universidade Tecnológica Federal do Paraná. Av. Sete de Setembro, 3165, Curitiba/PR – 80230-901
b COPEL - Companhia Paranaense de Energia. R. Cel. Dulcídio, 800, Curitiba/PR – 80420-170
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 286
1 INTRODUÇÃO
A crescente demanda de energia elétrica exige das concessionárias constantes medidas para
expansão e adequação dos sistemas de distribuição. Estas medidas visam não só melhorar a
qualidade da energia disponibilizada para o consumidor final, mas também permitir a operação
em conformidade com as especificações dos órgãos regulamentadores (ANEEL, 2004).
Entre os diversos procedimentos técnicos utilizados para melhoria dos sistemas de
distribuição enquadra-se o balanceamento de cargas, também conhecido por faseamento de
consumidores. Este procedimento consiste em alocar uniformemente as cargas de um circuito às
fases do sistema de distribuição, reduzindo assim o desequilíbrio de corrente entre estas. Esta
redução é essencial, dado que o desequilíbrio prejudica a fase mais carregada, comprometendo o
circuito como um todo pelo surgimento de corrente no condutor neutro, quedas de tensão mais
representativas e maior carregamento dos cabos e do transformador (Souza, 2002).
Embora a formulação do balanceamento seja não-linear (Souza, 2002), apresenta-se aqui
um modelo matemático que utiliza uma formulação linear aproximada para o problema. A fim de
considerar as alterações dos níveis de tensão, o modelo inclui também um cálculo linear
aproximado de fluxo de potência, baseado em uma iteração do método backward-forward sweep
(Cheng e Shirmohammadi, 1995; Chindris et al., 2007).
O modelo matemático proposto é desenvolvido em Programação Linear Inteira Mista,
tendo por objetivo sugerir uma reconfiguração das conexões dos consumidores ligados às fases
do circuito em análise. A sugestão deve atender a parâmetros técnicos especificados previamente,
apresentando não só a solução com menor custo operacional, representado aqui pelo menor
número de trocas de fases de consumidores, mas também aquela que forneça os melhores índices
de balanceamento ou de queda de tensão dentre as alternativas de mesmo custo.
Estrutura-se o artigo da seguinte forma: a Seção 2 apresenta breve revisão de abordagens
correlatas disponíveis na literatura; a Seção 3 discorre sobre o problema do balanceamento e a
formulação aproximada desenvolvida; a Seção 4 apresenta os conceitos de fluxo de potência e a
simplificação proposta; a Seção 5 apresenta o modelo matemático desenvolvido; resultados
obtidos com o modelo são discutidos na Seção 6; considerações finais são tecidas na Seção 7.
2 REVISÃO DE LITERATURA
Diversas abordagens para o problema de balanceamento de cargas são encontrados na
literatura, normalmente associados a outras operações em redes de distribuição de energia
elétrica. Uma compilação estruturada de modelos e técnicas dedicados ao planejamento destas
redes é apresentada por Khator e Leung (1997).
Knolseisen e Coelho (2003) apresentam um sistema computacional gráfico para simulação
de reconfigurações de distribuição de cargas em uma rede de distribuição trifásica. O sistema
contempla um módulo de otimização para sugestão da configuração de trocas que minimiza as
diferenças entre as potências totais das fases. Os autores ressaltam a influência direta desta
abordagem na melhora dos níveis de tensão, redução de corrente de neutro do transformador e
conservação da integridade dos condutores e transformadores.
Díaz-Dorado e Pidre (2004) propõem o uso de programação dinâmica para localização de
subestações e atribuição de condutores em redes de distribuição. O algoritmo considera, além da
queda de tensão e das perdas de potência, a distribuição individual de cargas dos consumidores
entre as fases.
Um modelo de programação inteira combinado com busca heurística para estudos de
reconfiguração de redes é apresentado por Wu e Baran (1989). O modelo objetiva a redução de
perdas pela redistribuição equilibrada de cargas entre circuitos e ramais através da comutação de
chaves de seccionamento de rede. Kashem, Ganapathy e Jasmon (1999) apresentam uma
abordagem similar, utilizando análise gráfica para selecionar trechos de maior influência e
reduzir o esforço computacional.
Souza (2002) apresenta diversos modelos PLIM para melhoria de redes de distribuição
secundária, entre eles o balanceamento de cargas minimizando o número de operações a realizar
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 287
(trocas de fase dos consumidores). Souza conceitua o índice de desvio do balanceamento como o
vetor resultante da somatória das demandas distribuídas nas fases do circuito. Este vetor é
decomposto em componentes ortogonais para evitar a não linearidade, utilizando otimização por
etapas para minimizá-lo através de execuções recursivas. Posteriormente, Souza, Neves Jr. e
Lopes (2006) tratam o problema de melhoria dos níveis de tensão através do balanceamento
aliado ao dimensionamento de condutores, fazendo uso de um modelo de otimização multiobjetivo sob uma abordagem evolucionária baseada no critério de Pareto.
Oliveira (2008) também apresenta diversos modelos PLIM aplicados à melhoria de redes
de distribuição. Entre estes inclui-se uma aproximação linear para o balanceamento de cargas,
baseada na formulação proposta por Souza (2002), que estabelece a limitação do vetor de desvio
através de seis restrições lineares.
O presente artigo apresenta um modelo PLIM de balanceamento de cargas que estende a
aproximação de Oliveira (2008), desenvolvendo uma formulação geral com um número ajustável
de restrições lineares para limitação do vetor de desvio. Assim, a exatidão da aproximação passa
a ser função do número de restrições consideradas. Incluem-se também considerações sobre
simplificações desta formulação. Acrescentam-se ainda restrições operacionais que permitem a
limitação do número de trocas de fase, bem como do número de postes em que se pode operar
estas trocas. A proposta implementa ainda uma formulação linear simplificada de fluxo de
potência, a fim de considerar na solução a melhoria dos níveis de tensão da rede, visando mantêlos dentro dos limites regulamentados (ANEEL, 2004).
3 BALANCEAMENTO DE CARGAS
O índice de balanceamento de um circuito estima a homogeneidade da distribuição das
cargas entre as fases da rede elétrica. Estas cargas são representadas, usualmente, pela potência
complexa estática referente à demanda das unidades consumidoras em horários de pico – dados
que podem ser estimados por seu perfil característico e consumo mensal. Em redes trifásicas
estas cargas podem ser conectadas de forma monofásica, bifásica ou trifásica ao longo dos postes
do circuito. Considerando uma rede trifásica ABC, este índice pode ser dado por (Souza, 2002):
(1)
Bal=1− R/  A BC 
Onde A, B e C representam a carga total de cada fase e R representa o desequilíbrio total,
ou seja, o módulo do vetor resultante da soma vetorial das cargas das fases. Assim, sendo θf o
ângulo correspondente à fase f, R pode ser obtido por:
R=

∑
f ∈{ A , B , C }
2
 
f cos   f  
∑
f ∈{ A , B , C }

f sin   f 
2
(2)
Para um dado balanceamento mínimo (BalMin) é possível estabelecer uma região máxima
aceitável para valores de R. Reescrevendo (1) obtêm-se a seguinte inequação:
(3)
R1− BalMin  ABC 
Esta região aceitável é representada pela área interna à circunferência exibida na Figura 1.
3.1 Aproximação Linear para o Balanceamento de Cargas
Uma possível representação linear de uma circunferência pode ser feita, de maneira
aproximada, por um polígono regular circunscrito a ela. Quanto maior o numero de arestas deste
polígono, mais acurada se torna esta aproximação. Propõe-se aqui o uso desta aproximação para
delimitar uma região circular aceitável para o vetor de desequilíbrio R, conforme ilustrado na
Figura 1, permitindo estabelecer o número de equações lineares (retas) utilizadas. A formulação
correspondente, a partir de (3), sendo nr o número de retas, é dada por:
∑
*
f cos  k /nr  2 f   1− BalMin  ABC  ∀ k∈ℕ ∣knr
f ∈{ A , B , C }
(4)
A Figura 1 ilustra a aproximação linear proposta em (4) para o caso de 6 retas. Observa-se
que esta aproximação permite pequenos desvios em situações limite, conforme representado pela
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 288
Figura 1 - Limites do vetor de desequilíbrio
Figura 2 - Erro na projeção aproximada da tensão
diferença entre os vetores R e R'. É possível verificar geometricamente que estes desvios podem
atingir até 15% utilizando a aproximação linear em 6 segmentos – podendo ser reduzidos para até
3,5% utilizando 12 segmentos (Oliveira, 2008). Assim, supondo que se solicite um desequilíbrio
máximo de 10% (BalMin=90%), as soluções calculadas podem admitir vetores similares a R' com
11,5% de desequilíbrio real para aproximação com 6 segmentos ou 10,35% para 12 segmentos.
3.1.1 Simplificações para Retas Paralelas e Perpendiculares aos Eixos das Fases
Para a aproximação linear sugerida em (4), é possível observar que o aumento do número
de retas melhora a exatidão, mas também causa o aumento de complexidade da modelagem
correspondente. Expandindo o lado esquerdo da inequação (4) obtém-se:
A cos  k  Bcos   k  4 /3  C cos  k  2  /3   1− BalMin  A BC 
(5)
Nesta inequação verifica-se que as cargas atribuídas a uma fase não sofrem limitação por
retas paralelas a esta fase. Exemplificando, para nr=8 e k=2, ou seja, θk=90°, o valor atribuído à
demanda da fase A é nulo, eliminando da formulação a análise das variáveis correspondentes.
Outro caso específico surge quando uma reta é perpendicular a uma das fases, como as
retas exibidas na Figura 1. Por exemplo, para nr=6 e k=6 , ou seja, θk=0°, tem-se a partir de (5):
A−  B /2 −  C/2 1−BalMin  A BC  ⇔  3A / 2−  ABC  / 2 1− BalMin  ABC 
(6)
Esta inequação foi reescrita de forma a isolar a demanda total, uma constante do circuito,
permitindo que a fórmula seja expressada apenas em função das variáveis referentes a uma fase.
Assim, a formulação simplificada para as retas perpendiculares às fases pode ser dada por:
∣ 3 f /2− A BC  /2∣1− BalMin  A BC  ∀ f ∈ { A , B ,C }
(7)
Assim, em função da redução do número de variáveis envolvidas, é interessante adotar
aproximações em que o número de retas contemple o maior número de inequações pertencentes a
estes dois casos particulares – por exemplo: nr∈ { 6,8,12, 16, 24, 36, 48 } .
4 FLUXO DE POTÊNCIA
O fluxo de potência, ou fluxo de carga, consiste na determinação das tensões fasoriais nas
barras (nós) de uma rede elétrica, em condições determinadas de geração e carga. Em função de
sua simplicidade e adequação a redes radiais fracamente malhadas (Denis e Padilha, 1999), o
método Backward-Forward Sweep (BFS) (Cheng e Shirmohammadi, 1995) e suas variações são
amplamente adotados para cálculo de fluxo com potências constantes em redes de distribuição
trifásicas. A formulação proposta neste trabalho parte de uma variação do BFS (Chindris et al.,
2007), adotada pela simplicidade e conformidade com o padrão da concessionária COPEL,
considerando apenas uma iteração do método. O algoritmo original pode ser descrito por:
1 ) Arbitrar os fasores das tensões nodais iniciais Vif para cada nó i e fase f, normalmente
utilizando a tensão da fonte (transformador).
2 ) Na iteração k, calcular as correntes nodais Iif através das demandas de potência Sif estimadas:
k
k *
I if = S if /V if 
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
(8)
Pág. 289
3 ) Backward: Partindo dos nós extremos, calcular as correntes acumuladas Iaif e de neutro IaiN .
Sendo F o conjunto de fases e J o conjunto de nós alimentados diretamente pelo nó i :
Ia if = I if ∑ Ia jf
(9)
Ia iN = ∑ Ia if
(10)
k
k
k
j ∈J
k
k
f ∈F
4 ) Calcular as quedas de tensões nodais ΔVif . Sendo Zif a impedância do trecho que alimenta o
poste i na fase f e ZiN a impedância do neutro, caracterizando (a) como a queda da fase e (b)
como a tensão do neutro:
k
k
k
V if = 
Z if Ia if 
Z iN IaiN
a
(11)
 b
5 ) Forward: Partindo do nó fonte, calcular as tensões nodais Vif :
k
k1
k
(12)
V if =V  i −1  f −V if
6 ) Teste de convergência: se a variação máxima de tensão (ou potência calculada) entre a
iteração atual e a anterior for menor que uma tolerância desejada, o algoritmo é encerrado
com convergência do cálculo. Caso contrário, inicia-se nova iteração retornando ao passo 2.
4.1 Formulação Linear Simplificada de Fluxo de Potência
Em função da natureza iterativa do método e da não-linearidade dada pelos produtos
vetoriais, propõe-se aqui uma simplificação do fluxo de potência que assume apenas uma iteração
do método BFS e adota uma formulação linear baseada nos módulos com ângulos estimados.
Embora estas simplificações impliquem perda de exatidão, sua adoção permite linearizar o
equacionamento e reduzir o número de variáveis e equações. Explicitando, considere-se a
expansão do passo 2 do algoritmo, a partir de (8), em notação polar X =∣X ∣∢ X :
I if =∣I if∣∢ I = ∣S if∣∢ S  / ∣V f∣∢V
if
if
f

*
⇔ I if =∣S if∣/∣V f∣ ∢ V − S
f
if
(13)

Adotando a consideração usual de que as demandas i do circuito possuem o mesmo fator
de potência obtêm-se ângulos idênticos nas correntes de cada nó i de uma mesma fase. Assim,
permite-se a soma modular na equação (9):
Iaif =∣Ia if∣∢ I =∣I if ∣∢I  ∑ ∣I jf ∣∢ I
f
if
j∈ J


⇔ Iaif = ∣I if ∣∑ ∣I jf ∣ ∢I
jf
j ∈J
(14)
f
E a queda de tensão individual das fases, partindo de (11), item (a), é obtida por:
V if  a  =∣ V if  a ∣∢ V =∣Z if∣∢ Z ⋅∣Iaif ∣∢ I ⇔ V if  a =∣Z if ∣⋅∣Ia if ∣∢   Z  I
if
if  a
f
if
f

(15)
Ou, decomposto em notação retangular:
V if  a  =∣V if  a∣cos   V
 V  ∣sin  

 i∣
if a
if a 
 a1 
 V if  a 
 a2 
(16)
Sendo que, por substituição de (13) em (15), obtém-se o ângulo estimado correspondente:
(17)
 V = Z  I =Z V −S
if  a 
if
f
if
f
A parcela da tensão do neutro, (11), item (b), considerando-se idênticos os condutores
instalados em cada trecho, pode ser calculada diretamente pela soma das tensões das fases:
V iN =
Z iN Ia iN = Z if
b 
Z if Iaif
∑ Iaif = ∑ 
f ∈F
f ∈F
⇔ V iN = ∑  V if  a
f ∈F
a
(18)
Ou, decomposto em notação retangular:
V iN = ∑ ∣V if  a ∣cos  V
i ∑ ∣ V  ∣sin   
 
f ∈F
if  a 
if a
f ∈F
 b1
V if a 
(19)
 b2 
O ângulo do fasor da tensão do neutro não pode ser estimado, uma vez que depende da
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 290
soma vetorial da tensão das fases. Para evitar o cálculo vetorial utiliza-se a projeção aproximada
das quedas de tensão, desprezando o erro ilustrado pela diferença entre (i) e (ii) na Figura 2.
Visando obter o valor absoluto da tensão, desprezam-se também pequenos desvios angulares na
tensão nominal da fonte. Assim, partindo de (11) e aplicando a projeção das quedas da fase (16) e
do neutro (19) nos fasores originais da tensão, tem-se o módulo da queda de tensão nodal:
 V if  a ∣cos   V  cos  V  ∣V if  a ∣sin   V  sin  V 
∣ V if ∣ ≈ ∣

f
if  a
f
if a 
a1
a2
 ∑ ∣V ig  a ∣ cos   V
 cos  V  ∑ ∣ V ig  a∣ sin   V  sin  V 

ig  a 
g∈F

f
ig a 
g ∈F
b1
(20)
f
b2
Que, considerando o cosseno da diferença entre dois ângulos, pode ser reduzido a:
V if  a∣cos   V −V   ∑ ∣ V ig  a ∣cos  V −V 
∣ V if ∣≈∣
if  a 
f
a
g ∈F

ig a 
f
(21)
b
Finalmente, a partir de (12), calculam-se os módulos das tensões efetivas para cada nó:
∣V if∣≈∣V i −1  f∣−∣V if ∣
(22)
Embora as equações descritas envolvam operações trigonométricas, a linearidade da
formulação adotada é garantida uma vez que os cossenos e senos envolvidos representam
coeficientes constantes, calculados previamente a partir dos ângulos estimados.
5 MODELAGEM MATEMÁTICA
Em essência, o modelo matemático PLIM desenvolvido pretende apresentar uma
reconfiguração da conexão dos consumidores às fases da rede, minimizando o número de trocas
necessárias para que o circuito atenda aos índices exigidos de balanceamento mínimo (BalMin) e
queda máxima de tensão (QtMax). O modelo desenvolvido assume as seguintes premissas:
• A configuração da rede é sistema trifásico ABC.
• O custo operacional é igual para troca de qualquer fase de qualquer consumidor do circuito.
• A troca de fase de um consumidor é operada pontualmente, ou seja, são mantidas as demais
características da conexão original, como o poste associado e a distância do transformador.
• As demandas individuais são estáticas e distribuídas de forma homogênea entre as fases.
• O valor do fator de potência é padrão, igual para todas as demandas do circuito.
• As características dos condutores, como comprimento e impedância, são iguais para as fases
e o neutro em um mesmo trecho.
A nomenclatura adotada no modelo matemático é apresentada no Quadro 1.
5.1 Função Objetivo
∑ ∑  CdFP −CdFP ⋅CdFI  1−0.99⋅PrQt −Bal /100 
 
min
uf
u ∈Cons f ∈Fase
uf
uf
b
a

0.010.99⋅PrQt  QtMaxC /10  vBalM
/0.005vQtMax /0.001

c
(23)
d
Objetiva-se primeiramente minimizar o custo representado pelo número de trocas de fase
de consumidores. O primeiro termo da função objetivo (a) contabiliza a troca de fase de cada
consumidor ao identificar a sua ligação em nova fase. Em segundo plano, objetiva-se a
maximização do balanceamento (b) e a minimização da queda máxima de tensão (c), termos
ponderados por coeficientes que mantém sua funcionalidade com mínima interferência nos
demais. O valor do parâmetro binário PrQt altera estes coeficientes, permitindo optar pela
priorização (aumento do peso na função objetivo) do termo do balanceamento ou da queda de
tensão. Por fim, no último termo (d) são minimizadas as variáveis de violação sob alto custo,
utilizadas a fim de evitar a infactibilidade do modelo.
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 291
Parâmetros
Conjuntos
Fase
Cons
ConF
Post
PsEx
Cond
{1, 2, 3}
{1..nc}
⊆ Cons
{1..np}
⊆ Post
{1..nd}
BalMin
QtMax
PrQt
Fases A, B e C
Consumidores
Consumidores Fixos
Postes
Postes Extremos
Condutores (cabos)
nR
Índices
f, g, h
u
i, j
c
nmPosT
∈Fase
∈Cons, ∈ConF
∈Post, ∈PsEx
∈Cond
Fase
Consumidor
Poste
Condutor
Variáveis
Bal
vBalM
QtMaxC
vQtMax
CdFPuf
PosTFCi
CrPFif
CrAcPFif
QtPFif
QtRePFif
≥0, ≤1
≥0, ≤1
Balanceamento
Violação do Balanceamento
Mínimo
≥0, ≤1
Queda de Tensão Máxima
Calculada
≥0, ≤1
Violação da Queda de Tensão
Máxima Calculada
∈{0,1}
Condição Proposta para o
Consumidor u na Fase f
∈{0,1}
Troca de Fase de algum
Consumidor no Poste i
Corrente individual do Poste i na Fase f
Corrente Acumulada do Poste i na Fase f
Queda Tensão individual do Poste i Fase f
Queda Tensão Relativa do Poste i na Fase f
nMxTro
nMnTro
TenTra
FatPt
DemTot
DemConu
PosConu
NatFCu
CdFIuf
CdInsij
LAdji
DistPij
ImpCdc
aICc
aTFf
Balanceamento Mínimo Exigido (%)
Queda de Tensão Máxima Exigida (%)
∈{0,1}
Seleção de Prioridade para
Queda de Tensão
≥6
Número de Retas da
Linearização do
Balanceamento
≥1
Número Máximo de Postes a
operar Trocas
≥0
Número Máximo de Trocas
≥0
Número Mínimo de Trocas
Tensão Fornecida pelo Transformador (V)
Fator de Potência das Cargas (%)
Demanda Total do Circuito (kVA)
Demanda do Consumidor u (kVA)
∈Post
Poste em que o Consumidor u
está ligado
∈{1,2,3}
Natureza de Faseamento do
Consumidor u
∈{0,1}
Condição de Inicial do
Consumidor u na Fase f
∈Cond
Condutor Instalado entre os
Postes i e j
∈Post
Lista de Adjacência: poste que
alimenta diretamente o poste i
Distância entre os Postes i e j (m)
Impedância do Condutor c (Ω/m)
Ângulo da Impedância do Condutor c (rad)
Ângulo do fasor de Tensão da Fase f (rad)
Quadro 1 - Nomenclatura adotada na modelagem matemática
5.2 Restrições
Balanceamento de Cargas
O vetor de desequilíbrio deve estar contido nos limites definidos pelas equações da
aproximação linear do balanceamento, conforme exposto na Seção 3.1. A partir de (4) tem-se:
∑ ∑
u ∈Cons f ∈Fase

CdFP uf
 
DemCon u
2
cos k
aTF f  1−Bal  DemTot
NatFC u
nR

∀ k∈ { 1nR }
(24)
No entanto, a simplificação sugerida em (7) pode ser utilizada em retas paralelas aos eixos
das fases, ou seja, quando k(6/nR)∈{1..6}. Nestes casos, sendo “÷” o operador “resto de divisão
inteira”, a formulação é substituída por:
 k nR6 
−1 




DemConu DemTot
3
CdFP uf
−
 1− Bal  DemTot
∑
2 u ∈Cons
NatFC u
2
(25)
∀ k ∈ {1nR } , ∀ f ∈Fase∣ f =1−k÷3
Relaxação do Balanceamento
O balanceamento do circuito deve ser maior que o balanceamento mínimo, permitindo-se o
relaxamento pela variável de violação.
BalvBalM  BalMin
(26)
Conservação da Natureza dos Consumidores
A natureza do faseamento dos consumidores tem que ser mantida, ou seja, consumidores nfásicos devem continuar sendo n-fásicos.
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 292
∑
f ∈Fase
CdFP uf = NatConu ∀ u∈Cons
(27)
Conservação dos Consumidores Fixos
Para os consumidores fixos, ou seja, consumidores para os quais não se permite que seja
modificada a conexão à rede, a condição de faseamento proposta deve ser igual à inicial.
CdFP uf =CdFI uf
∀ u∈ConF , ∀ f ∈ Fase
(28)
Limitação do Número de Postes com Trocas
Visando concentrar as atividades de manutenção em campo, o número de postes em que
são operadas trocas não deve ser superior ao número máximo definido.
∑
i∈Post
PosTFC i nmPosT
(29)
Os postes que sofrem trocas são identificados pela variável binária auxiliar PosTFCi. Duas
restrições, elaboradas de modo a representar implicações lógicas (Magatão, 2005), são utilizadas
para definir seu valor: se não houver trocas em i então PosTFCi assume valor 0.
∑ ∑ CdFP uf −CdFPuf ⋅CdFI uf −PosTFC i 0
u
∀ i∈ Post ,∀ u∈Cons∣PosConu=i
f ∈ Fase
(30)
E, sendo ε um número relativamente pequeno (como 1/nc, a fim de garantir que o termo
associado seja menor que 1): se houver uma ou mais trocas em i então PosTFCi assume valor 1.
∑ ∑ CdFP uf −CdFPuf ⋅CdFI uf ⋅−PosTFC i 0
u
∀ i∈ Post ,∀ u∈Cons∣PosConu =i
f ∈ Fase
(31)
Limitação do Número de Trocas
Se especificados, nMxTro define o limite máximo para o número de trocas e nMnTro o
limite mínimo. A combinação de ambos permite fixar um intervalo para o número de trocas.
∑ ∑  CdFP uf −CdFP uf⋅CdFI uf nMxTro
(32)
∑ ∑  CdFP uf −CdFP uf⋅CdFI uf nMnTro
(33)
u ∈Cons f ∈Fase
u ∈Cons f ∈Fase
Fluxo de Potência Simplificado
As quedas de tensão são calculadas usando o fluxo de potência simplificado exposto na
Seção 4.1. A partir de (13) obtêm-se os módulos das correntes individuais por poste e fase:

CrPF if = ∑ CdFPuf
u

DemConu 1000
NatFC u TenTra
∀ i∈ Post ,∀ f ∈ Fase ,∀ u∈Cons∣PosCon u=i
(34)
E, a partir de (14), obtém-se o módulo da corrente acumulada para cada poste e fase:
CrAcPF if =CrPF if ∑ CrAcPF jf
∀ i∈Post , ∀ f ∈Fase ,∀ j ∈Post∣LAdj j =i
(35)
j
De (15) obtém-se o módulo da queda de tensão individual por poste e fase:
QtPF if =CrAcPF if⋅DistP ji⋅ImpCd CdIns
ji
∀ i∈Post , ∀ f ∈Fase ,∀ j ∈Post∣ j= LAdj i
(36)
Finalmente, partindo de (21), utilizando o princípio de (22) para obter a queda normalizada
pela tensão do transformador, tem-se a queda de tensão relativa para cada poste e fase:
QtRePF if =QtRePF jf 
1
QtPF if cos  if −aTF f  ∑ QtPF ig cos  ig −aTF f 
TenTra
g∈ Fase


(37)
∀ i∈ Post , ∀ f ∈Fase ,∀ j ∈Post∣ j= LAdj i
Sendo que os ângulos estimados α são parâmetros dados por (17):
ih =aIC CdIns aTF h −acos  FatPt  ∀ i∈ Post ,∀ h∈Fase ,∀ j ∈Post∣ j= LAdj i
ji
(38)
Identificação da Queda Máxima Calculada
O valor da queda de tensão máxima calculada corresponde à maior dentre as quedas dos
postes/fases. O equacionamento é fundamentado na minimização de QtMaxC na função objetivo.
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 293
QtRePF if QtMaxC
∀ i ∈ Post , ∀ f ∈ Fase
(39)
Limitação da Queda de Tensão
A queda de tensão máxima calculada não deve superar o valor da queda máxima exigida,
permitindo o relaxamento pela variável de violação correspondente.
(40)
QtMaxCvQtMaxQtMax
5.3 Implementação
A implementação computacional do modelo foi realizada em GMPL, linguagem própria do
solver livre adotado, GLPK (2009). A execução faz uso de um arquivo de dados contendo os
parâmetros desejados e as informações do circuito em análise. Após a resolução apresenta-se o
relatório de saída do modelo, conforme exemplo ilustrado na Figura 3.
Execucao Priorizando Balanceamento Maximo
Queda de Tensao Maxima Final:
5.726 %
Balanceamento Inicial, nao Linear: 84.848 %
Balanceamento Solicitado:
90.000 %
Balanceamento Final, Linear:
91.077 %
Consum. A B C
Codigo
Poste Codigo
===========================================
52
D L
3291582
10 351419553
55
D
L
3498039
10 351419553
132
D L
3429610
22 351419992
===========================================
L eg en da : D= De sl ig am en to L =L ig aç ão
Limite de Postes com Trocas:
Numero de Trocas Sugeridas:
Custo Total Estimado (R$):
2
3
135.00
Queda Tensao
Fase A
Fase B
Fase C
===========================================
Poste 1
0.000
0.000
0.000
Poste 2
1.106
0.868
0.665
Poste 3
1.439
1.243
0.882
...
...
...
...
Poste 21
4.796
2.952
[ 2.687 ]
Poste 22
[ 5.726 ]
3.445
2.527
===========================================
Figura 3 - Exemplo sumarizado de relatório de saída do modelo
Foram implementadas ainda considerações adicionais: limites e valores padrão para alguns
parâmetros; arredondamento de coeficientes a fim de eliminar termos desprezíveis; geração
seletiva de restrições incorporada ao próprio modelo. Esta última reduz a complexidade do
problema ao considerar apenas formulações referentes aos parâmetros de entrada especificados.
6 SIMULAÇÃO E RESULTADOS
A validação do modelo proposto, bem como de sua implementação computacional, foi
realizada através de estudos de caso com diversos circuitos reais. Para execução dos testes foram
utilizados dados obtidos a partir da base geo-referenciada da COPEL – Companhia Paranaense de
Energia. Os circuitos mencionados nos testes seguintes são apresentados na Tabela 1.
A simulação foi realizada em três etapas: validação da linearização do balanceamento,
validação do fluxo de potência simplificado e análise das trocas sugeridas sob diversos cenários.
Tabela 1 - Informações gerais dos circuitos utilizados
Circuito
1
2
3
4
5
6
7
8
Bal. Inicial (%)
95,987
91,834
93,906
84,848
87,440
98,723
79,811
84,738
Nº Cons
167
129
126
133
80
141
48
32
Nº ConF
42
42
31
22
20
39
20
12
Nº Post
42
42
31
22
20
39
20
12
DemTot (kVA)
84,1
61,2
94,5
59,3
42,8
66,8
17,6
18,5
Q.T. Máxima (%)
8,947
9,091
10,515
8,085
5,731
5,454
2,117
0,735
A validação do balanceamento foi realizada comparando os valores do balanceamento
inicial, calculado usando a formulação não linear (equações (1) e (2)), aos fornecidos pelo
modelo linearizado com variados números de retas. Um sumário dos resultados selecionados é
exibido na Tabela 2. Observa-se que o erro apresentado pelos valores do balanceamento
linearizado não foi significativo. Percebe-se também que o aumento do número de retas implica
melhor aproximação quando o vetor de desequilíbrio encontra-se nas áreas próximas às
intersecções das retas – fato bem representado pelos circuitos 3 e 7 na aproximação com 12 retas
e por praticamente todos os circuitos na aproximação com 36 retas.
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 294
Tabela 2 - Sumário da validação do balanceamento
Circuito
1
2
3
4
5
6
7
8
I. Balanceamento
Não Linear (%)
95,987
91,834
93,906
84,848
87,440
98,723
79,811
84,738
Bal. Linear – 6 retas
II. Bal (%)
Erro (II-I)
96,009
0,022
91,940
0,106
94,224
0,318
84,858
0,010
87,447
0,007
98,734
0,011
80,187
0,376
84,740
0,002
Bal. Linear – 12 retas
III. Bal (%)
Erro (III-I)
96,009
0,022
91,940
0,106
94,027
0,121
84,858
0,010
87,447
0,007
98,734
0,011
79,960
0,149
84,740
0,002
Bal. Linear – 36 retas
IV. Bal (%)
Erro (IV-I)
95,997
0,010
91,834
0,000
93,908
0,002
84,858
0,010
87,447
0,007
98,724
0,001
79,815
0,004
84,740
0,002
Para a validação do fluxo de potência foram executados testes comparando os níveis de
queda de tensão obtidos pelo modelo proposto aos resultados de duas execuções distintas da
ferramenta de fluxo de potência adotada como referência (também baseada no método BFS),
proprietária da concessionária COPEL (Oliveira, 2008): a primeira com apenas uma iteração e a
segunda com múltiplas iterações (execução completa). Parte dos resultados é exibida na Tabela 3,
apresentando apenas os valores máximos para os circuitos com quedas mais significativas.
Tabela 3 - Sumário da validação do fluxo de potência
Circ.
Fase
1
2
3
4
5
6
1
Queda de Tensão Máxima (%)
I. Modelo Proposto II. COPEL 1 iteração III. COPEL Completo
A
B
C
A
B
C
A
B
C
5,140 8,947 5,797 5,129 8,941 5,795 5,397 9,798 6,106
2,721 5,456 9,091 2,717 5,416 9,083 2,654 5,708 10,053
7,726 10,566 6,564 7,683 10,554 6,559 8,309 11,952 7,032
8,085 3,978 2,536 8,080 3,976 2,520 8,849 4,145 2,518
2,815 5,731 3,774 2,813 5,728 3,771 2,852 6,078 3,868
5,454 3,921 5,292 5,436 3,908 5,290 5,766 3,986 5,573
Maior Erro Encontrado no Modelo1
P/ COPEL 1 iteração P/ COPEL Completo
A
B
C
A
B
C
0,026 0,010 0,016 -0,260 -0,851 -0,309
0,013 0,041 0,008 0,130 -0,258 -0,962
0,042 0,012 0,060 -0,583 -1,386 -0,468
0,020 0,015 0,040 -0,764 -0,167 0,115
0,003 0,003 0,006 -0,040 -0,347 -0,094
0,016 0,012 0,004 -0,312 -0,112 -0,281
Apresentam-se os maiores erros entre os valores tensão de todos os postes, não se restringindo apenas aos sumarizados em II e III.
Os erros de queda de tensão encontrados no modelo, comparando este ao fluxo de potência
de referência sob uma iteração, são consideravelmente reduzidos, indicando a validade das
simplificações adotadas na formulação. Contudo, como esperado, erros maiores surgem em
comparação ao fluxo de várias iterações, agravados principalmente nos pontos de maior queda.
A análise das trocas sugeridas pelo modelo de balanceamento foi realizada para 6
diferentes especificações (testes de A a F), considerando diferentes valores para os parâmetros de
execução. As soluções para circuitos selecionados são apresentadas na Tabela 4, incluindo dados
referentes ao custo computacional da resolução. A validação de cada resposta foi realizada pela
análise e simulação das implicações das trocas sugeridas. O teste A apresenta o estado inicial do
circuito, apenas para comparação.
Os testes B e C representam as exigências normalmente praticadas pelas concessionárias
para um circuito secundário: balanceamento mínimo de 90% e queda máxima de tensão de 8%.
Observa-se que para todos os circuitos foi sugerida ao menos uma troca a fim de atender a estas
exigências. O teste C difere do B apenas por adotar a priorização da queda de tensão (PrQt=1).
Ambos os testes sugerem o mesmo número de trocas para cada circuito, mas diferenciam a
solução por fornecer resultados com melhores valores para o respectivo critério priorizado.
Nos testes D, E e F as exigências são mais acentuadas, acarretando um maior número de
trocas. No caso do circuito 3 estes testes apresentam soluções com violação, uma vez que não é
possível obter solução com queda de tensão inferior à solicitação de 6%. No teste E as trocas são
limitadas a apenas 2 postes. Para o circuito 3 esta limitação implica em violação ainda maior da
queda de tensão, enquanto para o circuito 4 as exigências são atendidas sob o custo de mais uma
troca. O teste F diferencia-se do D pela linearização de balanceamento com apenas 6 retas,
ilustrando a diferença que o número de retas adotado na aproximação pode causar nos índices de
balanceamento (circuitos 1 e 2) e, eventualmente, de queda de tensão (circuito 4). A simplificação
obtida pela presença de retas perpendiculares na aproximação pode ser induzida ao comparar os
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 295
custos computacionais dos testes D e F – efeito corroborado por diversos outros ensaios
realizados, principalmente em circuitos de maior dimensão.
Tabela 4 - Sumário da execução1 dos testes com diferentes especificações
Teste
BalMin/QtMax/nR
Circ.
Parâmetros adicionais
1 Bal(%) / QtMaxC(%)
Trocas / Tempo(s)
Mem(Mb) / Iterações
2 Bal(%) / QtMaxC(%)
Trocas / Tempo(s)
Mem(Mb) / Iterações
3 Bal(%) / QtMaxC(%)
Trocas / Tempo(s)
Mem(Mb) / Iterações
4 Bal(%) / QtMaxC(%)
Trocas / Tempo(s)
Mem(Mb) / Iterações
5 Bal(%) / QtMaxC(%)
Trocas / Tempo(s)
Mem(Mb) / Iterações
7 Bal(%) / QtMaxC(%)
Trocas / Tempo(s)
Mem(Mb) / Iterações
1
A
- / - / 12
96,009 / 8,947
0 / 0,2
3 / 353
91,94 / 9,091
0 / 0,1
2,7 / 274
94,224 / 10,515
0 / 0,1
2,2 / 221
84,858 / 8,085
0 / 0,1
2 / 247
87,447 / 5,731
0/0
1,5 / 162
79,960 / 2,117
0/0
1,3 / 78
B
90% / 8% / 12
96,009 / 6,588
1 / 1,4
3,2 / 871
94,463 / 6,257
1 / 0,2
2,8 / 416
96,729 / 7,991
2 / 3,1
2,7 / 1546
91,405 / 5,71
3 / 0,3
2,1 / 316
91,269 / 5,731
1 / 0,1
1,5 / 180
95,663 / 2,472
2 / 0,4
1,4 / 240
C
90% / 8% / 12
PrQt:1
95,673 / 5,948
1 / 0,4
3,2 / 414
94,463 / 6,257
1 / 0,2
2,8 / 329
96,116 / 7,952
2 / 2,3
2,7 / 1846
91,001 / 4,951
3 / 0,8
2,3 / 599
90,428 / 4,174
1 / 0,1
1,5 / 182
94,501 / 1,898
2 / 0,4
1,4 / 342
D
94% / 6% / 12
95,673 / 5,948
1 / 0,3
3,1 / 440
96,747 / 6,075
2 / 2,8
3 / 716
98,588 / 6,748
6 / 786,3
75,7 / 307665
93,554 / 5,496
4 / 0,2
2,1 / 322
94,776 / 4,513
2 / 0,1
1,5 / 183
95,663 / 2,472
2 / 0,3
1,4 / 190
E
94% / 6% / 12
nmPosT:2
95,673 / 5,948
1 / 0,2
3,2 / 173
96,747 / 6,075
2 / 3,8
3,2 / 570
98,305 / 7,094
4 / 19,7
3,2 / 5097
94,692 / 4,815
5 / 2,5
2,8 / 1053
94,776 / 4,513
2/0
1,6 / 107
95,663 / 2,472
2 / 0,4
1,5 / 305
F
94% / 6% / 6
96,009 / 5,948
1 / 0,2
2,8 / 552
97,109 / 6,075
2 / 1,9
2,8 / 927
98,588 / 6,748
6 / 379,8
36,1 / 150237
93,554 / 5,491
4 / 0,2
1,9 / 429
94,776 / 4,513
2 / 0,1
1,4 / 229
95,663 / 2,472
2 / 0,5
1,4 / 422
Execução em IBM/PC, AMD Athlon64 2,2 GHz, 2 Gb RAM DDR400, MS Windows XP SP3 32 bits, GLPK v4.33.
98
Queda Máxima de
Tensão Calculada (% )
Balanceamento (% )
Por fim, a Figura 4 apresenta o comparativo dos resultados de 12 execuções do modelo
para um mesmo circuito arbitrado. Cada execução foi realizada estabelecendo apenas o número
fixo de trocas (progressivamente, de 1 a 6), priorizando ora o índice de balanceamento (PrQt=0),
ora o de queda de tensão (PrQt=1). Nota-se que, quanto maior o número de trocas sugerido, mais
significativa torna-se a diferença entre priorizar ou não um índice. Contudo, em âmbito geral as
melhorias nos índices do circuito tendem a ser menos representativas à cada troca adicional. Isto
ocorre pela aproximação dos resultados aos limites dos índices e pelo aumento de combinações
possíveis – fato que também implicou custos computacionais mais elevados nos testes.
96
94
92
90
88
86
84
0
1
2
3
4
PrQt=0
PrQt=1
5
6 Trocas
10
9
8
7
6
5
4
3
0
1
2
PrQt=0
3
4
5
6 Trocas
PrQt=1
Figura 4 - Resultados do circuito 4 para diversas trocas alternado a priorização da queda de tensão
7 CONCLUSÕES
No presente artigo foi desenvolvido um modelo matemático PLIM que aborda o problema
de balanceamento de cargas minimizando o número necessário de trocas de fases de
consumidores para atender a índices exigidos de balanceamento mínimo e queda máxima de
tensão. A fim de permitir o cálculo destes índices no modelo, apresentam-se aproximações
lineares para a formulação do balanceamento de cargas e do fluxo de potência – soluções de
escopo abrangente, passíveis de aplicação em problemas similares. Simplificações adicionais
foram implementadas, permitindo reduzir a complexidade do modelo matemático. Implementamse ainda considerações sobre priorização e relaxação dos índices exigidos, bem como restrições
operacionais que permitem aproximar as soluções obtidas às necessidades práticas.
Os estudos de caso demonstraram a validade das formulações adotadas, apresentando
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 296
desvios considerados adequados para atender aos propósitos do modelo. Os resultados do fluxo
de potência simplificado ficaram muito próximos dos de referência sob uma iteração, embora as
diferenças para os de várias iterações sugiram considerar a adoção de medidas corretivas como,
por exemplo, a solicitação de valores mais estritos (menores) de queda máxima de tensão. No
entanto, estas diferenças encontradas no fluxo de potência, bem como os erros da aproximação
linear do balanceamento frente à formulação original, são pouco significativos se considerada a
imprecisão inerente aos valores estáticos e estimados de demanda dos consumidores.
Os testes realizados demonstraram que simulações exigindo índices demasiadamente
rigorosos (alto balanceamento e baixa queda de tensão), principalmente nos circuitos com maior
número de consumidores, apresentam custo computacional elevado, demonstrando o aumento
progressivo da complexidade do problema combinatório quando o número de trocas necessárias
aumenta significativamente. Contudo, para os índices convencionais de balanceamento e queda
de tensão, ou seja, valores normalmente praticados pelas concessionárias, obtém-se soluções com
custo computacional baixo, favorecendo a adoção do modelo desenvolvido em ferramentas de
sugestão que visem a prática do balanceamento de cargas em redes de distribuição.
AGRADECIMENTOS
Apoio financeiro da Agência Nacional de Energia Elétrica (ANEEL) por meio do Projeto
de Pesquisa e Desenvolvimento 2866-018/2005.
BIBLIOGRAFIA
ANEEL, Agência Nacional de Energia Elétrica. Resolução Nº 505, de 26 de novembro de 2001. Diário
Oficial da União, Vol. 141, No. 147, pp. 73. Brasília, 2 ago. 2004.
Cheng C.S.; Shirmohammadi D. A Three-Phase Power Flow Method For Real-Time Distribution System
Analysis. IEEE Transactions on Power Systems, Vol. 10, No. 2, pp. 671-679, 1995.
Chindris, M.; Sudria i Anderu, A.; Bud, C.; Tomoiaga, B. The Load Flow Calculation In Unbalanced
Radial Electric Networks With Distributed Generation. Electrical Power Quality and Utilisation, 9th
International Conference on, pp. 1-5, 2007.
Denis, I. F. E. D.; Padilha, A. Three-Phase Fast Decoupled Power Flow. Electric Power Engineering,
International Conference on, pp. 283, 1999.
Díaz-Dorado, E.; Pidre, J. C. Optimal Planning of Unbalanced Networks Using Dynamic Programming
Optimization. IEEE Transactions on Power Systems, Vol. 19, No. 4, pp. 2077-2085, 2004.
GLPK (GNU Linear Programming Kit) v. 4.37. GNU Project, (http://www.gnu.org/software/glpk), 2009.
Kashem, M.A.; Ganapathy, V; Jasmon, G.B. Network Reconfiguration for Load Balancing in
Distribution Networks. IEE Proceedings Generation, Transmission and Distribution, Vol. 146, No. 6, pp.
563-567, 1999.
Khator, S. K.; Leung, L. C. Power Distribution Planning: A Review of Models and Issues. IEEE
Transactions on Power Systems, Vol. 12, No. 3, pp. 1151-1159, 1997.
Knolseisen, A. B.; Coelho, J. Balanceamento de Cargas em Sistemas de Distribuição de Baixa Tensão.
CTAI - Revista de Automação e Tecnologia da Informação. Vol. 2, No. 1, pp. 6-11, Florianópolis, 2003.
Magatão, L. Mixed Integer Linear Programming and Constraint Logic Programming: Towards a Unified
Modeling Framework. Tese de Doutorado, Centro Federal de Educação Tecnológica do Paraná. Curitiba,
2005.
Oliveira, R. P. de. Ferramenta de Otimização Para Projetos de Melhoria de Redes de Distribuição
Secundária de Energia Elétrica. Dissertação de Mestrado, Universidade Tecnológica Federal do Paraná,
Curitiba, 2008.
Souza, A. A. A. Modelos Matemáticos para Auxílio ao Planejamento e Projeto de Redes Secundárias de
Distribuição. Dissertação de Mestrado, Centro Federal de Educação Tecnológica do Paraná. Curitiba,
2002.
Souza, A. A. A.; Neves Jr., F.; Lopes, H. S. Sistema de Avaliação da Rede Secundária de Distribuição
Utilizando Algoritmos Genéticos. Espaço Energia, No. 5, pp. 34-41, 2006.
Wu, F. F.; Baran, M. E. Network reconfiguration in Distribution Systems for Loss Reduction and Load
Balancing. IEEE Transactions on Power Delivery, Vol. 4, No. 2, pp. 1401-1407, 1989.
XLI SBPO 2009 - Pesquisa Operacional na Gestão do Conhecimento
Pág. 297