Resumen: El método de flujo de carga basado en incrustación holomorfa, o Holomorphic Embedding Load-Flow Method (HELM) por sus siglas en inglés, es un método para la solución del flujo de potencia usando técnicas de análisis complejo, resultando en un método determinístico, no iterativo, e independiente de estimaciones iniciales; a diferencia de los conocidos algoritmos iterativos, los cuales son susceptibles a dar soluciones incorrectas en ausencia de aproximaciones iniciales adecuadas. Modelos para Barra Slack Distribuida (BSD), los cuales se emplean para plantear el flujo de potencia llevando a cabo un reparto realista de pérdidas de potencia activa entre los generadores del sistema en lugar de asignar dichas pérdidas a un solo generador, se encuentran formulados únicamente para métodos iterativos. Por lo tanto, el presente artículo propone la novedosa formulación e implementación de un modelo de BSD para HELM. Este modelo novedoso es implementado usando el sistema IEEE de 118 barras, arrojando resultados de perfiles de tensión similares a los obtenidos usando un método convencional de Newton-Raphson.
Palabras clave:Barra slack distribuidaBarra slack distribuida,Holomorphic Embedding Load-Flow MethodHolomorphic Embedding Load-Flow Method,Flujo de potenciaFlujo de potencia,Método de Newton-RaphsonMétodo de Newton-Raphson.
Abstract: The Holomorphic Embedding Load flow Method (HELM) is an approach aimed at solving the power flow problem based on complex analysis techniques, resulting into a deterministic, non-iterative, and initial-estimation independent method, as opposed to well-known iterative algorithms, which are prone to incorrect solutions in absence of appropriate initial estimations. Distributed Slack Bus (DSB) models, which are a realistic mean of distributing active power losses among system generators instead of assigning the losses to a single generator, are formulated for iterative methods only. Thus, the present paper proposes a novel approach for a DSB model suitable for HELM. This novel approach is tested by using an IEEE 118-bus benchmark system, yielding very similar results to those obtained using a conventional Newton Raphson method.
Keywords: Distributed slack bus, Holomorphic Embedding Load-Flow Method, Power flow, Newton-Raphson method.
Artículos
Distribución de pérdidas de potencia activa en el flujo de potencia basado en el método de incrustación holomorfa
Active power loss distribution in power flow based on the holomorphic embedding method
Universidad de los Andes
Esta obra está bajo una Licencia Creative Commons Atribución-NoComercial-CompartirIgual 3.0 Internacional.
Recepción: 20 Octubre 2018
Aprobación: 12 Febrero 2019
Las redes eléctricas están en constante crecimiento y expansión, alcanzando gigantescas dimensiones con gran cantidad de generadores, y en algunos casos operando al límite de su capacidad cerca del borde de colapso. Es importante tener herramientas computacionales que logren determinar el estado estable de operación del sistema y que funcionen acorde a las realidades del mismo. El problema de flujo de potencia radica en la computación de la magnitud y ángulo de fase del voltaje en cada barra en un sistema de potencia bajo condiciones de estado estable y balance de las tres fases (Glover y col., 2012). Esta computación se lleva a cabo mediante la resolución del sistema de ecuaciones no lineales que se plantea con las conocidas Ecuaciones de Balance de Potencia (EBP). Los algoritmos clásicos que se utilizan para dar solución a este problema se fundamentan en métodos iterativos, siendo el de mayor uso el de Newton-Raphson (NR), por sus propiedades de convergencia cuadráticas. Dicho método tiende a exhibir en algunos casos problemas convergencia a una solución inviable y divergencia debido a su propia naturaleza no determinística, especialmente para sistemas altamente cargados o al no establecerse una estimación inicial apropiada. Diferentes trabajos se han desarrollado para obtener una estimación inicial que ayude a disminuir los casos de divergencia (Leonidopoulos 1995), (Klump y Overbye 2000), sin embargo, no está garantizado que este inconveniente se solucione, o más aún, que el algoritmo converja a la solución físicamente posible (Trias 2012).
A. Trias propuso el método de flujo de carga basado en incrustación holomorfa, o Holomorphic Embedding Load-Flow Method (HELM) (Trias 2012) por sus siglas en inglés, cuya característica principal es ser determinístico. Este método se basa en la incrustación de un parámetro complejo en las EBP para convertir los voltajes en funciones holomorfas de este parámetro, es decir, funciones analíticas en el plano complejo. Estas funciones holomorfas de voltajes se representan en series de potencia para así calcular su valor en el punto de interés mediante técnicas de continuación analítica. La propuesta original (Trias 2012) solo planteaba una incrustación para la barra slack y barras PQ; trabajos posteriores realizaron incrustaciones para la integración de barras PV (Rao y col., 2016). Por otra parte, una técnica desarrollada para resolver el problema de flujo de potencia aplicando límites de control, y a su vez, aumentar la precisión de la solución obtenida para este problema, es formulada en (Trias y Marín 2018). La principal ventaja de HELM es que, independientemente de valores iniciales, encontrará la solución correcta al problema de flujo de potencia en sistemas altamente cargados, o en su defecto, señalará la inexistencia de la misma.
El modelo clásico de flujo de potencia idealiza que las pérdidas del sistema son absorbidas por un solo generador, pero en la realidad los sistemas de potencia se componen de una gran cantidad de generadores entre los cuales las pérdidas son repartidas. Por este motivo, surge el modelo de Barra Slack Distribuida (BSD), el cual resulta ser más realista en comparación al modelo clásico. En (Meisel 1993), la metodología aplicada para implementar el modelo de BSD se fundamenta en la asignación de factores de participación a cada uno de los generadores para llevar a cabo dicha repartición. Diversos estudios han sido formulados para calcular estos factores de participación, basándose en distintos criterios (Meisel 1993, Zobian y col., 1997, Calovic y col., 1981). El modelo de BSD solo ha sido implementado en algoritmos iterativos para resolver el flujo de potencia, es por ello que surge la necesidad de aplicar este modelo mediante un método que no dependa de aproximaciones iniciales acertadas para dar resultados correctos; HELM cubre esa necesidad.
La contribución de este artículo es la novedosa formulación e implementación de un modelo de BSD para HELM, empleando la propuesta de selección de factores de participación desarrollada en (Meisel 1993), los cuales se calculan dependiendo de la generación de potencia activa en cada generador. Dicha formulación se lleva a cabo modificando las ecuaciones incrustadas para HELM formuladas en (Trias 2012), y (Rao y col., 2016), agregando una nueva variable Ploss que representa las pérdidas de potencia activa totales del sistema, esta se expresa, al igual que los voltajes, como función holomorfa del parámetro complejo de incrustación.
El resto del artículo se organiza de la siguiente manera: la sección 2 describe HELM y presenta todas sus ecuaciones; la sección 3 introduce la teoría de BSD a emplear; la sección 4 propone la formulación del modelo de BSD para HELM; la sección 5 explica los criterios tomados en cuenta para llevar a cabo la implementación de la formulación propuesta; la sección 6 presenta resultados numéricos de pruebas realizadas a la implementación de la formulación del modelo de BSD para HELM; finalmente, las conclusiones se presentan en la sección 7.
El método propone la incrustación de las ecuaciones originales de balance de potencia en una extensión de funciones holomorfas de sí mismas, para aprovechar las propiedades otorgadas por la analiticidad compleja (Trias 2012).
Como es bien sabido, dependiendo de las incógnitas y valores fijos para cada tipo de barra, se definen distintas EBP tradicionales. A continuación se presentan las distintas EBP tradicionales, seguidas de su incrustación holomorfa respecto al parámetro complejo s para un sistema de N barras, una barra slack, un conjunto u de barras PV y uno h de barras PQ:
• Barra Slack:
(1)Donde Visp es la magnitud de voltaje especificado para esa barra, en este caso la slack, y Visu voltaje. Su incrustación holomorfa se presenta como sigue (Rao y col., 2016):
(2)• Barra PQ:
(3)Yik representa la admitancia serie entre las barras i y k, Sies la potencia compleja inyectada en la barra i, y Vies el voltaje en la barra i. Luego de aplicar la incrustación de (3) respecto al parámetro complejo s se obtiene (Trias 2012):
(4)donde el término Yi shuntcorresponde al elemento que modela las admitancias en paralelo a la barra i; por otra parte, Yik transrepresenta el elemento (i, k) de la matriz de admitancias Y sin tomar en consideración admitancias en paralelo.
• Barra PV:
(5)
(6)Donde Pies la potencia activa inyectada en la barra i, la cual es una variable conocida para este tipo de barras, al igual que la magnitud de voltaje, esta última debe ser fijada a su valor especificado y esto se logra al incrustar (6). Esta incrustación se muestra a continuación (Rao y col., 2016):
(7)Asimismo, se fija el valor de Piincrustando a (3) de la siguiente forma:
(8)Se observa que en todas las EBP incrustadas, la variable de voltaje V pasa a convertirse en función holomorfa del parámetro complejo de incrustación s como V (s), al igual que la potencia reactiva inyectada en la barra i, Qi(s), en (8). Además, es notable que el término V * es incrustado como v *(s*), dado que de esta manera se obtiene holomorficidad (Trias 2015). Por otra parte, la inversa de las funciones de voltaje de todas las barras se indica como sigue (Rao y col., 2016):
(9)Las funciones holomorfas poseen la capacidad de ser representadas por medio de series de potencia, utilizando como parámetro de expansión al parámetro complejo s (Subramanian y col., 2013). Como ejemplo, se tiene la serie de potencia de voltaje de la barra 
, donde Vi[n] denota el coeficiente de orden n de la serie.
Para la aplicación de HELM, es necesario el cálculo de "germ solution", solución germinal, en español, que representa la solución del flujo de potencia bajo condición sin carga ni generación en el sistema; y se logra al evaluar s = 0 en las ecuaciones incrustadas, obteniéndose V (0) = V [0] = 1 y W (0) = W [0] = 1 para todas las barras, y Q(0) = Q[0] = 0 para las barras PV (Rao y col., 2016).
Desarrollando las series de potencia en ambos lados de cada una de las ecuaciones incrustadas e igualando los coeficientes del mismo orden de s, se plantean relaciones de recurrencia necesarias para la formulación de expresiones que permitan obtener el coeficiente de orden n de dichas series. En este sentido, para hallar los coeficientes de la serie de potencia que representa a (9), se tiene (Rao y col., 2016):
(10)Para definir algunas de las expresiones siguientes se utiliza la notación especial δnb, presente en términos que aparecen para coeficientes de una potencia específica de s en las series de potencia (Rao y col., 2016).
(11)Según el tipo de barra, las expresiones para obtener los coeficientes, para n > 0, son las siguientes (Rao y col., 2016):
• Barra Slack:
(12)• Barra PQ:
(13)• Barra PV:
(14)
(15)Es importante resaltar que para las barras PV, (14) fija la magnitud de sus voltajes a través de la parte real de sus coeficientes de series de potencia de voltaje. Por otra parte, los coeficientes de la serie de potencia de Qi(s) en (8), además de ser incógnitas, solo pueden ser de valor real. Lo mismo ocurre con la parte imaginaria de los coeficientes de las series de potencia de voltaje para las barras PV, y con la parte real e imaginaria de los coeficientes de las series de potencia de voltaje para las barras PQ, es decir, son incógnitas de valor real. Por lo tanto, para poder calcular estas variables desconocidas se hace necesario separar en parte real e imaginaria a los coeficientes de las series de potencia de voltaje, así como los demás elementos que conforman a (12), (13) y (15). De este modo, se plantea una ecuación matricial lineal de valores reales que se resuelve recursivamente con el fin de hallar la cantidad de coeficientes de series de potencia necesarios (Rao y col., 2016). Dicha ecuación está formada por (12), (13), (14) y (15).
A manera ilustrativa, en (19) se muestra la ecuación matricial resultante de un sistema de múltiples barras, donde .. corresponde al voltaje de una barra PQ y su ecuación ocupa las dos primeras filas del sistema, las dos siguientes filas corresponden a la ecuación de una barra PV y su voltaje correspondiente es V2, y donde los términos Gik y Bik corresponden a la parte real e imaginaria del elemento (i, k) de la matriz Ytrans, y se encuentran en la matriz denominada como matriz de recursión, la cual es cuadrada de tamaño 2N para un sistema de N barras.
Los valores de voltajes finales en el punto de interés se obtienen al evaluar s = 1 en las series de potencia de las funciones holomorfas correspondientes, sin embargo, el radio de convergencia de dichas series es normalmente mucho menor que 1, por lo que se hace necesario la aplicación de técnicas de continuación analítica con el fin de extender la región de convergencia de las mismas (Trias 2012). Por lo tanto, para garantizar la obtención de la máxima continuación analítica, (Trias 2012) propone emplear los aproximantes de Padé, una función racional que se denota como (Rao y col., 2016):
(16)donde L es el grado del polinomio del numerador, y M el grado del polinomio del denominador. Para aplicar este método, la serie de potencia debe tener un largo de L+M +1 términos; además, se obtiene una aproximación más precisa si se calcula el aproximante de Padé diagonal (L = M ) o cercano al diagonal (|L − M | = 1) (Baker y Graves-Morris 1996). En este artículo se emplea el método directo o de la matriz para encontrar los polinomios que conforman a (16), presente en (Baker y Graves-Morris 1996), del mismo modo, se evalúan solo los aproximantes de Padé diagonales.
La formulación de flujo de potencia por factores de participación, planteada en (Meisel 1993) propone que los generadores deben absorber una porción de las pérdidas de potencia activa totales del sistema. A esta porción se le denomina factor de participación, se asigna dependiendo de la capacidad de generación de potencia activa del generador, y se define como:
(17)donde e es el conjunto de generadores conformado por la barra slack y por las barras PV que tengan una potencia activa generada mayor a 0, y Pg icorresponde a la potencia activa generada en la barra i, la cual se calcula para la barra slack como sigue (Tong y Miu 2005):
(18)siendo Pd kla carga de potencia activa en la barra k.
Es importante resaltar que (18) introduce la variable conocida Pg slack al problema del flujo de potencia.
Un modelo de BSD para NR orientado a sistemas de distribución, es decir, con fases desbalanceadas, se encuentra definido en (Tong y col., 2005). En el presente artículo se implementa dicho modelo, sin tomar en consideración desbalances de fases, y utilizando los factores de participación definidos en (17).
(19)Para llevar a cabo dicha implementación, las pérdidas de potencia activa totales del sistema son convertidas en una nueva incógnita, denotada como Ploss (Tong y Miu 2005), por lo que se hace necesaria la inclusión de la ecuación de balance de potencia activa en la barra slack para poder resolver el problema de flujo de potencia. La ecuación de balance de potencia activa para las barras PV y la slack toman entonces la siguiente forma:
(20)donde Pise define como:
(21)Con el fin de representar el reparto de pérdidas de potencia activa basado en factores de participación, calculados a través de (17) y al definir a Pg slack, se modifica la ecuación incrustada que define a las barras PV:
(22)donde se introduce la incógnita Ploss convertida en función holomorfa del parámetro complejo s como Ploss(s), siguiendo el mismo tratamiento que se le da a las variables de voltaje y potencia reactiva para barras PV en la propuesta de HELM. Es importante notar que, al igual que Qi(s), Ploss(s) solo puede tener coeficientes de serie de potencia de valor real.
Para calcular la solución germinal, considerando la inclusión de la nueva función holomorfa, se resuelve el sistema de ecuaciones incrustadas que describe a los distintos tipos de barras, evaluando en s = 0. Dicho sistema de ecuaciones es:
(23a)
(23b)
(23c)
(23d)Obteniendo de (23a), (23c) y (23d) que para todas las barras V (0) = V [0] = 1, resulta:
(24)Es por ello que para satisfacer (24) se tiene que Ploss(0) = Ploss[0] = 0 y Qi(0) = Qi[0] = 0. Este es el resultado esperado de pérdidas de potencia activa considerando un sistema sin carga ni generación, condición planteada para determinar la solución germinal. Por otra parte, para todas las barras se tiene que W (0) = W [0] = 1.
Con el planteamiento de la nueva ecuación incrustada (22) para representar a las barras PV, se debe formular una nueva expresión con el fin de calcular el coeficiente de orden . de las series de potencia presentes en dicha ecuación. Por lo tanto, al desarrollar estas series en ambos lados de (22), e igualando los coeficientes del mismo orden de s, se obtiene que para n > 0:
(25)Debido a que ninguna otra ecuación incrustada se modifica al formular el modelo de BSD para HELM, las expresiones empleadas para obtener el coeficiente de orden n de las series de potencia en dichas ecuaciones no cambian.
El hecho de agregar la función holomorfa Ploss(s) para representar las pérdidas de potencia activa totales del sistema implica que una variable libre ha sido introducida en (22). Por consecuente, se hace necesaria la adición de una ecuación para poder resolver el flujo de potencia, tomando en cuenta que esta ecuación debe fijar el valor de Pg slack, ya que esta pasa a ser conocida, y que la función holomorfa Vslack(s) ya se encuentra definida por (2). Es por ello que se plantea utilizar (22) evaluada para la barra slack, por lo tanto, se genera la función holomorfa Qslack(s) como variable libre, además de Ploss(s). En ese sentido, se introduce (25) evaluada para dicha barra en el sistema de ecuaciones lineal para el cálculo de coeficientes de orden n de las series de potencia. Se puede entonces representar por medio de (26) a la ecuación matricial a ser resuelta recursivamente para hallar la cantidad deseada de coeficientes de series de potencia para un sistema de múltiples barras, donde la primera barra es PQ, la segunda es PV, y t es la barra slack.
De este modo, con respecto a (19), se observa que en (26) son agregadas dos filas en las que se realiza la inclusión de (25) evaluada para la barra slack; lo que produce que una columna, ubicada en la última posición, sea introducida a la matriz de recursión con el único propósito de añadir al elemento Qslack[n] presente en la parte imaginaria de (25). Asimismo, para agregar al elemento Ploss[n], se incluye una segunda columna adicional ubicada en la penúltima posición de dicha matriz, conformada por el negativo de cada uno de los factores de participación de los generadores del sistema, para representar al término Fi Ploss[n], el cual aparece en la parte real de (25). Por lo tanto, para un sistema conformado por . barras, la matriz de recursión es cuadrada de tamaño 2N + 2.
Adicionalmente, es importante destacar que la solución germinal mantiene sus mismos valores ya que la nueva ecuación agregada al sistema es igual a la de una barra PV.
(26)Para calcular el perfil de tensiones del sistema se lleva a cabo la continuación analítica de cada una de las series de potencia de los voltajes, la cual se puede realizar cada vez que se calculan dos coeficientes de las series. El criterio de convergencia manejado corresponde a la diferencia entre dos voltajes calculados sucesivamente mediante la continuación analítica; si el valor absoluto de la diferencia para cada barra es menor que una tolerancia especificada se ha obtenido convergencia. Para el manejo del cambio de tipo de barra PV a PQ se calcula la generación de potencia reactiva para cada barra PV después de haber obtenido un perfil de tensiones una vez comprobada la convergencia; luego se verifican los límites de potencia reactiva, y si existe violación de dichos límites se realiza el cambio de tipo de barra PV a PQ de manera tradicional para todos los casos en que esto ocurra, y se replantea el problema de flujo de potencia aplicando los cambios recién hechos.
En esta sección, empleando precisión aritmética de 64 bits para todos los cálculos realizados, se presentan resultados de pruebas realizadas a la implementación de la formulación propuesta, para el sistema IEEE de 118 barras, obtenido de (Christie 1993), el cual cuenta con 54 generadores, barras con susceptancias en paralelo y transformadores con posición de toma no nominal. Por otra parte, el perfil de voltajes obtenido de las pruebas aplicadas a dicha implementación es comparado, para la misma tolerancia, con el perfil de voltajes resultante de una implementación del modelo de BSD para NR como la descrita en la sección 3 de este artículo, siendo obtenidos resultados idénticos en ambos casos.
El primer experimento numérico consiste en obtener la cantidad de coeficientes de serie de potencia necesarios para alcanzar convergencia, a medida que se escala tanto el nivel de carga base de potencia activa y reactiva, como la generación de potencia activa de cada uno de los generadores del sistema de prueba, hasta llegar al punto de colapso. En la figura 1 se puede observar el resultado de la prueba recién explicada, para una tolerancia de 10−8, donde se hace notable el aumento de la cantidad de términos de serie de potencia en proporción al incremento de cargabilidad del sistema. El punto de colapso obtenido mediante la formulación propuesta se encuentra para una escala de cargabilidad de 1.92, siendo necesario el cálculo de 29 coeficientes de serie de potencia.
Fig. 1
Número de coeficientes de serie potencial vs. escala de cargabilidad del sistema.
Por otra parte, en la figura 2 se grafica el resultado de la segunda evaluación numérica, donde se busca determinar la cantidad de términos de serie de potencia en función de la tolerancia, considerando la carga base del sistema. Como se puede observar en la figura 2, se requiere calcular mayores cantidades de coeficientes de serie de potencia en la medida en que se incremente la precisión requerida.
Fig. 2
Número de coeficientes de serie potencial vs. tolerancia.
Es importante mencionar que al intentar llegar a una escala de cargabilidad mayor al punto de colapso determinado, la cantidad de coeficientes de serie de potencia requeridos aumenta de manera indefinida sin encontrar convergencia, esto puede deberse a que la solución al problema planteado para tal nivel de carga no es físicamente posible, o al hecho de que se excede la precisión utilizada para los términos de serie de potencia requeridos, por lo tanto se introduce error de redondeo en las operaciones numéricas siguientes, lo que imposibilita la obtención de convergencia para el problema. Este fenómeno provocado por el error de redondeo puede afectar a cualquier problema planteado que requiera el cálculo de una alta cantidad de coeficientes de serie de potencia, bien sea por la escala de cargabilidad establecida, por la tolerancia requerida, o aún más, por una combinación de ambos factores. Asimismo, es importante destacar que la implementación presentada no se encuentra optimizada.
Una formulación para la aplicación del modelo de BSD para HELM ha sido propuesta en este artículo, la cual fue implementada y probada en un sistema clasico de la IEEE. Es importante destacar que el perfil de voltajes obtenido de las pruebas numéricas aplicadas para dicha implementación fue comparado, para la misma tolerancia, con el perfil de voltajes resultante de una implementación del modelo de BSD para NR como la descrita en la sección 3 de este artículo, siendo obtenidos resultados idénticos en ambos casos; por lo tanto se verifica la validez de la formulación propuesta.
Pese a determinar un punto de colapso alto para el sistema de prueba, al intentar aumentar la escala de cargabilidad mas allá de dicho punto, o al exigir tolerancias muy bajas para condiciones de carga mayores a la carga base, se presentan problemas de convergencia producidos por errores de redondeo. Como trabajo a futuro se plantea dar solución a estos problemas, una de las opciones es acoplarse a trabajos desarrollados en este aspecto como los ya mencionados. También se recomienda llevar a cabo la optimización de la implementación de la formulación propuesta.
Fig. 1
Número de coeficientes de serie potencial vs. escala de cargabilidad del sistema.
Fig. 2
Número de coeficientes de serie potencial vs. tolerancia.
