Resumos

SISTEMAS ELÉTRICOS DE POTÊNCIA

Técnicas de parametrização global para o fluxo de carga continuado

Global parameterization techniques to the continuation power flow

Alfredo Bonini Neto; Dilson Amancio Alves

UNESP - Departamento de Engenharia Elétrica, Av. Brasil, 56, 15385-000 Ilha Solteira - Brasil

RESUMO

Este trabalho analisa a utilização de técnicas de parametrização global para o fluxo de carga continuado. Essas técnicas são consideradas inadequadas para a obtenção da margem de carregamento de sistemas com problemas de estabilidade de tensão com características fortemente locais. Isto se deve ao fato de que no ponto de máximo carregamento a singularidade da matriz Jacobiana do método de parametrização global coincide com a da matriz Jacobiana do fluxo de carga. Nesses casos, a parametrização local é considerada como a única forma de se eliminar a singularidade. Entretanto, este trabalho mostra que a singularidade também pode ser eficientemente eliminada não só para estes sistemas, mas para qualquer outro, através de uma nova técnica de parametrização. A técnica utiliza a equação de uma reta que passa através de um ponto no plano determinado pelas variáveis fator de carregamento e a somatória das magnitudes, ou dos ângulos, das tensões nodais de todas as barras do sistema, que são as variáveis comumente usadas pelas técnicas de parametrização global. Os resultados obtidos para diversos sistemas confirmam o aumento da eficiência dos métodos propostos e mostram sua viabilidade para aplicações no planejamento da operação nos atuais sistemas de gerenciamento de energia.

Palavras-chave: Método da Continuação, Fluxo de Carga, Técnicas de Parametrização, Ponto de Máximo Carregamento.

ABSTRACT

This paper presents an analysis of the use of global parameterization techniques to the continuation power flow. Those techniques are considered inadequate for computation of the loading margin of power systems characterized by strong local static voltage stability. In such systems, at maximum loading point, the singularity of the Jacobian matrices of global parameterization techniques coincide with the one of the power flow Jacobian matrix. In those cases, the local parameterization is considered as the only way to overcome the singularity. However, this paper shows that this kind of singularity can be efficiently eliminated not only for these systems, but also for all others, by a new parameterization technique. This technique uses the addition of a line equation, which passes through a point in the plane determined by the sum of all the bus voltage magnitudes, or angles, and loading factor variables, that are variables commonly used by global parameterization techniques. The obtained results for several systems confirm the efficiency increased of the proposed methods and show its viability for applications in the operating planning in a modern energy management system.

Keywords: Continuation Method, Load Flow, Parameterization Techniques, Maximum Loading Point.

1 INTRODUÇÃO

Nos estudos relacionados com a estabilidade estática de tensão de sistemas elétricos de potência (SEP), as empresas do setor elétrico internacional (WSCC, 1998) e nacional (ONS, 2002; FTCT, 1999) recomendam o levantamento do perfil de tensão das barras em função de seu carregamento (curvas P-V e Q-V). Essas curvas possibilitam avaliação qualitativa de diversas condições operativas do sistema em diferentes condições de carregamento e contingências. Estes perfis são utilizados, entre outras, para: determinar os limites de transferência de potência entre as áreas de um sistema; ajustar margens; observar o comportamento das tensões das barras do sistema em análise; e comparar estratégias de planejamento visando à proposição adequada de ampliações e reforços da rede com o intuito de se evitar o corte de carga. Nestes estudos o levantamento da curva P-V é considerado o mais adequado para a determinação das margens de estabilidade, enquanto que as metodologias Q-V e análise modal são consideradas como complementares. As curvas Q-V possibilitam a avaliação das margens de potência reativa enquanto que a análise modal fornece informações para a determinação de locais mais apropriados para o reforço do sistema. Por se tratar de uma análise linear, deve-se proceder a análise modal num ponto de operação o mais próximo possível do ponto de máximo carregamento (PMC). Assim, um dos principais objetivos desses estudos é a obtenção do PMC do sistema. A determinação exata do PMC é importante para identificar se ocorrem mudanças significativas de informação. Apesar dos métodos de FC convencionais possibilitarem o cálculo de pontos de operação muito próximos ao PMC, a sua determinação exata fica limitada. Isso ocorre porque além dos problemas relacionados com a singularidade da matriz Jacobiana (J) nesse ponto, a solução também dependerá das características comuns aos processos de solução de equações algébricas não-lineares, tais como do método utilizado na resolução, da existência da solução, das múltiplas soluções existentes e da estimativa inicial. Assim, sempre será necessário ponderar se os problemas de não convergência são devidos a problemas numéricos ou a limitações físicas do sistema. Em geral, as diferenças não são óbvias.

Métodos diretos foram propostos e são assim denominados por possibilitarem o cálculo direto do PMC, sem a determinação dos demais pontos existentes entre o caso base e o PMC (Alvarado e Jung 1988; Cañizares et al., 1992; Cañizares e Alvarado, 1993; Ferreira et al., 2004). O método direto foi proposto inicialmente em (Alvarado e Jung 1988) e posteriormente adaptado em (Cañizares et al., 1992; Cañizares e Alvarado, 1993). Este método é baseado na teoria da bifurcação e na singularidade da matriz Jacobiana, e é usado para detectar uma bifurcação sela-nó para certa direção de crescimento de carga. Uma bifurcação do tipo sela-nó (ou em particular PMC) é caracterizada por um Jacobiano que tenha um autovalor simples e único igual à zero, com os autovetores direito (v) e esquerdo (w) não nulos. Neste método o sistema de equações não-lineares resultante é duas vezes maior do que o sistema de equações do fluxo de carga convencional. Estes métodos têm-se mostrado eficientes na localização do ponto de PMC quando o sistema se encontra muito perto deste ponto, e quando o sistema se encontra longe do PMC, e especialmente quando todos os limites e controles são considerados, o método tende a falhar (Cañizares e Alvarado, 1993). Os autovalores e autovetores são responsáveis pelos problemas de convergência quando o sistema se encontra longe do PMC, uma vez que mudam muito conforme o sistema se aproxima do ponto de colapso, principalmente quando os limites são encontrados ao longo da trajetória de solução. Se os limites são desconsiderados, o processo de solução é rápido, porém não realista. Uma alternativa para superar a dificuldade e obter uma boa estimativa para v ou w é calcular esses vetores a partir de um ponto de operação inicial, ou seja, obter o autovetor correspondente ao menor autovalor de J no caso base. De acordo com a referência (Cañizares e Alvarado, 1993), das experiências obtidas com o uso deste método, constatou-se que o uso das equações correspondentes ao autovetor esquerdo produzem resultados melhores que os alcançados com o uso das equações correspondentes ao autovetor direito. O fluxo de potência continuado obtém todo o perfil de tensão das barras enquanto que o método direto concentra-se apenas na obtenção do PMC do SEP. Assim, o método direto não provê informações a respeito do mecanismo do colapso de tensão, nem do modo como o sistema se aproxima do PMC, mas apenas possibilita o cálculo do mesmo.

Nos métodos de fluxo de carga continuados (FCC) as equações do FC são reformuladas visando à eliminação da singularidade da matriz Jacobiana no PMC e, consequentemente, dos problemas numéricos que ocorrem em torno deste. Em geral isto é efetuado por meio da adição de equações parametrizadas (Cañizares et al., 1992; Alves et al., 2003, Ajjarapu e Christy, 1992; Chiang et al., 1995; Garbelini et al., 2006, Seydel, 1994; Leite e Da Costa, 2003). Assim, o traçado completo do perfil de tensão pode ser efetuado variando automaticamente o valor de um determinado parâmetro do sistema. A diferença entre os métodos de FCC está no modo como o novo parâmetro é escolhido e em como a singularidade é eliminada. As técnicas de parametrização mais utilizadas pelos FCC para eliminar a singularidade de Jacobiana são a geométrica (Cañizares et al., 1992; Chiang et al., 1995; Li e Chiang, 2008) e a local (Ajjarapu e Christy, 1992).

Iba et al., (1991) apresentaram uma técnica para contornar a singularidade de Jacobiana sem a necessidade de parametrização. A técnica consiste em definir um vetor perpendicular ao vetor tangente ao ponto da curva, e que passe pelos pontos previsto subseqüente e um outro que se encontra sobre a curva da trajetória de soluções. Esta técnica exige um bom controle de passo nas proximidades do PMC (Cañizares et al., 1992).

Chiang et al. (1995) propuseram uma parametrização onde o comprimento de arco (s) é utilizado como parâmetro. Entretanto, como o sistema de equação formado no passo preditor é não linear, a sua solução exige um método especial, o que pode implicar num tempo computacional muito alto. A opção sugerida pelos autores foi o uso do preditor secante logo após a obtenção de dois pontos da curva. De acordo com os autores esta técnica de parametrização é mais robusta, possibilitando que sejam dados passos maiores do que a técnica utilizando parametrização local.

A técnica de parametrização local (Ajjarapu e Christy, 1992; Seydel 1994; Leite e Da Costa, 2003) consiste na troca de parâmetro próximo ao PMC. Leite e Da Costa (2003) propõem para o fluxo de potência continuado, o uso da formulação de injeção de corrente em coordenadas retangulares ao invés da formulação polar convencional (Ajjarapu e Christy, 1992). Os resultados mostram que esse método apresenta as mesmas características de convergência da formulação polar convencional, tanto no passo preditor, quanto no corretor, e que a sua utilização propicia uma grande redução no esforço computacional, devido principalmente ao fato dos elementos fora da diagonal da matriz Jacobiana serem elementos da matriz admitância nodal, os quais permanecem constantes durante todo processo iterativo. Ambos os métodos fazem uso do vetor tangente no qual a variável escolhida é aquela que apresentar a maior variação, sendo que o fator de carregamento (λ) passa a ser tratado como variável dependente, enquanto que a variável escolhida passa a ser o novo parâmetro. Essa técnica de parametrização tem demonstrado que ao aproximar-se do PMC, o parâmetro muda de λ para a tensão que apresenta a maior variação, retornando novamente para λ após alguns pontos. Embora em geral o uso desta técnica para a escolha automática do parâmetro não tenha apresentado dificuldades, o conjunto das barras cuja magnitude de tensão pode ser utilizada como parâmetro da continuação fica muito restrito, particularmente nos casos dos sistemas com problemas de estabilidade de tensão com características fortemente locais.

Souza (1996) apresenta a técnica de extrapolação quadrática que é baseada na análise do comportamento do vetor tangente como função do carregamento, e realiza a busca pelo PMC, através de tentativa e erro, a partir de dois pontos de operação conhecidos. De acordo com Souza et al. (1996) o método apresenta resultados tão precisos quanto os métodos da continuação, porém com um baixo esforço computacional. Silva (2007) utiliza a busca binária como estratégia de redução de passo durante o procedimento de obtenção do PMC. No método deseja-se encontrar um valor dentro de um determinado intervalo. No caso de não convergência, inicia-se uma busca binária a partir do último ponto convergido, reduzindo o passo preditor pela metade. O processo prossegue com uma seqüência de avanços, no caso de convergência, ou recuos em caso contrário, sempre reduzindo o intervalo de busca pela metade. A finalização do processo se dá quando o entre dois pontos consecutivos, um correspondente ao caso convergente e o outro ao não convergente, for inferior a um valor prefixado. O último ponto de operação convergido é considerado como o PMC. O método de busca binária alcança resultados próximos dos demais métodos, e o tempo computacional requerido e relativamente reduzido.

Em Bonini e Alves (2008) foi proposto um novo esquema de parametrização geométrica para o fluxo de carga continuado obtida a partir da observação das trajetórias de solução do fluxo de carga. Neste método a singularidade da matriz Jacobiana é eliminada pela adição da equação de uma reta que passa por um ponto no plano formado pelas variáveis fator de carregamento e a magnitude da tensão nodal. Observa-se que em vários métodos existentes na literatura pode ser necessário se efetuar a mudança de parâmetro ao longo do traçado da curva P-V, o que poderá acarretar mudanças na estrutura da matriz Jacobiana modificada. Nesse método, ao contrário do proposto por Ajjarapu e Christy, 1992, não se necessita realizar a troca de parâmetro ao longo de todo o traçado da curva P-V. Da mesma forma que no método proposto por Ajjarapu e Christy (1992), quando for necessário efetuar a mudança de coordenadas do centro do feixe de retas, esta também não implicará em mudanças na estrutura da nova matriz, mas apenas do valor do elemento correspondente a derivada da equação de reta em relação ao fator de carregamento, ou seja, no valor do coeficiente angular da reta (α). Uma vantagem importante acrescida pelo uso deste novo método foi à ampliação do conjunto de variáveis de tensão que podem ser adotadas como parâmetro da continuação. Este conjunto passou agora a incluir as barras cujas magnitudes de tensão permanecem constantes ao longo de uma faixa da curva P-V, e dessa forma não podem ser utilizadas como parâmetro para se obter essa parte da curva P-V. Também incluiu as barras cujas magnitudes da tensão apresentam uma inversão na sua tendência de variação simultaneamente com o fator de carregamento, i.e., os pontos de inflexão ("narizes das curvas") são coincidentes. Nesses casos, há coincidência da singularidade de ambas as matrizes Jacobianas no PMC (Bonini e Alves, 2008; Ajjarapu e Christy, 1992 ). Adicionalmente, conforme foram detalhadas em Bonini e Alves (2008), várias outras condições, p.ex. os limites de potência reativa e as contingências, podem modificar sensivelmente a curvatura da trajetória de soluções (curva P-V). Esse fato dificulta em muito a identificação de qual magnitude de tensão é mais apropriada para ser usada para a obtenção de todos os pontos da curva. Assim, apesar de todas as vantagens apresentadas pelo uso desta técnica, em alguns poucos casos, da mesma forma que no método proposto por Ajjarapu e Christy (1992), ainda permanece a necessidade do uso da técnica de parametrização local próximo do PMC, particularmente em sistemas com problemas de estabilidade de tensão com características fortemente locais.

Em Alves et alii (2000) foi proposto a adição da equação da perda total de potência ativa às equações do FC e assim, ao invés de especificar o carregamento e obter o estado convergido, especifica-se o valor da perda total de potência ativa, e obtém-se o estado convergido (ou o ponto de operação), incluindo o nível de carregamento para o qual o valor da perda total de potência ativa especificada ocorre. Adotando-se um passo fixo para o valor do novo parâmetro podem-se determinar, através de sucessivas soluções do novo sistema de equações, os demais pontos da curva P-V. A vantagem apresentada para o uso desta técnica era a de que na maioria dos casos analisados a parametrização local só se fazia necessária para pontos localizados pouco depois do PMC. Posteriormente verificou-se que para muitos sistemas, a singularidade da matriz Jacobiana aumentada praticamente coincidia com a da matriz Jacobiana do FC. Com isso, ainda permanecia a dificuldade em discernir se a divergência era conseqüente das limitações físicas do sistema ou de problemas numéricos. Com o intuito de se superar as limitações dessa técnica, no método apresentado em Garbelini et al. (2006) propôs-se acrescentar às equações do FC a equação da reta que passa por um ponto escolhido no plano formado pelas variáveis perda total de potência ativa e o fator de carregamento. No procedimento apresentado, inicia-se a obtenção dos pontos da curva P-V incrementando gradualmente o valor do coeficiente angular da reta que passa pela origem e pelo ponto correspondente ao caso base, obtido por um FC convencional. Na primeira divergência do método, retorna-se à solução anterior e efetua-se a redução do passo. Após uma nova divergência, adota-se a equação da reta pertencente ao feixe que passa pela última solução obtida e pelo ponto situado no eixo das abscissas (PMA), cujo valor da abscissa corresponde ao valor médio entre o fator de carregamento do caso base e o maior valor obtido antes do processo divergir novamente. Retoma-se o passo inicial e calculam-se apenas alguns pontos. Observa-se que o uso das retas pertencentes ao feixe que passa pelo ponto PMA é importante do ponto de vista da robustez do método, uma vez que o mesmo é necessário para eliminar a singularidade da matriz Jacobiana. O procedimento mostrou-se eficiente no traçado da curva P-V de todos os sistemas até então analisados.

Em Zhao e Zhang (2006) são apresentadas as limitações das técnicas de parametrização global para o caso de sistemas com estabilidade de tensão com característica predominantemente local, i.e., para sistemas cujo perfil de tensão de uma pequena área, ou magnitude de tensão de algumas poucas barras, não permanecem dentro da faixa normal de operação, conforme se pode ver nas figuras 1(a) e 1(b). A conseqüência direta dessa característica é que as curvas P-V da maioria das barras destes sistemas apresentam um "nariz agudo", como aquele apresentado na figura 2. As partes superior e inferior da curva apresentam praticamente a mesma inclinação, ao invés de um sinal oposto como no caso da barra crítica. Nesses casos ambos o fator de carregamento quanto a magnitude da tensão apresentam uma inversão simultânea na sua tendência de variação, i.e., os "narizes" são coincidentes. Essa coincidência, conforme já comentado, implica que a singularidade da matriz Jacobiana, a qual ocorre no PMC quando λ é usado como parâmetro, é coincidente com a singularidade da matriz Jacobiana modificada quando a magnitude de tensão de uma dessas barras é utilizada como parâmetro (Bonini e Alves, 2008). Como se pode observar nas figuras 1(c) e 1(d), onde se apresentam as curvas de perdas de potência ativa versus λ , correspondentes às das figuras 1(a) e 1(b) respectivamente, essa característica se reflete nos métodos que utilizarão técnicas de parametrização global, como os apresentados em (Chiang et al., 1995; Garbelini et al., 2006; Zhao e Zhang, 2006). Assim, no ponto de máximo carregamento, a singularidade da matriz Jacobiana do fluxo de carga também coincide com a da nova matriz. Nesses casos, conforme afirmado em Zhao e Zhang (2006), a parametrização local é considerada como a única forma de se eliminar a singularidade. Por outro lado, o método apresentado em Garbelini et al. (2006) não apresenta dificuldades numéricas, relacionadas com a singularidade da matriz Jacobiana, para sistemas com as características apresentadas na figura 1(d), entretanto, falha para os sistemas com as características apresentadas na figura 1(c).

Em Bonini e Alves (2008) mostrou-se que a mudança das coordenadas do centro do feixe de retas para o ponto médio (PM), situado entre os dois últimos pontos obtidos antes da primeira divergência e, portanto, próximo ao PMC, possibilitava a determinação do PMC utilizando a magnitude de tensão de barras cuja curva P-V apresentavam a característica de "nariz agudo", como as apresentadas na figura 2.

Neste trabalho são propostas técnicas mais eficientes de parametrização geométrica para o FCC. Nestas técnicas são utilizadas as equações de retas que passam através de pontos nos planos determinados pelas variáveis fator de carregamento e a somatória das magnitudes, ou dos ângulos, das tensões nodais de todas as barras do sistema. Apresenta-se o algoritmo para o traçado completo da curva P-V. O critério utilizado para se efetuar a mudança das coordenadas do centro do feixe de retas é baseado na análise da evolução do mismatch total de potência, a qual possibilita a identificação mais eficiente de condições de mal-condicionamento.

Também é proposta uma reformulação do método apresentado em Garbelini et al. (2007) de forma a não só simplificar o procedimento apresentado, reduzindo com isso o número de iterações necessárias para o traçado completo da curva P-V, mas também torná-lo eficiente para o traçado das curvas P-V de sistemas com estabilidade de tensão com característica predominantemente local. Estas vantagens são obtidas com duas alterações. A primeira consiste na mudança das coordenadas do centro do feixe de retas para o ponto médio (PM) ao invés do PMA. A segunda consiste no uso do critério baseado na análise da evolução do mismatch total de potência e não mais no número máximo de iterações prefixado. Uma outra vantagem é a simplificação do algoritmo proposto uma vez que agora o tamanho do passo pode ser mantido fixo ao longo de todo o procedimento, tornando-se desnecessário a sua redução na vizinhança do PMC.

Os métodos são aplicados para a obtenção da curva P-V do sistema teste IEEE 300-barras e de um sistema de 904 barras do sudoeste Americano, o qual apresenta problemas de estabilidade de tensão com características fortemente local. Os resultados mostram que o uso do PM próximo ao PMC é a principal estratégia do ponto de vista do aumento da robustez dos métodos.

2 FLUXO DE CARGA CONTINUADO PROPOSTO

O levantamento da curva P-V tem por objetivo determinar o quanto a demanda poderá aumentar antes que o sistema entre em colapso, ou seja, qual é a margem de carregamento (MC) para as condições operativas preestabelecidas. O traçado da curva P-V é realizado por meio de sucessivas soluções do sistema de equações do FC, G(θ,V)=0, considerando um crescimento da carga numa direção predefinida. De forma a se automatizar o processo de levantamento da curva P-V, o sistema de equações do FC é reescrito na seguinte forma

onde G é um vetor composto pelas equações dos balanços de potências ativa e reativa nodais, P^esp(λ) - P(θ,V)=0 e Q^esp(λ) - Q(θ,V)=0, respectivamente. Os vetores V e θ são os respectivos vetores das magnitudes e ângulos das tensões nodais. O vetor P^esp(λ) é o vetor da diferença entre as potências ativas gerada, P_G(λ) = , e consumida, P_C(λ) = , para as barras de carga (PQ) e de geração (PV). O vetor Q^esp(λ) é o vetor da diferença entre as potências reativas gerada, Q_G, e consumida, Q_C(λ) = , para as barras de carga PQ. Os vetores k_PG, k_PC e k_QC são parâmetros prefixados usados para caracterizar um cenário de carga específico. Eles descrevem as taxas de variação de potência ativa (P_G) nas barras de geração (barras PV), e das potências ativa (P_C) e reativa (Q_C) nas barras de carga (barras PQ). Assim, é possível realizar uma variação de carregamento individual, isto é, para cada barra do sistema, considerando para cada uma, um crescimento de carga com fatores de potência diferentes aos do caso base. Tradicionalmente, entretanto, assume-se que o aumento de carga de uma determinada área é feito com fator de potência constante e proporcional ao carregamento do caso base com modelo de carga de potência constante (nesse caso k_PG, k_PC e k_QC são todos iguais a um), visto que este fornece a condição operacional mais segura para o sistema (WSCC, 1998). Esta condição será adotada para todos os casos apresentados neste trabalho.

Nesse procedimento, P^esp e Q^esp são as variáveis independentes, enquanto que as magnitudes de tensão (V) e os ângulos de fase (θ), excetuando os da barra referência, são as variáveis dependentes. Com a inclusão de λ como variável na equação (1), o sistema resultará em n equações e n + 1 incógnitas. Assim, qualquer uma das n + 1 incógnitas pode ser definida como parâmetro. A diferença entre os métodos da continuação está na forma de tratar esta nova variável e em como contornar a singularidade da matriz Jacobiana (J). Dos quatro elementos básicos dos métodos da continuação, a parametrização, o passo preditor, o passo corretor e o controle de passo, a parametrização é o mais importante, uma vez que é através desta que se procura garantir a não singularidade da matriz J no PMC. Para a obtenção de curvas P-V, a adição de equações parametrizadas à equação (1) tem sido um procedimento padrão (Seydel, 1994). De uma forma geral a equação a ser acrescida pode ser colocada na seguinte forma:

onde α e β são coeficiente angulares que definem a reta a ser utilizada, e [y] é a medida escalar do vetor y=[y₁,...,y_n]^T, a qual pode ser escolhida dentre várias formas (Seydel 1994):

As duas primeiras formas encontram-se dentre as técnicas de parametrização local, enquanto que a última é um exemplo da técnica de parametrização global. As técnicas que utilizam o comprimento de arco (s) (Chiang et al., 1995) e a perda total de potência ativa (Alves et al., 2003) como parâmetro são outros exemplos de técnicas de parametrização global.

No método proposto por Ajjarapu e Christy, (1992), quando λ é usado como parâmetro α =1 e β =0, e quando [y]= ||y||_∞ = ||t||_∞, α =0 e β =1. O vetor t=[dθ, dV, dλ ]^T é o vetor tangente. No caso do uso de λ , o seu valor é incrementado gradualmente, a partir do caso base (λ⁰=1), até um valor para o qual não mais se obtenha solução, o processo iterativo do FC não converge ou diverge. Próximo ao PMC é efetuado a troca de parâmetro, e a variável com a maior variação no vetor t é a escolhida como o novo parâmetro, sendo que λ passa a ser, a partir daí, tratado como variável dependente. Esta característica de troca de variáveis corresponde à rotação de 90º do diagrama θ_k versus λ , ou V_k versus λ (Seydel, 1994). No caso, por exemplo, de V_k ser o novo parâmetro, seu último valor calculado é diminuído gradualmente e as respectivas soluções determinadas. O processo é continuado até que λ comece a decrescer. Nesse ponto λ pode ser novamente adotado como parâmetro, e seu valor diminuído gradualmente até que a curva esteja completamente traçada.

No presente trabalho considera-se [y]= Σy_k, sendo y o vetor das magnitudes (V) ou dos ângulos (θ) das tensões nodais. O parâmetro α é o coeficiente angular da reta (figura 3) que passa por um ponto escolhido "O" (λ⁰, ) no plano formado pelas variáveis λ e Σy_k. No presente caso escolheu-se β=1. Com a adição desta equação, λ passa a ser tratado como uma variável dependente e a como uma variável independente, ou seja, escolhida como parâmetro da continuação. Assim, o número de incógnitas é igual ao de equações, isto é, a condição necessária para que se tenha solução é atendida, desde que a matriz tenha posto máximo, isto é, seja não singular.

A partir da solução do sistema de equações (1) para o caso base (θ¹, V¹, λ¹=1), um passo preditor é executado para encontrar uma estimativa para a próxima solução. Os preditores mais utilizados são o tangente e o secante. No preditor tangente, a estimativa é encontrada dando um passo de tamanho apropriado na direção do vetor tangente à curva P-V, no ponto correspondente à solução atual (Ajjarapu e Christy, 1992). Os dois métodos secantes mais utilizados são (Seydel, 1994): o de primeira ordem, que usa as soluções atual e anterior, para estimar a próxima, e o de ordem zero, que usa a solução atual e um incremento fixo no parâmetro (θ_k, V_k ou λ , e α caso do método proposto) como uma estimativa para a próxima solução. O preditor trivial será a técnica adotada neste trabalho. Assim, o cálculo das demais soluções é realizado pelo fluxo de carga continuado proposto (FCCP), considerando sucessivos incrementos (Δα) no valor de α.

Finalmente, após se efetuar a previsão, torna-se necessário realizar a correção da solução aproximada para se obter a solução final. Na maioria das vezes o ponto obtido pelo passo preditor está próximo da solução correta e assim, poucas iterações são necessárias no passo corretor para a obtenção da solução correta, dentro da precisão desejada. O método de Newton é o mais usado no passo corretor. Neste passo a equação do tipo z0z^est=0, onde z e z^est correspondem à variável escolhida como parâmetro de continuação e seu respectivo valor estimado, obtido pelo passo preditor, pode ser acrescentada ao sistema de equações formado pelas equações ( 1 ) e (2). A solução também poderia ser obtida simplesmente fixando o valor do parâmetro em z^est, como é feito no caso do preditor de ordem zero. Assim, com a solução do caso base (θ¹, V¹ e λ¹ calcula-se o valor de α a partir do ponto inicial escolhido "O" (λ⁰, ) e dos seus respectivos valores obtidos no caso base "P" (λ¹, ):

Para α = α¹- Δα , a solução do sistema de equações formado pelas equações (1) e (2) fornecerá o novo ponto de operação (θ², V² e λ²) correspondente à interseção da trajetória de soluções (curva λ-[y]) com a reta cujo novo valor de coeficiente angular (α¹+ Δα) foi especificado. Para α = α¹, a solução convergida deverá resultar em λ =1. A expansão do sistema de equações formado pelas equações (1) e (2) em série de Taylor, incluindo somente os termos de primeira ordem, considerando o valor prefixado no valor do parâmetro α calculado para o caso base, resulta em:

onde x=[θ^T V^T]^T, J é a matriz Jacobiana do FC, e G_λ corresponde à derivada de G em relação a λ . ΔG e ΔR representam os fatores de correção (mismatches) das respectivas funções no sistema de equações. Deve-se observar que estes, serão iguais a zero (ou praticamente nulos, isto é, inferior a tolerância adotada) para o caso base convergido. Assim, somente ΔR será diferente de zero devido à variação de α , i.e., devido ao seu incremento Δα .

2.1 Procedimento Geral Adotado para Traçado da Curva P-V

Em função das análises realizadas definiu-se um procedimento geral para o traçado da curva λ-[y]. A curva P-V de qualquer barra do sistema é obtida armazenando-se, durante a obtenção da curva λ-[y], os correspondentes valores desejados. O procedimento utilizado para o traçado da curva λ-[y] é o seguinte:

As coordenadas iniciais do centro do feixe de retas, ponto "O", foram escolhidas levando-se em conta o fato de que as magnitudes de tensões e, portanto, a sua soma, diminuem com o aumento gradual de λ . O mesmo comportamento é apresentado pelos ângulos das tensões nodais, que tendem a se tornar mais negativos à medida que λ aumenta. As curvaturas das trajetórias de soluções de todos os sistemas analisados também apresentam esse mesmo comportamento, o que simplifica o estabelecimento do passos necessários para o sucesso do método. Assim, escolheu-se para a ordenada o mesmo valor do caso base, enquanto que para a abscissa, o valor zero e assim, os valores do coeficiente angular inicial será nulo (α¹=0). Observa-se que a escolha de uma outra coordenada, por exemplo (0,0), praticamente não afeta o desempenho do método durante o traçado completo da curva. Por outro lado, observa-se que o a mudança das coordenadas do centro do feixe de retas para o ponto médio (PM) é importante do ponto de vista da robustez do método, uma vez que é necessária para eliminar a singularidade da matriz Jacobiana. No caso da curva apresentada na figura 3(b), a singularidade da matriz J_m, definida na equação (4), ocorre no ponto em que a reta R tangencia a curva, ou seja, num ponto da curva situado entre os pontos "e" e "f". Voltando para a figura 3(a) verifica-se que no caso em questão, não há intersecção entre a reta tracejada (R) e a curva P-V e assim, o problema na realidade não apresenta solução. Voltando ao ponto "b" e efetuando a mudança do centro do feixe de retas, o problema volta a ter solução e com isso se obtém o ponto "c" e os demais, sem problemas com a singularidade da matriz J_m. Assim, vê-se que é desnecessário se efetuar sucessivas reduções de passo e dar continuidade ao traçado da curva a partir do ponto anterior (no caso, o ponto "b"). Observe também que a troca de feixe de retas implicará apenas na alteração do valor de α na matriz Jacobiana, ver equação (5), e não em mudanças na sua estrutura ou na criação de novos elementos. Da mesma forma que no caso do método apresentado em Garbelini et al. (2006), a curvatura dessas trajetórias é previamente conhecida assim, a mudança das coordenadas do centro do feixe de retas, quando se fizer necessária, é previamente estabelecida.

Antes de se efetuar a mudança das coordenadas do centro do feixe de retas de PM para as do ponto "O", logo após se verificar que o valor da ΣV_k (ou da Σθ_k) do ponto atual é maior que a do ponto anterior, recomenda-se o cálculo de mais alguns pontos. Essa recomendação se deve ao fato de se estar ainda muito próximo ao PMC e à característica muito acentuada apresentadas pelas curvas λ - ΣV_k (ou λ - Σθ_k) nessa região.

É importante ressaltar que o critério adotado pelo método para a mudança das coordenadas do centro do feixe de retas não se baseia apenas no número máximo de iterações adotado, que no caso é 10, mas também no critério baseado na análise do comportamento do mismatch total de potência apresentado em (Bonini e Alves, 2008). Esse mismatch é definido como sendo a soma dos valores absolutos dos desbalanços de potência ativa e reativa. Na maioria dos casos as parametrizações resultam numa convergência mais rápida, entretanto algumas vezes isso não ocorre. Nestes casos a evolução dos respectivos mismatches indica a possibilidade de mau-condicionamento. Com o uso desta análise obtêm-se uma redução do número de iterações, isto é, as mudanças sempre ocorreram antes de atingir-se o número máximo de iterações estipulado, mostrando-se assim, mais vantajoso para o traçado das curvas.

3 RESULTADOS

Para todos os testes realizados, a tolerância adotada para os mismatches foi de 10^-4 p.u. O primeiro ponto de cada curva é obtido com o método de FC convencional. Os limites de potência reativa (Q) nas barras PV's são os mesmos utilizados no método convencional de FC. Em cada iteração a geração de reativos de cada uma essas barras é comparado com seus respectivos limites. No caso de violação, ela é alterada para tipo PQ. Estas barras podem voltar a ser PV nas iterações futuras. As cargas são modeladas como de potência constante e o parâmetro λ é usado para simular incrementos de carga ativa e reativa, considerando fator de potência constante. Cada aumento de carga é seguido por um aumento de geração equivalente usando λ . O objetivo dos testes é mostrar a robustez do método proposto na obtenção da curva P-V de sistemas elétricos de potência.

3.1 Desempenho do Método Proposto

As figuras 4 e 5 apresentam os respectivos resultados da aplicação do método proposto ao sistema teste IEEE-300, i.e., considerando a parametrização pelo coeficiente angular (α) da reta situada no plano λ - ΣV_k, ou da reta situada no plano λ - Σθ_k. As figuras 4(a) e 5(a) mostram as respectivas curvas da somatória das magnitudes (ΣV_k) e dos ângulos (Σθ_k) das de tensões nodais, em função do fator de carregamento (λ). Observe a característica muito acentuada dessas curvas. Conforme afirmado em Zhao e Zhang (2006), nestes casos a parametrização local é considerada como a única forma de se eliminar a singularidade porque todas as demais técnicas de parametrização tais como a global que utiliza o comprimento de arco ou a apresentada em Iba et al., (1991), que usa um vetor perpendicular ao vetor tangente ao ponto da curva, falham na obtenção da solução no PMC. Com isso, esses métodos impossibilitariam se concluir que a divergência é devida à limitação física do sistema e não apenas restrições puramente numéricas. Os casos abaixo mostram que a nova técnica aqui proposta também possibilita o cálculo do PMC.

O processo se inicia com a obtenção, por meio de um FC convencional, das coordenadas do ponto "P", no caso (λ¹=0,8, =302,84 p.u.) para curva λ - ΣV_k e ( λ¹=0,8, = -32,01 rad.) para curva λ - Σθ_k. O processo se inicia com o valor de 0,8 para λ apenas para proporcionar uma melhor visualização da curva, uma vez que o ponto de máximo carregamento (1,0553) é muito próximo do valor do caso base ( λ =1,0), ou seja, a margem de carregamento desse sistema é de apenas 5,53%.

O valor escolhido para a ordenada do ponto "O" foi o mesmo do ponto P, ou seja, as coordenadas do ponto inicial escolhido são λ⁰=0,0 e =302,84 p.u. no plano λ - ΣV_k, e λ⁰=0,0 e = -32,01 rad. no plano λ - Σθ_k. Assim, os valores de ambos os coeficientes angulares iniciais serão nulos (α¹=0). O desempenho de qualquer método da continuação dependerá do tamanho do passo e do sistema em estudo, sua escolha não será, portanto, uma limitação exclusiva do método proposto. Maiores detalhes sobre as dificuldades relacionadas com o controle do passo podem ser encontrados em (Chiang et alii, 1995). Os passos (Δα) adotados para a obtenção dos demais pontos nas curvas λ - ΣV_k e λ - Σθ_k foram 1,0 e 6,0, respectivamente. Esses valores foram adotados visando o cálculo de poucos pontos ao longo da curva. Seu ajuste foi feito para o sistema IEEE-300. Partindo de um valor inicial pequeno 0,1 incrementou-se gradativamente o seu valor até um valor que se julgou como o mais apropriado foi Δα = 1 para coeficiente angular da reta localizada no plano λ - ΣV_k, e Δα = 6 para coeficiente angular da reta localizada no plano λ - Σθ_k. Já no caso do sistema de 904 barras adotou-se um valor proporcional a esse, i.e., Δα =NB/k sendo NB o número de barras e k a constante de proporção. Assim os seus respectivos valores foram 3 e 18. Possibilitando assim, um controle de passo proporcional para qualquer sistema a ser estudado.

As figuras 4(b) e 5(b) apresentam os detalhes da região do PMC, de onde se observam que tanto λ quanto ΣV_k e Σθ_k apresentam pontos de inflexão ("narizes") coincidentes, ou seja, a singularidade da matriz Jacobiana, a qual ocorre no PMC quando λ é usado como parâmetro, é coincidente com as singularidades das matrizes Jacobianas modificadas quando da parametrização global considerando ΣV_k ou Σθ_k. Dessa forma, essas variáveis não podem ser utilizadas como parâmetro para a obtenção do PMC porque o método apresentará dificuldades numéricas em sua vizinhança. Por outro lado, como se pode verificar das figuras, com o FCCP possibilita a determinação do PMC sem os problemas numéricos relacionados com a singularidade da matriz Jacobiana.

Os pontos da parte de cima das curvas apresentadas nas figuras 4 e 5 foram obtidos considerando os passos 1 e 2 do procedimento geral. Quando o FCCP não encontrar solução, tomam-se às coordenadas dos dois últimos pontos convergidos na curva λ - ΣV_k (ou λ - Σθ_k) e calculam-se as coordenadas do novo centro do feixe de retas, ponto médio (PM) nas figuras 4(b) e 5(b). Observe nestas figuras que o último ponto "b" do primeiro feixe de retas também é o primeiro ponto da primeira reta do segundo feixe de retas.

Utilizando os critérios do segundo feixe de retas e mantendo-se o mesmo passo inicialmente adotado para α , determinam-se os pontos pertencentes à região do PMC. Quando o valor da ΣV_k (ou da Σθ_k) do ponto atual for maior que o do ponto anterior, mudam-se novamente as coordenadas do centro do feixe de retas, de acordo com o passo 4 do procedimento geral. Assim, os demais pontos pertencentes à parte inferior das respectivas curvas λ - ΣV_k e λ - Σθ_k, foram obtidos considerando a equação da reta que passa pelas coordenadas do centro do feixe de retas inicial (ponto "O") e do último ponto calculado, ponto "e". O passo adotado é mantido. A troca do sinal se deve à mudança no sentido do traçado, que agora é anti-horário. Assim completa-se o traçado da parte de baixo da curva. As figuras 4(c) e 5(c) apresentam as curvas P-V da barra crítica (barra 526). Elas foram traçadas com os valores de λ e da magnitude da tensão da barra armazenados pelo FCCP durante o traçado das respectivas curvas λ - ΣV_k ou λ - Σθ_k. No PMC, os respectivos valores encontrados para λ e a magnitude de tensão (V₅₂₆) foram de 1,0553 e 0,735 p.u. Estes valores correspondem ao ponto "c" (1,0553, 294,07) da curva λ - ΣV_k apresentada na figura 4(b). No caso da parametrização por λ - Σθ_k, os valores encontrados para λ e para o ângulo da tensão da barra crítica (θ₅₂₆) foram 1,0553 e -77,652 rad., respectivamente. Observe nas figuras 4(d) e 5(d), que o método proposto possibilita o traçado de toda a curva P-V com um número reduzido de iterações.

Observe também nas figuras 4(c) e 5(c), que são computados um maior número de pontos na região do PMC. Isto ocorre naturalmente em virtude da proximidade das coordenadas do ponto PM com as do PMC. Assim, apesar do uso do mesmo tamanho de passo durante todo o traçado da curva, o método proporciona automaticamente um maior número de pontos na região do PMC, i.e., a mudança do centro de feixe de retas para o PM acarreta um controle automático de passo nessa região. Esta característica é vantajosa porque indiretamente possibilita o cálculo mais preciso do valor de λ do PMC sem a necessidade de uma estratégia de redução no passo de α.

As figuras 6 e 7 apresentam os respectivos resultados da aplicação do método proposto ao sistema teste IEEE-904, i.e., considerando a parametrização pelo coeficiente angular (α) da reta situada no plano λ - ΣV_k, ou da reta situada no plano λ - Σθ_k. Conforme se pode observar das figuras 6(a), 7(a) e 7(b), para este sistema as curvas λ - ΣV_k e λ - Σθ_k também são muito acentuadas. Com a utilização da técnica aqui proposta, é possível efetuar-se o traçado eficiente da curva e obter o PMC com uma boa precisão. Isso pode ser confirmado nas figuras 6(c) e 7(d), onde se verifica o baixo número de iterações necessários para que o FCCP, utilizando o procedimento geral apresentado no item 2.1, obtenha o PMC. O método apresentou um bom desempenho durante o traçado completo da curva. A principal razão para a melhoria no desempenho pode ser visto nas figuras 8(a) e 8(b), onde se vêem os comportamentos de ΣV_k e Σθ_k em função do novo parâmetro α. O ponto A corresponde ao PMC onde λ atinge seu valor máximo. No mesmo ponto ΣV_k e Σθ_k atingem seus valores mínimos e assim, no caso do uso de uma destas variáveis como parâmetro, a matriz Jacobiana modificada apresentará uma singularidade nesse mesmo ponto. Por outro lado, o PMC pode ser facilmente computado com o método proposto porque a singularidade da matriz Jacobiana é removida quando α é usado como parâmetro. Esta característica é a responsável pelo bom desempenho do FCCP na vizinhança do PMC.

A figura 9 mostra uma comparação entre o desempenho do FCCP e o método que utiliza a técnica de parametrização local apresentado em (Ajjarapu e Christy, 1992). A figura apresenta os resultados considerando um valor de 0,097 para o tamanho do passo (σ), e o critério para troca de parâmetro baseado na componente que apresentar a maior variação no vetor tangente. Observa-se que para valores de σ menores do que 0,097, o método apresentado em (Ajjarapu e Christy, 1992) não apresentou problemas numéricos para obter o PMC. Das figuras 9(a) e 9(b) verifica-se que, ao contrário do que ocorre com o FCCP, o método que utiliza a técnica de parametrização local falha após a previsão do ponto A, ou seja, antes da obtenção do PMC. Estes resultados confirmam que para ter sucesso na determinação do PMC, o método apresentado em (Ajjarapu e Christy, 1992) exige passos menores que o proposto.

A figura 10 apresenta o desempenho do método (FCCPa) apresentado em Garbelini et al. (2006), o qual utiliza a equação da reta que passa por um ponto escolhido no plano formado pelas variáveis perda total de potência ativa (Pa) e λ. Esse método permite a obtenção do PMC do sistema IEEE-300 sem os problemas numéricos relacionados com a singularidade da matriz Jacobiana. Entretanto, a aplicação do procedimento como foi proposto falha na determinação do PMC do sistema de 904 barras. A obtenção do PMC somente se torna possível com uma redução sensível do tamanho do passo (Δα), de 10 para 150. Entretanto, isso acarretará num cálculo de um elevado número de pontos na região do PMC e consequentemente, num grande número de iterações, conforme apresentado na figura 10(e). Por outro lado, é possível reformulá-lo (FCCPaR) de forma a não só simplificar o procedimento apresentado, reduzindo com isso o número de iterações necessárias para o traçado completo da curva P-V, mas também torná-lo eficiente para o traçado das curvas P-V de sistemas com estabilidade de tensão com característica predominantemente local. Estas vantagens são obtidas com duas alterações. A primeira consiste na mudança das coordenadas do centro do feixe de retas para o ponto médio (PM) ao invés do PMA, cujo valor da abscissa correspondia ao valor médio entre o fator de carregamento do caso base e o maior valor obtido antes do processo divergir novamente. Agora os valores das coordenadas do PM((ΣPaª+ ΣPa^b)/2, ( λª+ λ^b)/2) são calculadas a partir das coordenadas de dois pontos, "a" e "b" próximos ao PMC. A segunda consiste no uso do critério baseado na análise da evolução do mismatch total de potência apresentado em (Bonini e Alves, 2008), ao invés de adotar apenas um número máximo de iterações prefixado como limite. Uma outra vantagem é a simplificação do algoritmo proposto uma vez que agora o tamanho do passo pode ser mantido fixo ao longo de todo o procedimento, tornando-se desnecessário a sua redução na vizinhança do PMC. Todos os pontos da curva λ-Pa foram obtidos considerando-se um passo fixo (Δα) de 0,1; para passos maiores o método falha. Como se pode verificar agora, essa mudança de estratégia é suficiente para tornar esse método eficiente também na determinação do PMC desse sistema. Entretanto, observa-se que nesse caso, conforme se pode ver na figura 10(d), em relação ao método parametrizado pelo coeficiente angular (α) da reta situada ou no plano λ - ΣV_k, ou no plano λ - Σθ_k, ainda é necessário um maior número de iterações para se traçar o mesmo trecho da curva.

4 CONCLUSÕES

Neste trabalho apresentam-se propostas de parametrização que possibilitam o traçado completo das curvas P-V para qualquer sistema, inclusive para o caso de sistemas com estabilidade de tensão com característica predominantemente local. Para remoção da singularidade da matriz Jacobiana do fluxo de carga no PMC propõem-se a utilização da equação de uma reta que passa através de um ponto no plano determinado pelas variáveis fator de carregamento e a somatória das magnitudes, ou dos ângulos, das tensões nodais de todas as barras do sistema.

A reformulação proposta para o método apresentado em Garbelini et al. (2007) proporcionou um aumento da eficiência e uma simplificação do procedimento geral apresentado. O método tornou-se mais eficiente para o traçado completo da curva P-V de qualquer sistema, incluindo os sistemas com estabilidade de tensão com característica predominantemente local. Com as alterações propostas o tamanho do passo pode ser mantido fixo ao longo de todo o procedimento, tornando-se desnecessário a sua redução na vizinhança do PMC, reduzindo com isso, o número total de iterações.

Mostra-se também que a principal razão da eficiência dos métodos propostos para a determinação do PMC está na mudança das coordenadas do centro do feixe de retas para o ponto médio (PM), situado entre os dois últimos pontos obtidos antes da primeira divergência e, portanto, próximo ao PMC. Observa-se que em todas as parametrizações propostas, essa mudança do centro de feixe de retas para o ponto médio acarreta um controle automático de passo na região do PMC. Assim, os métodos propostos não só obtêm êxito em encontrar, com a precisão desejada, as soluções na região do PMC e no próprio ponto, mas também permitem a obtenção de soluções além deste (isto é, pontos da parte inferior da curva P-V) com um número baixo de iterações.

Observa-se que também que em vários métodos existentes na literatura pode ser necessário se efetuar a mudança de parâmetro ao longo do traçado da curva P-V, o que poderá acarretar mudanças na estrutura da matriz Jacobiana modificada. No caso dos métodos propostos não se necessita realizar a troca de parâmetro ao longo de todo o traçado da curva P-V, sendo que algumas vezes apenas se faz uma mudança de coordenadas do centro do feixe de retas. Essa mudança não implicará em alterações na estrutura da nova matriz, mas apenas no valor do elemento correspondente à derivada da equação da reta (R) em relação ao fator de carregamento (λ), ou seja, no valor do coeficiente angular da reta (α).

Os resultados apresentados no artigo contradizem a afirmação de que em sistemas com problemas de estabilidade de tensão com característica predominantemente local, a parametrização local é a única forma de se eliminar a singularidade da matriz Jacobiana. Os resultados obtidos para diversos sistemas confirmam, portanto, a maior eficiência dos métodos propostos e mostram sua viabilidade para aplicações no planejamento da operação nos atuais sistemas de gerenciamento de energia.

5 AGRADECIMENTOS

Os autores agradecem CNPq e a FAPESP pelo apoio financeiro.

Artigo submetido em 05/02/2009 (Id.: 00946)

Revisado em 01/04/2009, 15/10/2009, 15/03/2010

Aceito sob recomendação do Editor Associado Prof. Alexandre Bazanella

Datas de Publicação

Histórico