Optimización de carga de vehículos eléctricos con asistencia solar fotovoltaica y almacenamiento de energía

Proyecto Fin de Carrera Ingeniería de Telecomunicación Formato de Publicación de la Escuela Técnica Superior de Ingeniería Autor: F. Javier Payán Somet Tutor: Juan José Murillo Fuentes Dep. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2013 Trabajo Fin de Grado Ingeniería Electrónica, Robótica y Mecatrónica Optimización de carga de vehículos eléctricos con asistencia solar fotovoltaica y almacenamiento de energía Autor: Antonia López Llamas Tutor: Francisco Rodríguez Rubio y Carlos Vivas Venegas Dpto. Ingeniería de Sistemas y Automática Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2022 Trabajo Fin de Grado Ingeniería Electrónica, Robótica y Mecatrónica Optimización de carga de vehículos eléctricos con asistencia solar fotovoltaica y almacenamiento de energía Autor: Antonia López Llamas Tutor: Francisco Rodríguez Rubio y Carlos Vivas Venegas Catedrático de Universidad y Profesor Titular Dpto. Ingeniería de Sistemas y Automática Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2022 Trabajo Fin de Grado: Optimización de carga de vehículos eléctricos con asistencia solar fotovoltaica y almacenamiento de energía Autor: Antonia López Llamas Tutor: Francisco Rodríguez Rubio y Carlos Vivas Venegas El tribunal nombrado para juzgar el trabajo arriba indicado, compuesto por los siguientes profesores: Presidente: Vocal/es: Secretario: acuerdan otorgarle la calificación de: El Secretario del Tribunal Fecha: Agradecimientos En primer lugar, quiero dar las gracias a mi familia por el apoyo y por el esfuerzo que han realizado durante estos años, por haber estado ahí en momentos de dificultad. Me gustaría agradecer también la ayuda y confianza que he recibido por parte de mis profesores de la Universidad, especialmente a Francisco Rodríguez Rubio, por darme la oportunidad de trabajar en este proyecto y a Carlos Vivas Venegas, por guiarme durante todo el proceso. Finalmente, quiero mencionar a mi profesora de tecnología del instituto IES Fray Luis de Granada, Purificación Entrena Núñez, ya que sin ella no hubiera llegado a donde estoy ahora. Me abrió los ojos al mundo de la ingeniería y creyó en mis capacidades para dedicarme a esto. Simplemente, gracias. En definitiva, me gustaría dar las gracias a todas las personas que me han acompañado durante esta etapa. Antonia López Llamas Sevilla, 2022 I Resumen Este trabajo aborda el problema del modelado y gestión de una estación de carga inteligente para vehículos eléctricos que sea adecuada para la carga rápida de CC al tiempo que proporciona un estrés mínimo en la red eléctrica. Se ha diseñado un prototipo de estación de carga con asistencia solar fotovoltaica y un sistema de almacenamiento para poder acumular el excedente y disponer de él en momentos de baja irradiancia. Teniendo como objetivo la minimización del coste monetario de la instalación, se ha empleado la técnica de optimización matemática Programación Cuadrática de Enteros Mixtos (MIQP) y un método basado en reglas heurísticas para gestionar de forma eficiente el flujo de potencia. Para mitigar la fuente de incertidumbre sobre la demanda de energía y la producción solar se han utilizado estrategias de predicción basadas en deep learning. III XÍndice 4.2 Estrategia de predicción 22 Modelo neuronal simple ANN 22 Modelo neuronal recurrente RNN 22 Modelo neuronal con memoria a largo plazo LSTM 22 Modelo neuronal con memoria a largo plazo bidireccional Bi-LSTM 22 Modelo neuronal convolucional CNN 22 4.3 Predicción de la potencia demandada por los vehículos 23 4.4 Predicción de la potencia suministrada por los paneles fotovoltaicos 23 5 Problema de optimización para determinar la acción de control 25 5.1 Técnica de optimización 25 5.2 Optimización de una estación de carga de vehículos eléctricos 25 5.3 Condiciones de factibilidad de la solución 27 5.3.1 Restricciones operativas y técnicas 27 Potencia de la red eléctrica 27 Potencia de la batería 28 5.3.2 Restricciones físicas 28 5.4 Planteamiento del problema de optimización 29 5.5 Adaptación del modelo 30 5.5.1 Transformación del producto cruzado de variables continuas y lógicas 30 5.5.2 Transformación de la expresión bidireccional entre variables continuas y lógicas 31 5.6 Programación 32 5.6.1 Programación de la función objetivo 32 Variable a optimizar x 32 Término cuadrático Q 32 Término proporcional c 32 Término independiente α33 5.6.2 Programación de las restricciones de la planta 33 Matriz de restricciones A 34 Vector de restricciones b 35 6 Algoritmo de gestión de la energía basado en reglas heurísticas 37 6.1 Funcionamiento en situación de sobrecarga 38 6.2 Funcionamiento en situación de subcarga 39 7 Resultados obtenidos y conclusiones 41 7.1 Gestión de la energía basado en técnica de optimización matemática 41 7.1.1 Escenario 1: Irradiancia solar alta 41 Potencia absorbida o cedida por el sistema de almacenamiento de energía 42 Potencia absorbida o cedida por la red eléctrica 43 Balance de potencia 44 Función de optimización 44 7.1.2 Escenario 2: Irradiancia solar baja 45 Potencia absorbida o cedida por el sistema de almacenamiento de energía 46 Potencia absorbida o cedida por la red eléctrica 47 Balance de potencia 47 Función de optimización 48 7.1.3 Comparativa con predicciones óptimas 48 7.2 Gestión de la energía basado en reglas heurísticas 49 Potencia absorbida o cedida por el sistema de almacenamiento de energía 49 Potencia absorbida o cedida por la red eléctrica 50 Balance de potencia 50 Función de optimización 51 7.3 Conclusiones 52 Índice XI 8 Lineas futuras 53 8.1 Modelado de los vehículos eléctricos 53 8.2 Definición del problema de optimización 54 8.3 Condiciones de factibilidad de la solución 55 8.4 Planteamiento del problema de optimización 55 Índice de Figuras 57 Índice de Tablas 59 Bibliografía 61 1 Introducción 1.1 Motivación En la actualidad, la sociedad es totalmente dependiente de la energía, ya sea en forma de calor, electricidad o para el transporte. A pesar de esto, el modelo de producción de energía actual es insostenible, ya que en nuestro país, más de la mitad de esa energía procede de fuentes no renovables [21]. Figura 1.1 Estructura de la generación española 2021. Dentro del 46.7% de la energía procedente de fuentes renovables[ 21 ], es decir, energía presente en la naturaleza de manera continua y prácticamente inagotable, sorprende que solo el 9.9% fuera obtenida de la energía solar, siendo España uno de los países más soleados de Europa. Por otro lado, la energía no renovable es un recurso limitado que además, contribuye al calentamiento global y a la contaminación de la atmósfera. Por ejemplo, cuando se quema cualquier combustible, se produce dióxido de nitrógeno, un contaminante del aire nocivo al que la Agencia Europea de Medio Ambiente (AEMA) en 2018 atribuyó 54.000 muertes prematuras [11]. El NO2 , uno de los gases responsables del agujero de la capa de ozono, es producido a partir del NO expulsado por los coches a motor. Por estas razones, son muchas las personas que han decidido adquirir un coche eléctrico de cero emisiones, una alternativa sostenible pero con otras desventajas respecto al vehículo de combustión. El principal inconveniente a tener en cuenta es el tiempo de carga de las grandes baterías de litio, que puede durar entre 6 y 8 horas. Muchos conductores se inclinan por instalar un cargador en su cochera y realizar la carga durante la noche. Sin embargo, este método solo es válido si una carga completa es suficiente para los desplazamientos realizados a lo largo del día, ya que la autonomía de estos vehículos es más limitada. 1 2Capítulo 1. Introducción 1.2 Descripción del problema Las electrolineras son puntos de carga de vehículos eléctricos situados en lugares públicos, áreas de servicio, entre otros. Surgen de la necesidad de recargar los vehículos cuando se realizan largos trayectos. Una estación convencional se alimenta de la red eléctrica, la cual no fue diseñada para abastecer tanta energía. Se estima que para 2035 habrá un total de 130 millones de vehículos eléctricos circulando por Europa [ 7 ], cifra que no dejará de crecer. Este aumento en la demanda puede provocar desviaciones en la tensión y reducción en la calidad del suministro eléctrico, llegando en casos extremos a un corte del mismo. Por esta razón aparece el concepto de "fotolinera". Se trata de una estación de carga de vehículos eléctricos alimentada por energía renovable, en este caso, con energía solar captada a partir de placas fotovoltaicas. Se trata de una alternativa sostenible y económica que suaviza el problema mencionado anteriormente. Figura 1.2 Estación de carga de vehículos eléctricos con soporte fotovoltaico[9]. Este tipo de estaciones de carga suelen incluir un sistema de almacenamiento de energía para acumularla en momentos donde haya excedente procedente del sistema fotovoltaico y utilizarla para cargar los vehículos en momentos de baja irradiancia. Además, usualmente van conectadas a la red eléctrica, pudiendo en algún momento concreto intercambiar energía en caso de necesidad por alguna de las partes. Para optimizar la carga de los vehículos, será necesario diseñar un controlador para gestionar la energía de la electrolinera. Esta no es una tarea sencilla, dado que la energía procedente del sol es intermitente y la demanda por parte de los vehículos eléctricos tiene una naturaleza cíclica con una componente estocástica importante. A fin de mitigar esta fuente de incertidumbre, se pueden emplear técnicas de control predictivo basadas en estimaciones de la demanda de energía y de la producción solar. El control predictivo es una estrategia de optimización que hace uso de modelos que describen el comportamiento del sistema para prever el efecto que tendrán sobre la misma las acciones de control futuras. Dichas técnicas pueden ser de naturaleza heurística o matemática y se describirán con más detenimiento en el apartado 2. Como se ha mencionado anteriormente, uno de los principales inconvenientes de los vehículos eléctricos es el tiempo de carga de sus baterías. Para aminorar este problema, se han diseñado cargadores de carga rápida que pueden suministrar una potencia igual o superior a 40 kW [ 28 ], lo que reduce el tiempo de espera considerablemente. Por ejemplo, los supercargadores de la marca Tesla suministran una potencia de 250kW, pudiendo cargar la batería en 15 minutos aproximadamente [24]. Por todo lo mencionado anteriormente, para alentar a la sociedad a adquirir vehículos eléctricos se deben instalar más puntos de carga que, además, permitan recargar las baterías de forma rápida. Se estima que para poder hacer la transición a la movilidad eléctrica, se requerirán más de 65 millones de cargadores para 2035 [ 7 ], siendo aproximadamente 220.000 los instalados en Europa el pasado mes de julio. España solo cuenta con 7.407 puntos, es decir, el 3,3% de la totalidad [1]. 1.3 Propósito del trabajo 3 1.3 Propósito del trabajo En este documento se pretende desarrollar una estrategia de gestión de la energía para maximizar el uso de los recursos renovable y reducir el coste monetario de una estación de carga de vehículos eléctricos. Para ello, se parte de un modelo formado por un sistema fotovoltaico (PV), un sistema de almacenamiento de energía (BES) y la red eléctrica (GRID), como se muestra en la figura. Figura 1.3 Esquema estación de carga de vehículos eléctricos. El modo de carga de los vehículos eléctricos (EV) puede realizarse mediante corriente alterna (AC), en cuyo caso el vehículo debe traer instalado un convertidor AC/DC, o corriente continua (DC), instalando el convertidor en la estación de carga. En este modelo se ha optado por la segunda opción, utilizando un bus DC para guiar la energía hacia los vehículos. La principal ventaja es el aumento en la velocidad de carga de los vehículos, ya que admite tecnología de carga rápida. Con la finalidad de diseñar un controlador robusto para gestionar la energía, se propone emplear técnicas de predicción de la demanda y de la producción solar mediante algoritmos de aprendizaje automático (machine learning). Partiendo de estos datos, se obtendrá una acción de control más precisa y eficaz cuanto mejor sean las estimaciones realizadas. 1.4 Objetivos Los objetivos que se pretenden alcanzar en este documento son los siguientes: •Documentación y estudio de modelos matemáticos de tecnologías de carga de vehículos eléctricos. •Modelado matemático de los elementos constitutivos de la instalación y su interrelación. • Ensayo y desarrollo de algoritmos inteligentes de gestión de flujos de potencia y validación para un conjunto de escenarios tipo. •Análisis comparativo de resultados. 4Capítulo 1. Introducción 1.5 Programas utilizados 1.5.1 MATLAB Figura 1.4 Logo Matlab. MATLAB (abreviatura de MATrix LABoratory) es una plataforma de programación y cálculo numérico utilizada por cantidad de ingenieros y científicos para analizar datos, desarrollar algoritmos y crear modelos. Se trata de una plataforma que ofrece un entorno de desarrollo con un lenguaje propio (M). Entre sus prestaciones se encuentran la implementación de algoritmos, la manipulación de matrices y el cómputo iterativo, herramientas necesarias para el desarrollo que se explicará a continuación. Además, se puede comunicar con otros programas como GUROBI aunque el lenguaje de programación sea diferente, lo que facilita la resolución del modelo. 1.5.2 GUROBI Figura 1.5 Logo Gurobi. GUROBI es una herramienta de optimización utilizada por más de 2500 empresas. Permite expresar problemas complejos mediante modelos matemáticos y los resuelve de forma iterativa, encontrando la solución óptima entre millones de posibilidades. En este caso se va a utilizar para resolver el problema de optimización una programación cuadrática de enteros mixtos (MIQP) con restricciones. 2 Estado del arte El término optimización hace referencia a la asignación de los recursos disponibles de forma óptima con el fin de maximizar o minimizar la función objetivo, cumpliendo las restricciones impuestas. La integración de energía renovable y la demanda irregular por parte de los vehículos eléctricos plantea varios retos a la hora de determinar la acción de control a aplicar, debido a la incertidumbre sobre la oferta y demanda de energía, pudiendo poner el peligro la estabilidad y la fiabilidad de la red [27]. Por ejemplo, según [ 8 ], si el 60% del total de los vehículos eléctricos se carga durante las horas de menor consumo, las pérdidas de energía pueden llegar al 40% en función de la estrategia de carga utilizada. Por lo tanto, se utilizan técnicas de optimización para superar estos problemas y conseguir el máximo beneficio en la estación de carga. En este apartado se van a mostrar técnicas de optimización heurística y técnicas de optimización matemática, dividiendo este último en técnicas lineales y no lineales, incluyendo varios ejemplos de cada una de las opciones. 2.1 Diferentes métodos de optimización 2.1.1 Técnicas de optimización heurísticas Un método heurístico es un procedimiento en el cual se utilizan una serie de aproximaciones intuitivas para resolver un problema complejo. Hay problemas de optimización donde no se puede encontrar soluciones óptimas usando métodos numéricos exactos con esfuerzos computacionales aceptables, ya que el esfuerzo puede crecer de forma exponencial respecto al número de variables de decisión (fenómeno explosión combinatorial [ 12 ]). En estos casos puede ser conveniente recurrir a un modelo heurístico. En este tipo de técnicas no se puede garantizar la obtención del resultado óptimo, pero sí de una solución aceptable. Además, son métodos no deterministas, lo que indica que los resultados pueden variar aunque la condiciones sean las mismas. La mayoría se basan en la observación y análisis de la naturaleza, como se comportan los animales, insectos, incluso la evolución de las especies. Se ha realizado un breve resumen de las más conocidas y utilizadas, sobre todo en el campo de los vehículos eléctricos. 5 6Capítulo 2. Estado del arte Genetic Algorithm El algoritmo genético [ 10 ] es un método adaptativo con un comportamiento intrínseco no determinista que utiliza analogías inspiradas en la biología para obtener soluciones cuasi-óptimas, garantizando mayor eficiencia que si se realiza una búsqueda completamente aleatoria. Se emplean técnicas metaheurísticas para guiar el proceso de búsqueda de soluciones, alcanzando, en general, soluciones aceptables próximas al óptimo buscado. Se parte de una población de individuos, la cual evoluciona siguiendo dos mecanismos principalmente; la mutación genética y la selección natural, en los cuales se descarta a los que tengan peores aptitudes para el problema concreto. Este algoritmo no sigue ninguna teoría matemática estricta y es el investigador el que debe determinar el mecanismo de evolución del desplazamiento en función del problema. Cada individuo tiene su información codificada en una cadena, normalmente binaria, llamada cromosoma. Los cromosomas evolucionan en cada generación y se evalúa su aptitud siguiendo algún criterio concreto. La implementación clásica de este algoritmo sería: •Iniciación: Se genera una población inicial de individuos. •Clasificación: A cada cromosoma se le evalúa según una función de aptitud y se clasifica. • Selección: Se decide qué individuos se reproducen, siendo los más aptos los que tienen más probabilidad de ser seleccionados. • Recombinación: Se generan dos descendientes procedentes de los individuos seleccionados, combinando sus características. • Mutación: Algunos de los cromosomas pueden verse alterados de forma aleatoria, permitiendo alcanzar nuevos resultados que no estaban cubiertos anteriormente. •Reemplazo: Finalmente, la generación siguiente estará formada por los individuos más aptos. El proceso se repite hasta alcanzar el número máximo de generaciones o cuando no se produzcan cambios en la población, ya que se desconoce la solución óptima al problema. Particle Swarm Optimization El método de optimización por enjambre de partículas fue desarrollado por James Kennedy y Russell Eberhart en 1995 [ 15 ] y se inspira en la evolución del comportamiento colectivo, tratando de imitar las conductas sociales de varios grupos de animales. Se puede ver un ejemplo de esta técnica en la migración de las aves en busca de climas más cálidos. Estos animales parten de una posición inicial aleatoria y durante el vuelo, la bandada va comprobando las condiciones meteorológicas del lugar. En función de la medida del conjunto, se orienta en una dirección u otra. El proceso se repite con el objetivo no garantizado de hallar en algún momento una solución lo suficientemente satisfactoria. Se trata de un método adaptativo que utiliza los siguientes principios [3]: • Evaluación: El aprendizaje de los organismos vivos se basa en la observación y evaluación del lugar, con el fin de distinguir si las características son atrayentes. • Comparación: Se comparan las características evaluadas con las obtenidas por el resto de partículas. • Imitación: Lorenz asegura que solo los seres humanos y algunas aves son capaces de imitar. La imitación es central para la adquisición y mantenimiento de las habilidades mentales. Este método permite resolver problemas con funciones no lineales. Comprende un concepto sencillo y puede ser implementado fácilmente utilizando operadores matemáticos básicos. Además, es computacionalmente barato en términos de requisitos de memoria y velocidad. 2.1 Diferentes métodos de optimización 7 Simulated Annealing El algoritmo de recocido simulado se basa en una analogía con el proceso de recocido del acero, en el que el material se calienta y luego se enfría lentamente para variar sus propiedades físicas. Se debe controlar la velocidad de enfriamiento para conseguir la resistencia o flexibilidad adecuadas, teniendo como objetivo alcanzar el estado de mínima energía del metal. Se utiliza para resolver problemas de optimización combinatoria de gran escala. Para definir el algoritmo se deben especificar siete componentes: •Temperatura inicial. •Criterio de parada: Determina cuando se termina la ejecución del algoritmo. •Criterio de explotación: Regla para seleccionar una solución del entorno de vecindad. •Criterio de aceptación: Indica si la nueva solución sustituye a la anterior. •Longitud de la temperatura: Determina si la temperatura debe ser actualizada. • Esquema de enfriamiento: Actualiza la temperatura a un valor inferior si la longitud de la temperatura lo indica. •Reinicio de la temperatura: Restablece la temperatura a un valor determinado. Es una técnica de búsqueda local estocástica que, partiendo de una solución inicial, se va explorando su vecindad de forma iterativa, seleccionando una solución candidata según el criterio de explotación. La nueva solución se compara con la actual utilizando el criterio de aceptación, dando lugar a dos situaciones: • Si la nueva solución tiene una energía o una evaluación menor respecto al anterior, se sustituye por la que se tenía anteriormente. • En caso contrario, se le asigna una probabilidad para ser seleccionado, que depende de las diferencias en las evaluaciones y la temperatura del sistema. Conforme decrece la temperatura, la solución con peor evaluación es menos probable de ser seleccionada. Está demostrado que si la temperatura se va disminuyendo con la suficiente lentitud, el proceso converge a la solución óptima [4]. Cuando se alcance el criterio de parada, la solución final será la seleccionada. Figura 2.1 Gráfico ejemplo simulated annealing [17]. 14 Capítulo 3. Modelo matemático 3.2 Paneles solares 3.2.1 Potencia teórica Los paneles fotovoltaicos son un conjunto de celdas fotovoltaicas que producen energía eléctrica a partir de la energía solar incidente mediante el efecto fotoeléctrico. La cantidad de energía que absorben las celdas depende del ángulo de incidencia de los rayos solares, ya que solo la componente perpendicular a la superficie de los paneles será transformada en energía eléctrica. El resto, en cambio, será reflejada. Figura 3.1 Rayo incidente sobre las placas fotovoltaicas [13]. Suponiendo que los paneles solares tienen un ángulo de inclinación fijo, se obtiene la máxima potencia durante el año si se instalan con un ángulo de inclinación θPV igual a la latitud del lugar [13]. En las tablas se muestran todos los datos utilizados para los cálculos de la potencia obtenida del sistema fotovoltaico, empleando las coordenadas geográficas de la Escuela Técnica Superior de Ingeniería de Sevilla. Tabla 3.2 Datos de inicialización del sistema fotovoltaico. Nombre Unidades Valor Definición Sur f acePV m250 Superficie de paneles λPV grados 37.41136 Latitud de la ubicación de la planta φPV grados -6.00055 Longitud de la ubicación de la planta Tabla 3.3 Datos de inicialización de las celdas de las placas fotovoltaicas. Nombre Unidades Valor Definición θPV grados 37.41136 Inclinación de las placas con respecto al suelo GPV re f W/m21000 Radiación solar de referencia PwPV nom W 800e3 Potencia nominal de la planta solar instalada PwPV stc W 1000 Potencia nominal del panel en STC a ηPV % 0.85 Eficiencia del panel TcPV stc ºC 25 Referencia de temperatura de la celda en STC βPV %/ºC 0.45 Coeficiente de temperatura per ºC de celda Tamb ºC 25 Temperatura ambiente NOCTPV ºC 55 Normal Operation Cell Temperature aSTC: standard test conditions 3.2 Paneles solares 15 Inicialmente, haciendo uso de la Calculadora de Posición Solar de la NOAA [ 22 ], se calcula la declinación solar en grados para la fecha actual. Los cálculos se basan en ecuaciones de Algoritmos Astronómicos de Jean Meeus. Los resultados son teóricamente precisos con una exactitud de un minuto para ubicaciones entre ±72◦ de latitud, y con una exactitud de 10 minutos fuera de esas latitudes. Sin embargo, los valores pueden variar respecto a los reales debido a cambios de presión, temperatura o condiciones atmosféricas. •Ángulo de declinación solar: Es el ángulo que forma el ecuador y una línea trazada del centro de la Tierra hasta el centro del sol. En los equinoccios de primavera y otoño su valor es 0◦ y debido a la inclinación de la Tierra sobre su eje de rotación y a la traslación de la Tierra alrededor del sol, su valor varía hasta ±23,45◦. [13] Figura 3.2 Ángulo de declinación solar [13]. Tras conocer el ángulo de declinación y utilizando las coordenadas geográficas, se obtienen: • Ángulo cenital solar: Posición del sol relativa a la normal, es decir, es el ángulo que forma el sol respecto a la vertical. Por tanto, será el ángulo complementario a la elevación (α)[13]. • Ángulo de Azimuth: Es el ángulo desde donde proviene la luz del sol medido en sentido horario respecto al norte geográfico [ 13 ]. En este documento se ha supuesto que si las coordenadas geográficas se sitúan en el hemisferio norte, el ángulo de Azimuth se medirá en sentido antihorario respecto al sur geográfico. (a) Ángulo cenital. (b) Ángulo Azimuth. Figura 3.3 Ángulo cenital y Azimuth [13]. Estos cálculos se realizan dentro de las funciones EarthEphemeris.m, sunsag.m y sunslope.m, extraídas del repositorio de GitHub SunPosition. [5] 16 Capítulo 3. Modelo matemático Previo al cálculo de la potencia obtenida por los paneles fotovoltaicos en condiciones de operación reales, se deben consultar los siguientes parámetros en la ficha de características técnicas del fabricante [14]. • Normal Operation Cell Temperature ( NOCTPV ): Es la temperatura que alcanzan las celdas solares en función de la temperatura exterior y la irradiación solar. • Coeficiente de temperatura per ◦C de celda ( βPV ): Es un coeficiente de pérdida de potencia por cada grado superior a los 25◦C. Normalmente su valor está comprendido entre 0,43 y0,47%/◦C. Además, debe tenerse en cuenta: • Irradiancia solar: Se define como el caudal de energía solar que llega al panel tras atravesar la atmósfera terrestre, la cual actúa como filtro para el espectro lumínico. Se supone una irradiancia solar de referencia de 1000 W/m2 , que es el valor utilizado en las pruebas estándar. Sin embargo, como se puede observar en el siguiente gráfico, su valor depende de la estación y la hora del día en la que se encuentre. Figura 3.4 Gráfica de la irradiación solar durante el día [25]. A partir de estos datos, se calcula la radiación solar perpendicular al panel fotovoltaico ( GPV t ) ya que solo esta componente será la que se transforme en energía eléctrica. GPV t=GPV re f ∗(cos(α)∗sin(θPV )∗cos(azimuth)+sin(α)∗cos(θPV )) (3.1) La temperatura a la que están sometidas las celdas es un dato clave a la hora de conocer la potencia obtenida de los paneles fotovoltaicos, ya que estos trabajan de forma más eficiente, es decir, convierten una mayor proporción de radiciación solar en energía eléctrica, si se encuentran a menor temperatura. Precisamente, el término βPV indica las pérdidas ocasionadas debido al aumento de temperatura de la celda. TPV c=Tamb + (NOCTPV −20)∗GPV t 800 (3.2) A partir de este punto, se puede obtener la potencia extraída por un panel fotovoltaico con la siguiente expresión. PPV c=GPV t GPV re f ∗PwPV stc ∗ηPV ∗(1−βPV ∗(T cPV −TcPV stc )) (3.3) Finalmente, se multiplica el resultado por la superficie total de paneles y se comprueba que el resultado sea un valor positivo, ya que en términos matemáticos podría ocurrir que la potencia obtenida tuviera un valor negativo. En cambio, en la realidad no puede ocurrir que las placas fotovoltaicas absorban potencia del sistema. PPV =Sur f acePV ∗max([0PPV c]) (3.4) 3.2 Paneles solares 17 3.2.2 Potencia real En el modelo teórico anterior se ha supuesto que los paneles fotovoltaicos reciben energía solar durante todas las horas de sol del día. En cambio, si ese día llueve o simplemente se sitúa una nube sobre los paneles, la energía solar recibida será mucho menor o incluso nula. Para que el modelo fuera más realista, los datos utilizados son valores históricos de una planta fotovoltaico experimental situada en el tejado de los laboratorios de electrónica de la Escuela Técnica Superior de Ingeniería de Sevilla. Se han ordenado en forma de tabla para tener un registro por horas de la potencia generada por el sistema fotovoltaico. La estructura de la base de datos es la siguiente: Tabla 3.4 Base de datos de la demanda de potencia de vehículos eléctricos. Index Irradiancia (W/m2) ... ... 2019-03-31 07:24:39+00:00 0.62 2019-03-31 08:24:39+00:00 30.77 2019-03-31 09:24:39+00:00 268.78 2019-03-31 10:24:39+00:00 474.04 2019-03-31 11:24:39+00:00 480.30 2019-03-31 12:24:39+00:00 418.59 ... ... En la siguiente figura se muestra los diferentes valores de potencia generada por el sistema fotovoltaico teóricamente y de forma experimental. Feb 01 Feb 02 Feb 03 Feb 04 tiempo (seg) 2019 0 5 10 15 20 25 Potencia (KW) Potencia PV Experimental Teórica Figura 3.5 Comparativa potencia experimental vs potencia teórica. Como se puede apreciar, el día 1 de febrero la curva teórica se asemeja a la experimental, lo que indica que fue un día bastante soleado con algunos intervalos nubosos. En cambio, el día 3 de febrero fue un día encapotado, provocando que la generación de energía sea prácticamente nula durante todo el día. 18 Capítulo 3. Modelo matemático 3.3 Sistema de almacenamiento de energía Se trata de un conjunto de baterías de litio cuya estructura consiste en un electrodo positivo (cátodo de litio) y un electrodo negativo (formado por un ánodo de carbono) [ 18 ]. Mediante una serie de reacciones químicas se consigue la acumulación y posterior devolución de la energía, en este caso generada por el sistema fotovoltaico. Inicialmente, se deben consultar los siguiente parámetros en la ficha de características técnicas del fabricante: Tabla 3.5 Datos de inicialización del sistema de almacenamiento de energía. Nombre Unidades Valor Definición PBES min W -15000 Potencia mínima de/hacia las baterías PBES max W 15000 Potencia máxima de/hacia las baterías CapBES max W 17600 Capacidad de almacenamiento máxima SOCBES min % 0.05 Porcentaje mínimo de carga SOCBES max % 0.95 Porcentaje máximo de carga µBES ch % 0.8 Rendimiento de carga µBES disc % 0.8 Rendimiento de descarga ηBES ch % 0.90 Eficiencia de carga ηBES disc % 0.95 Eficiencia de descarga CCostBES €/Wh 0.125 Coste capital de las baterías CostBES degr,ch €/W2h 1e-9 Coste degradación de carga CostBES degr,disc €/W2h 1e-9 Coste degradación de descarga CyclesBES ciclos 3000 Ciclos de vida de las baterías La potencia suministrada o absorbida por las baterías situadas en la instalación se decide según la acción de control calculada mediante un algoritmo de optimización que se explicará en el apartado 5. Se ha decidido adoptar un modelo simple en el que se tiene en cuenta una eficiencia asimétrica para carga y descarga. Si la acción de control obtenida en un instante determinado ( PBES control(tk) ) es negativa, indica que se ha decidido cargar las baterías. PBES ch =−PBES control(tk)∗µBES ch ∗ηBES ch (3.5) En caso contrario, se procede a la descarga de las mismas. PBES disc =PBES control(tk)∗µBES disc ηBES disc (3.6) Además, se debe comprobar que no se supera el porcentaje máximo o mínimo de carga. Esto se debe a que si se supera alguno de estos límites se está forzando el sistema de almacenamiento, degradando las celdas de iones y por tanto, reduciendo el rendimiento y la energía máxima que son capaces de almacenar. Hay que tener presente que el número de ciclos de carga es limitado y viene indicado en ficha de características técnicas del fabricante ( CyclesBES ). Entonces, si se evita realizar cargas o descargas completas, se prolongará su tiempo de vida útil. 3.4 Vehículos eléctricos 19 Finalmente, se comprueba que no se supera la potencia mínima o máxima a la que pueden cargarse/descargarse. Si ocurriera alguno de estos casos se limitaría la potencia a ese valor. Potencia sumistrada/absorbida por las baterías: PBES(k) = −PBES ch +PBES disc (3.7) Porcentaje de carga actual de las baterías: SOCBES(k) = SOCBES(k−1)∗CapBES max +PBES ch −PBES disc CapBES max ;(3.8) 3.4 Vehículos eléctricos Se ha consultado la base de datos ACN-Data [ 30 ] para simular la llegada de vehículos eléctricos a la planta. Se trata de una base de datos pública que incluye registros reales de varias electrolineras de EEUU. Ha sido adaptada para tener un registro por horas de la potencia demanda por los vehículos. La estructura de la base de datos es la siguiente: Tabla 3.6 Base de datos de la demanda de potencia de vehículos eléctricos. Index val_sessions val_power_median (kW) 2019-01-01 00:24:39+00:00 3.0 2.854 2019-01-01 01:24:39+00:00 3.0 1.065 2019-12-10 02:24:39+00:00 2.0 4.782 2019-12-10 03:24:39+00:00 1.0 0.088 ... ... ... En cada instante, se consulta la potencia media (val_power_median) de la fecha y hora que coincida con la fecha y hora actual, sin importar los minutos y los segundos especificados en la tabla. Este dato se utiliza como potencia demandada por todos los vehículos durante esa hora para realizar el resto de cálculos. 3.5 Red eléctrica La potencia absorbida o cedida a la red será la que se obtenga haciendo el balance de potencia, ya que se debe cumplir que la potencia consumida por los vehículos debe ser igual a la suma de la potencia suministrada por las diferentes fuentes. PGRID(k) = PEV (k)−PBES(k)−PPV (k)(3.9) Se han definido unos límites de potencia que se desea consumir o ceder de la red utilizados para el cálculo de la acción de control, pero estos límites son orientativos y se podrían modificar si se necesitara. Este caso podría ocurrir si, por ejemplo, se predice que la demanda de vehículos es mucho menor a la que se tiene en realidad y no se puede abastecer con la potencia obtenida del resto de fuentes. Si ocurriera, se debe aumentar este límite para poder determinar una acción de control que cumpla las restricciones de diseño. Tabla 3.7 Datos de inicialización de la red eléctrica. Nombre Unidades Valor Definición PGRID max W 30e3 Máxima potencia que se desea consumir PGRID min W -30e3 Máxima potencia que se desea ceder 4 Estrategias de predicción de la demanda y producción solar La energía que generará el sistema fotovoltaico y la demanda de energía procedente de los vehículos eléctricos son datos que no se pueden conocer a priori y serán elementos decisivos para seleccionar la mejor acción de control a aplicar sobre el sistema. Por consiguiente, se requiere efectuar una predicción sobre lo que es más probable que ocurra, la cual se basa en datos pasados y en otros factores como puede ser la temperatura o humedad en el ambiente. Se ha decidido utilizar técnicas de Deep Learning para aprender el comportamiento de la demanda y de la generación fotovoltaica. El desarrollo completo de este apartado ha sido realizado en el Trabajo Fin de Grado Predicción de demanda y generación renovables con Deep Learning: Aplicación a la optimización de estaciones de carga de vehículos eléctricos [19]. 4.1 Red neuronal Las redes neuronales son modelos que tratan de emular el funcionamiento del sistema nervioso, es decir, la forma en la que el cerebro procesa la información. La unidad básica es la neurona, que consiste en un número que suele ir de 0 a 1. Este valor indica si la neurona va a ser activada (valor 1) o no (valor 0). La neurona se conecta al resto mediante enlaces con una fuerza de conexión variable (también llamada ponderación), formando la red neuronal. Generalmente, esta compuesta por 3 capas: • Capa de entrada: En esta capa, cada neurona se corresponde con una variable de entrada. Reciben los datos y los envía al resto para ser procesada. • Capa oculta: Pueden ser una o varias según la complejidad del modelo y se encargan de procesar la información para luego propagarla hasta la capa siguiente. • Capa de salida: Las neuronas de esta capa se corresponden con posibles resultados. La que obtenga un valor de activación más alto será el resultado elegido por el algoritmo. Inicialmente, las ponderaciones son aleatorias, por lo que es preciso entrenar a la red neuronal para ajustar estos parámetros y obtener un resultado lo más preciso posible. Para ello, se hace uso de datos de entrenamiento cuyo resultado final es conocido a priori. En cada prueba se compara el resultado de la red con el esperado y esta información se envía desde la capa de salida hasta la de entrada para modificar las ponderaciones de forma gradual. El proceso se repite hasta alcanzar alguno de los criterios de parada. En ese momento, la red se considera entrenada y lista para realizar predicciones suficientemente precisas. 21 22 Capítulo 4. Estrategias de predicción de la demanda y producción solar Figura 4.1 Esquema red neuronal básica [29]. 4.2 Estrategia de predicción Existen gran variedad de modelos de redes neuronales, por lo que se ha decidido testear los que se listan a continuación y seleccionar el más preciso, es decir, el que tenga menor error absoluto medio (MAE) y menor error cuadrático medio (RMSE). Modelo neuronal simple ANN Es el modelo básico de perceptrón multicapa o red neuronal artificial con neuronas "totalmente conectadas" donde la información fluye en una dirección, desde la capa de entrada hasta la de salida. Modelo neuronal recurrente RNN El modelo recurrente utiliza neuronas que incluyen retroalimentaciones, es decir, en cada instante de tiempo, las neuronas reciben como entrada las señales procedentes de la capa anterior y su propia salida del instante anterior. Es apropiado para el tratamiento de secuencias temporales. Modelo neuronal con memoria a largo plazo LSTM Las LSTM es una variante del modelo recurrente en el que la información puede persistir durante un periodo de tiempo más prolongado, "recordando" estados previos y utilizando esa información para decidir el valor del estado siguiente. Modelo neuronal con memoria a largo plazo bidireccional Bi-LSTM Se trata del modelo LSTM con una modificación. En este caso, no solo se tiene en cuenta el valor de estados previos sino que además se puede recibir información de estados futuros, fluyendo esta en ambas direcciones. Modelo neuronal convolucional CNN Es un modelo utilizado para el procesamiento de imágenes principalmente. Se ha utilizado en este problema para realizar convoluciones unidimensionales sobre la secuencia de entrada con el fin de extraer características relevantes. 4.3 Predicción de la potencia demandada por los vehículos 23 4.3 Predicción de la potencia demandada por los vehículos Con el objetivo de conseguir una predicción acerca de la demanda de energía procedente de los vehículos eléctricos se han utilizado los datos explicados en el apartado 3.4 como datos de entrenamiento. En primer lugar, se han procesado y filtrado para que tengan la estructura necesaria y puedan ser usados como datos de entrenamiento. Tras entrenar los modelos, se han validado los resultados mediante su comparación con los datos reales, calculando el error absoluto medio y el error cuadrático medio. MAE =1 n n ∑ j=1 |yj−ˆyj| RMSE =v u u t n ∑ j=1 (ˆyj−y2 j) n Finalmente, se ha seleccionado el modelo con redes LSTM con ventana temporal de 4 días, con el que se han obtenido un MAE de 4.41 kW y un RMSE de 4.10 kW sobre los datos de testeo. Figura 4.2 Comparativa predicciones demanda vs datos de testeo. 4.4 Predicción de la potencia suministrada por los paneles fotovoltaicos Se ha hecho uso de los datos explicados en el apartado 3.2.2 para entrenar a los modelos, con el objetivo de conseguir una predicción sobre la generación de energía del sistema fotovoltaico. Tras realizar el mismo proceso explicado en el caso anterior, se ha seleccionado la combinación de dos modelos CNN+LSTM con ventana temporal de 1 día, con el que se han obtenido un MAE de 55.60 kW y un RMSE de 104.61 kW sobre los datos de testeo. Figura 4.3 Comparativa predicciones producción solar vs datos de testeo. 30 Capítulo 5. Problema de optimización para determinar la acción de control 5.5 Adaptación del modelo Como se ha podido observar, las ecuaciones que describen el comportamiento de la planta entremezclan variables binarias y variables continuas, apareciendo productos cruzados como, por ejemplo, en la ecuación 5.17. Para poder resolver el problema utilizando la herramienta GUROBI, es necesario expresar las ecuaciones como relaciones lineales, por lo que habrá que adaptarlas a un formato resoluble. Siguiendo las transformaciones MLD (Mixed logic dynamics) [ 2 ], se puede linealizar este tipo de ecuaciones mediante unas reglas sencillas, incluyendo en ciertos casos alguna variable auxiliar. 5.5.1 Transformación del producto cruzado de variables continuas y lógicas Dada la expresión y=f(x)δ donde f(x) es una función continua, lineal y δ∈ {0,1} es una variable binaria, es fácil comprobar que es equivalente a estas dos expresiones: δ=0=⇒y=0 δ=1=⇒y=f(x) Definiendo un valor máximo (M) y mínimo (m) para la función f(x) , estas expresiones se pueden reescribir de la siguiente forma, consiguiendo sustituir el producto cruzado por cuatro ecuaciones en términos de desigualdades lineales: y≤Mδ y≥mδ y≤f(x)−m(1−δ) y≥f(x)−M(1−δ) Las ecuaciones 5.10 y 5.17, pueden ser transformadas siguiendo esta regla. •Restricción (5.10): PGRID pur (tk) = PGRID(tk)δGRID pur (tk) Se sustituye por: PGRID pur (tk)≤PGRID max δGRID pur (tk) PGRID pur (tk)≥PGRID min δGRID pur (tk) PGRID pur (tk)≤PGRID(tk)−PGRID min (1−δGRID pur (tk)) PGRID pur (tk)≥PGRID(tk)−PGRID max (1−δGRID pur (tk)) •Restricción (5.17): PBES disc (tk) = PBES(tk)δBES disc (tk) Se sustituye por: PBES disc (tk)≤PBES max δBES disc (tk) PBES disc (tk)≥PBES min δBES disc (tk) PBES disc (tk)≤PBES(tk)−PBES min (1−δBES disc (tk)) PBES disc (tk)≥PBES(tk)−PBES max (1−δBES disc (tk)) 5.5 Adaptación del modelo 31 5.5.2 Transformación de la expresión bidireccional entre variables continuas y lógicas Dada la expresión f(x)≤0←→ δ=1 donde f(x)∈[m,M] es una función continua, lineal y δ∈ {0,1} es una variable binaria, es fácil comprobar que es equivalente a estas dos desigualdades lineales: f(x)≤M(1−δ) f(x)≥ε+(m−ε)δ εes un valor suficientemente pequeño pero no cero. Se suele tomar ε=10−3. Las ecuaciones 5.12 y 5.19, pueden ser transformadas siguiendo esta regla: •Restricción (5.12): PGRID(tk)≤0←→ δGRID sale (tk) = 1 Se sustituye por: PGRID(tk)≤PGRID max (1−δGRID sale (tk)) PGRID(tk)≥ε+(PGRID min −ε)δGRID sale (tk) •Restricción (5.19): PBES(tk)≤0←→ δBES ch (tk) = 1 Se sustituye por: PBES(tk)≤PBES max (1−δBES ch (tk)) PBES(tk)≥ε+ (PBES min −ε)δBES ch (tk) En esta tabla se resumen todas las transformaciones MLD para linealizar ecuaciones mixtas. Figura 5.2 Tabla de transformaciones MLD [16]. 32 Capítulo 5. Problema de optimización para determinar la acción de control 5.6 Programación 5.6.1 Programación de la función objetivo Recordando el formato que debe tener la función objetivo (Ec. 5.24) m´ ınJ(t) = xTQx +cTx+α , se procede expresar el problema de optimización siguiendo ese formato. Variable a optimizar x Engloba todas las variables cuyo valor se puede controlar. Se debe indicar si la variable es binaria o continua. Se trata de un vector columna de 7*SH componentes, donde k toma valores desde 1 hasta SH: x= [x1x2... xSH ]T Donde: xk= [PBES disc (tk|t)PBES ch (tk|t)PGRID sale (tk|t)δBES disc (tk|t)δBES ch (tk|t)δGRID pur (tk|t)δGRID sale (tk|t)]T Dentro de este vector se ha excluido el término PGRID pur (tk|t) , ya que es una variable dependiente del resto que se puede calcular realizando balance de potencia (Ec. 5.23). Sustituyendo este término por su expresión en la función objetivo (5.7), el problema de optimización se ve modificado de la siguiente forma: m´ ın uJ(t) = SH ∑ k=1 −Γpred sale (tk|t)Psale(tk|t)Ts+ Γpred pur (tk|t)(Psale(tk|t)−(PBES disc (tk|t)−PBES ch (tk|t))−Ppv(tk|t)+Pev(tk|t))Ts CCostBES 2CyclesBES (PBES ch (tk|t)+PBES disc (tk|t))Ts+ (CostBES degr,ch(PBES ch (tk|t))2+CostBES degr,disc(PBES disc (tk|t))2)Ts) (5.25) Término cuadrático Q El término cuadrático de la ecuación de optimización es el relacionado con el desgaste de la batería debido a la carga/descarga. Se trata de una matriz cuadrada de 7*SH elementos. Q=diag(Q1Q2... QSH ) Donde: Qk=               Costdegr,disc 2Ts0 00000 0Costdegr,ch 2Ts00000 0 0 00000 0 0 00000 0 0 00000 0 0 00000 0 0 00000               Término proporcional c Este vector debe incluir el resto de términos de la expresión. Será un vector columna de 7*SH elementos. c= [c1c2... cSH ]T Donde: ck= [( CCostBES 2CyclesBES −Γpred pur (tk|t))Ts(CCostBES 2CyclesBES +Γpred pur (tk|t))Ts(Γpred pur −Γpred sale )Ts0000]T 5.6 Programación 33 Término independiente α Al ser un término constante que no depende de las acciones de control, su valor es irrelevante para el cálculo de la minimización de la función J(t), por lo que no será necesario definirlo. 5.6.2 Programación de las restricciones de la planta Inicialmente, se debe hacer una recopilación de todas las restricciones del modelo, sustituyendo las que al comienzo eran no lineales por las adaptadas explicadas en el apartado 5.5. Además, se ha sustituido las variables PGRID , PBES y PGRID pur por sus respectivas ecuaciones (5.3, 5.4 y 5.23), para utilizar en las expresiones únicamente variables que aparecen en el vector x. •Restricciones de la red eléctrica: PGRID sale (tk|t)−(PBES disc (tk|t)−PBES ch (tk|t)≥PPV (tk|t)−PEV (tk|t) PGRID sale (tk)≥0 δGRID pur (tk)+ δGRID sale (tk) = 1 −PBES disc (tk|t) +PBES ch (tk|t)−(PGRID min −ε)δGRID sale (tk)≥ε+PPV (tk|t)−PEV (tk|t) PBES disc (tk|t)−PBES ch (tk|t)−PGRID max δGRID sale (tk)≥ −PGRID max −PPV (tk|t) +PEV (tk|t) −PGRID sale (tk|t) +PBES disc (tk|t)−PBES ch (tk|t)+PGRID max δGRID pur (tk)≥ −PPV (tk|t)+PEV (tk|t) PGRID sale (tk|t)−PBES disc (tk|t) +PBES ch (tk|t)−PGRID min δGRID pur (tk)≥PPV (tk|t)−PEV (tk|t) −PGRID sale (tk)+ PGRID min δGRID pur (tk)≥PGRID min PGRID sale (tk)−PGRID max δGRID pur (tk)≥ −PGRID max •Restricciones del sistema de almacenamiento: PBES disc (tk)≥0 PBES ch (tk)≥0 δBES disc (tk)+δBES ch (tk) = 1 PBES disc (tk)−PBES ch −(PBES min −ε)δBES ch (tk)≥ε −PBES disc (tk)+PBES ch −PBES max δBES ch (tk)) ≥ −PBES max −PBES disc (tk)+PBES max δBES disc (tk)≥0 PBES disc (tk)−PBES min δBES disc (tk)≥0 −PBES ch (tk) +PBES min δBES disc (tk)≥PBES min PBES ch (tk)−PBES max δBES disc (tk)≥ −PBES max • Otras restricciones: Las limitaciones del flujo de potencia, tanto de las baterías (5.21) como de la red eléctrica (5.14), fueron incluidas al adaptar las ecuaciones en el apartado 5.5, por lo que solo falta incluir la restricción referente al estado de carga de las baterías (5.22). SOCBES(tk)≥SOCBES min −SOCBES(tk)≥ −SOCBES max Al ser todas las restricciones lineales se va a aplicar el siguiente formato Ax =b o Ax ≥b , según sean ecuaciones o inecuaciones. El vector xes la variable a optimizar, explicada en el la sección anterior. 34 Capítulo 5. Problema de optimización para determinar la acción de control Como se puede observar, las dos últimas restricciones, las relacionadas con el estado de carga de las baterías, no están expresadas en función de x, por lo que para poder imponerlas habrá que adaptarlas. SOCBES(tk+i) = SOCBES(tk+i−1) +KBES ch PBES ch (tk+i−1)Ts−KBES disc PBES disc (tk+i−1)Ts(5.26) Siendo: KBES ch =µBES ch CapBES max KBES disc =1 µBES disc ∗CapBES max Se trata de una ecuación recursiva, debido a que el estado de carga en cada instante depende de los estados anteriores. El estado de carga del instante actual ( SOCBES(tk) ) es medible y el del instante posterior ( SOCBES(tk+1) ) se puede calcular conociendo las acciones de control del instante actual, por lo tanto, también es conocido. El resto se deberán ser calculando conforme se decidan las acciones de control de sus estados anteriores. SOCBES(tk+1) = SOCBES(tk) +KBES ch PBES ch (tk)Ts−KBES disc PBES disc (tk)Ts SOCBES(tk+2) = SOCBES(tk+1) +KBES ch PBES ch (tk+1)Ts−KBES disc PBES disc (tk+1)Ts SOCBES(tk+3) = SOCBES(tk+2) +KBES ch PBES ch (tk+2)Ts−KBES disc PBES disc (tk+2)Ts =SOCBES(tk+1)+KBES ch PBES ch (tk+1)Ts−KBES disc PBES disc (tk+1)Ts +KBES ch PBES ch (tk+2)Ts−KBES disc PBES disc (tk+2)Ts ... (5.27) Matriz de restricciones A Se trata de una matriz de Nrest ∗SH filas, siendo Nrest el número de restricciones, y 7∗SH columnas, que incluye la parte de las restricciones que depende de las variables a optimizar. Se ha construido a trozos separando las restricciones referentes al estado de carga de las baterías (matriz ASOC ), ya que estas depende de los estados anteriores. A=                              A1000000... 0 0 0000000... 0 0 0A200000... 0 0 ASOC 000000... 0 0 0 0 A30000... 0 0 ASOC ASOC 00000... 0 0 000A4000... 0 0 ASOC ASOC ASOC 0000... 0 0 0000A50 0 ... 0 0 ASOC ASOC ASOC ASOC 000... 0 0 ... ... ... ... ... ... ... ... ... ... 0000000... 0ASH ASOC ASOC ASOC ASOC ASOC ASOC ASOC ... ASOC 0                              Donde: ASOC ="−Kdisc Kch 00000 Kdisc −Kch 00000# 5.6 Programación 35 Ak=                                           −1 1 1 0 0 0 0 0 0 1 0 0 0 0 −1 1 0 0 0 −(Pgrid min −ε)0 1−1 0 0 −Pgrid max 0 0 1−1−1 0 0 0 Pgrid max −1 1 1 0 0 0 −Pgrid min 0 0 −1 0 0 0 Pgrid min 0 0 1 0 0 0 −Pgrid max 1 0 0 0 0 0 0 0 1 0 0 0 0 0 1−1 0 0 −(Pbat min −ε)0 0 −1 1 0 0 −Pbat max 0 0 −1 0 0 Pbat max 0 0 0 100−Pbat min 0 0 0 0−1 0 Pbat min 0 0 0 010−Pbat max 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 1 1                                           Vector de restricciones b Se trata de un vector columna de Nrest ∗SH elementos que incluye la parte de las restricciones que no depende de las variables de interés. Como en la matriz A, se han separado las restricciones del estado de carga. b= [b1bSOC b2bSOC ... bSH bSOC]T Donde: bk=                                           Ppv(tk|t)−Pev(tk|t) 0 ε+Ppv(tk|t)−Pev(tk|t) −Pgrid max −Ppv(tk|t)+Pev(tk|t) −Ppv(tk|t)+Pev(tk|t) Ppv(tk|t)−Pev(tk|t) Pgrid min −Pgrid max 0 0 ε −Pbat max 0 0 Pbat min −Pbat max 1 1                                           36 Capítulo 5. Problema de optimización para determinar la acción de control bSOC ="SOCbat min −SOCbat (tk+1) −SOCbat max +SOCbat (tk+1)# Finalmente, habrá que indicar cuales de las restricciones son de igualdad ( = ) o de desigualdad. En este caso, todas las restricciones de desigualdad son del tipo (≥). 6 Algoritmo de gestión de la energía basado en reglas heurísticas A modo de comparativa, se ha construido otro algoritmo para calcular la acción de control basado en reglas heurísticas, siguiendo el propuesto en la tesis Charging Station for Electric Vehicles Using Hybrid Sources (Eltoumi,2020,p84) [6]. El objetivo es maximizar el uso de la energía renovable, como ocurría en el caso anterior, y reducir el consumo de energía de la red eléctrica, reduciendo con ello la inestabilidad y la contaminación. Los pasos seguidos por el algoritmo son: 1. Se comprueba la carga inicial que tiene el vehículo, la potencia fotovoltaica disponible y el estado energético del sistema de almacenamiento. 2. Las baterías proporcionan energía al vehículo o demandan energía al sistema para cargarse, dependiendo del SoC de la misma y de la producción solar. 3. La red se utiliza para mitigar el déficit de energía para la carga del vehículo y la carga del sistema de almacenamiento. 4. En caso de que haya un excedente de energía procedente del sistema fotovoltaico, la red o las baterías podrán absorber la energía sobrante. Los datos de los que se parte para determinar la acción de control en un instante determinado son las predicciones de demanda de energía de los vehículos y de generación del sistema fotovoltaico para ese instante, además de las medidas realizadas al sistema en el instante anterior. El algoritmo propuesto utiliza estrategias basadas en reglas, principalmente formulando escenarios de "si entonces", es decir, en función del escenario en que se sitúe el sistema, se realizarán uno o varios modos de funcionamiento. • Modo 1 (PV2EV): La potencia generada por las placas fotovoltaicas será utilizada para abastecer la demanda de los vehículos. •Modo 2 (BES2EV): Las baterías de la planta se descargan para suministrar energía a los vehículos. •Modo 3 (GRID2EV): Los vehículos se cargan comprando energía a la red eléctrica. • Modo 4 (PV2BES): Las baterías de la planta se cargan con la energía procedente del sistema fotovoltaico. •Modo 5 (PV2GRID): Se vende la energía sobrante de la producción solar a la red eléctrica. • Modo 6 (GRID2BESS): Las baterías de la planta se cargan con energía procedente de la red eléctrica. Una vez determinados los modos de funcionamiento, se decide la acción de control a aplicar sobre las baterías ( Pbat ). La energía que se consume o vende a la red eléctrica, se decide posteriormente utilizando la ecuación de balance de potencia (5.1). 37 38 Capítulo 6. Algoritmo de gestión de la energía basado en reglas heurísticas 6.1 Funcionamiento en situación de sobrecarga Figura 6.1 Esquema funcionamiento en situación de sobrecarga. Se trata de la situación en la que la potencia suministrada por las placas fotovoltaicas no sea suficiente para abastecer la demanda de los vehículos. En este caso, se debe decidir de donde se obtiene la energía restante para cargar los vehículos. Para ello, se comienza comprobando el precio de la energía de la red eléctrica en ese instante de tiempo. → Precio por debajo del umbral: Si ocurre esta situación, la energía restante se obtendrá de la red eléctrica. Además, si el sistema de almacenamiento no está totalmente cargado ( SOCmax ), se carga. Con esto se aprovecha el bajo coste de la energía y se almacena para un posible momento de necesidad posterior. Pbat =−(SOCbat max −SOC(tk))CCbat (6.1) → Precio por encima del umbral: En este caso se intenta obtener la energía del sistema de almacenamiento. →Si las baterías tienen carga: (SOC(tk)>SOCmin) → Si la energía de las baterías no es suficiente para cargar los vehículos: La energía restante se obtendrá de la red eléctrica, descargando las baterías de la planta hasta su valor mínimo (SOCmin). Pbat = (SOC(tk)−SOCbat min)CCbat (6.2) → Si la energía de las baterías es suficiente para cargar los vehículos: Los vehículos se cargarán por completo utilizando solo la energía de las baterías. Pbat =Pev −Ppv (6.3) → Si las baterías están descargadas: ( SOC(tk)≤SOCmin ) Los vehículos se cargarán comprando energía a la red eléctrica. Pbat =0(6.4) 6.2 Funcionamiento en situación de subcarga 39 6.2 Funcionamiento en situación de subcarga Figura 6.2 Esquema funcionamiento en situación de subcarga. Se trata de la situación en la que la potencia suministrada por las placas fotovoltaicas es suficiente para abastecer la demanda de los vehículos. En este caso, se debe decidir qué se realiza con la energía sobrante. Para ello, primero se comprueba si el sistema de almacenamiento está cargado. →Si las baterías no están cargadas completamente: Se utiliza el excedente para cargarlas. →Energía procedente de las placas fotovoltaicas no es suficiente para cargar las baterías: → Precio por debajo del umbral: Se cargarán las baterías hasta su valor máximo, aprovechando el bajo coste de la red. Pbat = (SOC(tk)−SOCbat min)CCbat (6.5) → Precio por encima del umbral: Se cargarán las baterías con la energía sobrante procedente del sistema fotovoltaico. Pbat =Pev −Ppv (6.6) → Energía procedente de las placas fotovoltaicas suficiente para cargar las baterías: Se cargarán hasta su valor máximo y la energía sobrante se venderá a la red eléctrica. Pbat =−(SOCbat max −SOC(tk))CCbat (6.7) → → Si las baterías están totalmente cargadas: Si ocurre esta situación, la energía sobrante se vende a la red eléctrica, obteniendo un beneficio económico. Pbat =0; (6.8) Una vez decidida la acción de control a aplicar sobre la potencia de las baterías, se debe comprobar que este valor se encuentra dentro de los límites establecidos (restricción 5.21). En caso contrario, se satura al valor máximo o mínimo, según corresponda. Finalmente, se calcula Pgrid realizando balance de potencia. Pgrid =Pev −Ppv −Pbat 46 Capítulo 7. Resultados obtenidos y conclusiones Potencia absorbida o cedida por el sistema de almacenamiento de energía La potencia positiva se corresponde con la potencia de descarga de las baterías, mientras que la potencia negativa con la carga de las mismas, como se puede ver en la figura donde se analiza el estado de carga del sistema de almacenamiento (figura 7.11). Feb 01 Feb 02 Feb 03 Feb 04 Tiempo 2019 -4 -3 -2 -1 0 1 2 3 4 5 Potencia (kW) Potencia BES Figura 7.10 Intercambio de potencia del sistema de almacenamiento de energía. En este caso, el sistema recarga las baterías aprovechando el precio reducido de la red eléctrica durante la madrugada y las descarga cuando hay picos de precio en la red. El resto del tiempo, es más económico comprar la energía directamente a cargar y descargar las baterías debido a dos razones: • Al no haber excedente de energía procedente del sistema fotovoltaico, la energía para cargar las baterías también debería comprarse a la red eléctrica. • En cada ciclo de carga-descarga, se produce una degradación y depreciación del sistema de almacenamiento que ha sido tenida en cuenta en la función de optimización que se trata de minimizar. Feb 01 Feb 02 Feb 03 Feb 04 Tiempo 2019 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Porcentaje (%) SoC BES Figura 7.11 Evolución del estado de carga del sistema de almacenamiento de energía. 7.1 Gestión de la energía basado en técnica de optimización matemática 47 Potencia absorbida o cedida por la red eléctrica La potencia positiva se corresponde con la compra a la red eléctrica, mientras que la negativa se corresponde con la venta. Feb 01 Feb 02 Feb 03 Feb 04 Tiempo 2019 -20 -15 -10 -5 0 5 10 15 Potencia (kW) Potencia GRID Figura 7.12 Intercambio de potencia con la red eléctrica. En este escenario, la energía que demandan los vehículos eléctricos es superior a la que es producida por las placas fotovoltaicas, por lo que la mayoría del tiempo se está comprando energía a la red eléctrica. Balance de potencia Se han representado los intercambios de energía producidos en la estación de carga de vehículos eléctricos para comprobar el balance de potencia. Feb 01 Feb 02 Feb 03 Feb 04 Tiempo 2019 -20 -15 -10 -5 0 5 10 15 20 25 Potencia (kW) Balance de potencia Potencia PV Potencia EV Potencia BES Potencia GRID Figura 7.13 Balance de potencia. Se puede observar que los picos de precio de la energía eléctrica suelen ocurrir en las últimas horas del día, ya que es cuando se suelen descargar las baterías. 48 Capítulo 7. Resultados obtenidos y conclusiones Función de optimización Finalmente, se muestra el cálculo de la función de optimización a minimizar (Ec. 5.7). El resultado positivo es muestra de que no se obtiene beneficio económico de la instalación ya que no se genera suficiente energía en el sistema fotovoltaico para abastecer la demanda de los vehículos eléctricos. Feb 01 Feb 02 Feb 03 Feb 04 Tiempo 2019 0 5 10 15 20 25 Precio (kWh) Función objetivo Figura 7.14 Representación del cálculo de la función de optimización. 7.1.3 Comparativa con predicciones óptimas Con el objetivo de conocer cuál es el impacto que provoca el uso de predicciones para determinar la acción de control, se ha repetido el segundo escenario suponiendo que se conoce la potencia que va a ser generada por el sistema fotovoltaico y la demanda por parte de los vehículos con exactitud. Para ello, se han sustituido los datos de las predicciones por los datos reales, obteniendo el siguiente resultado. Feb 01 Feb 02 Feb 03 Feb 04 Tiempo 2019 0 5 10 15 20 25 Precio (kWh) Función objetivo Predicciones Datos reales Figura 7.15 Comparativa función de optimización con predicciones óptimas. Como se puede observar, conforme más precisas sean las predicciones, menor será el coste de la función objetivo debido a que se decidirá una acción de control más precisa. 7.2 Gestión de la energía basado en reglas heurísticas 49 7.2 Gestión de la energía basado en reglas heurísticas Para poder comparar la técnica de optimización matemática explicada en el apartado 5 y la técnica heurística definida en el apartado 6, se ha repetido el escenario 1, pero ahora se decide la acción de control utilizando las reglas heurísticas. Finalmente, se compara la función de optimización a minimizar (Ec. 5.7) obtenida en cada caso. La gráfica de la potencia generada por el sistema fotovoltaico es la correspondiente a la figura 7.1 y la potencia demandada por los vehículos eléctricos la figura 7.2. Para que los resultados sean similares, en este caso tampoco se aplica acción de control sobre el sistema hasta pasadas 24 horas. Potencia absorbida o cedida por el sistema de almacenamiento de energía La potencia positiva se corresponde con la potencia de descarga de las baterías, mientras que la potencia negativa con la carga de las mismas. May 01 May 02 May 03 May 04 Tiempo 2019 -6 -4 -2 0 2 4 6 Potencia (kW) Potencia BES Figura 7.16 Intercambio de potencia del sistema de almacenamiento de energía. May 01 May 02 May 03 May 04 Tiempo 2019 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Porcentaje (%) SoC BES Figura 7.17 Evolución del estado de carga del sistema de almacenamiento de energía. 50 Capítulo 7. Resultados obtenidos y conclusiones Al utilizar las reglas heurísticas no se tiene en cuenta la degradación y depreciación que sufre el sistema de almacenamiento en cada ciclo de carga-descarga. Esto se traduce en un mayor uso del sistema de almacenamiento, como se ve en la figura 7.16. Además, este algoritmo solo tiene en cuenta el estado actual para decidir la acción de control a aplicar, por lo que se intenta tener las baterías cargada para anticiparse a un posible contratiempo. Potencia absorbida o cedida por la red eléctrica La potencia positiva se corresponde con la compra de potencia a la red eléctrica, mientras que la negativa se corresponde con la venta. Los resultados son bastante similares a los obtenidos en el caso anterior. May 01 May 02 May 03 May 04 Tiempo 2019 -25 -20 -15 -10 -5 0 5 10 15 Potencia (kW) Potencia GRID Figura 7.18 Intercambio de potencia con la red eléctrica. Balance de potencia La principal diferencia apreciable en el balance de potencia respecto al anterior es un mayor intercambio de energía procedente del sistema de almacenamiento. May 01 May 02 May 03 May 04 Tiempo 2019 -25 -20 -15 -10 -5 0 5 10 15 20 25 Potencia (kW) Balance de potencia Potencia PV Potencia EV Potencia BES Potencia GRID Figura 7.19 Balance de potencia. 7.2 Gestión de la energía basado en reglas heurísticas 51 Función de optimización Finalmente, se muestra el cálculo de la función de optimización a minimizar (Ec.5.7) obtenido en el caso de utilizar la técnica de optimización matemática y la técnica de optimización basada en reglas heurísticas. El resultado es bastante similar aunque se observa una pequeña diferencia al principio del día 2 de mayo, donde el coste es mayor si se utiliza la técnica heurística. May 01 May 02 May 03 May 04 Tiempo 2019 -50 -45 -40 -35 -30 -25 -20 -15 -10 -5 0 Precio (kWh) Función objetivo Optimización MIQP Técnica heurística Figura 7.20 Comparativa función de optimización con reglas heurísticas. Se ha decidido hacer una última prueba sustituyendo los datos de las predicciones por los datos reales para ver como cambia el comportamiento de cada técnica. May 01 May 02 May 03 May 04 Tiempo 2019 -60 -50 -40 -30 -20 -10 0 Precio (kWh) Función objetivo Optimización MIQP Técnica heurística Figura 7.21 Comparativa función de optimización con reglas heurísticas y datos reales. Se puede comprobar que conociendo los datos reales, la técnica matemática da mejores resultados que la heurística, ya que selecciona las acciones de control que consiguen extraer mayor beneficio a la planta. 52 Capítulo 7. Resultados obtenidos y conclusiones 7.3 Conclusiones En este trabajo se ha realizado un diseño a nivel matemático basado en el balance de potencia de una estación de carga de vehículos eléctricos formada por un sistema fotovoltaico, un sistema de almacenamiento de energía y la red eléctrica. Posteriormente, se han explicado dos algoritmos inteligentes de gestión de flujo de potencia con el objetivo de minimizar el coste monetario de la estación de carga. El primero empleaba técnicas de optimización matemática, utilizando la programación cuadrática de enteros mixtos. El segundo se basaba en una serie de reglas heurísticas. Como se mencionó en apartados anteriores, los algoritmos heurísticos no garantizan la obtención de la solución óptima, sino una solución buena que cumple las especificaciones. En cambio, los algoritmos matemáticos sí garantizan un resultado óptimo, aunque sean más lentos y tengan un coste computacional mayor. Tras realizar pruebas con ambos modelos, se cumplieron las expectativas. • Coste computacional: La simulación de tres días en el primer modelo podía demorarse varios minutos, mientras que en el segundo caso los resultados se calculaban en cuestión de segundos. • Mejor solución: Si se observa la figura 7.21, se puede apreciar cómo el método matemático consigue un coste monetario menor respecto al método heurístico, aunque la precisión en las predicciones puede alterar este resultado. En definitiva, la técnica de optimización matemática ofrece mejores resultados que la técnica heurística, a expensas de aumentar el coste computacional. Si se desea realizar pruebas de varios meses de duración, el tiempo de ejecución puede alargarse durante horas. En contraposición, el método heurístico devuelve buenos resultados de forma considerablemente más rápida. 8 Lineas futuras La tecnología V2G (vehicle to grid) permite un flujo de energía bidireccional entre el vehículo y la red eléctrica, convirtiendo la batería del vehículo en una fuente de energía controlable con la finalidad de satisfacer picos de demanda de la red. Para hacer uso de esta tecnología, el propietario del vehículo debe estar de acuerdo con el aprovechamiento de su batería como suministro de energía, ya que existen algunas desventajas desde el punto de vista del usuario. •El tiempo de carga de vehículo será más extenso de lo habitual. • La batería del vehículo se desgasta con cada ciclo de carga/descarga, por lo que su vida útil se verá reducida. En consecuencia, se necesitan incentivos y descuentos para motivar a los conductores a aceptar el uso de esta tecnología. Se presupone que la estación de carga cuenta con un gestor que hace de intermediario entre los proveedores de la red eléctrica y los usuarios, tratando de negociar unas tarifas que beneficien a ambas partes. El mecanismo de tarificación se divide en 3 etapas: 1. Existe un mercado encargado de gestionar las desviaciones no programadas con el fin de buscar un suministro de energía que ayude a estabilizar la red eléctrica. Cuando ocurre algún pico en la demanda, se envía un evento de demanda/respuesta (DR). 2. Cada electrolinera gestiona su conjunto de vehículos eléctricos cuya tecnología V2G esté habilitada y su sistema de almacenamiento, para ser capaz de reaccionar ante un evento DR. Cuando sucede, debe proporcionar una respuesta con el objetivo de maximizar el beneficio económico de la instalación. 3. Es aconsejable que los propietarios de los vehículos tengan acceso a una plataforma donde puedan indicar la hora de salida y el estado de carga deseado del vehículo. A partir de esta información será más sencillo diseñar una estrategia de control de la energía para minimizar costes y satisfacer al usuario. Se plantea la posibilidad de incluir esta tecnología en la estación de carga explicada en este documento. Para ello, sería necesario modificar el modelo matemático (apartado 3) y el problema de optimización (apartado 5) para tener en consideración la nueva fuente de energía de la electrolinera [26]. 8.1 Modelado de los vehículos eléctricos Se define la potencia intercambiada con los vehículos como potencia de carga menos potencia de descarga, ya que el flujo de potencia se ha definido positivo cuando se desplaza hacía el vehículo (5.1). PEV (t) = PEV ch (t)−PEV disc(t)(8.1) Además, el estado de carga de los vehículos en cada instante de tiempo depende del estado de carga del instante anterior y de la potencia consumida en el instante actual. CEV (t+1) = CEV (t)−ηEV ch PEV ch (t)Ts−PEV disc(t) ηEV disc Ts(8.2) 53 54 Capítulo 8. Lineas futuras Se requiere definir nuevas variables para expresar de forma matemática el comportamiento de los vehículos con la tecnología V2G. Tabla 8.1 Datos tecnología V2G. Nombre Unidades Definición ηEV ch % Eficiencia de carga de la batería ηEV disc % Eficiencia de descarga de la batería ρDR % Probabilidad de que ocurra un evento DR ra- Factor de ponderación al precio de venta de la instalación ru- Factor de ponderación al precio de venta del dueño del vehículo πsale €/Wh Precio unitario de venta de energía que recibe la instalación πuser €/Wh Precio unitario de venta de energía que recibe el dueño del vehículo PEV max W Máxima potencia a la que se puede carga el vehículo PEV min W Máxima potencia a la que se desea vender a la red CEV max % Porcentaje máximo de carga de los vehículos CEV min % Porcentaje mínimo de carga de los vehículos El esquema de la estación de carga se vería modificado como se muestra en la figura, incluyendo el intercambio de energía bidireccional entre los vehículos y la red eléctrica. Figura 8.1 Esquema estación de carga de vehículos eléctricos con tecnología V2G. 8.2 Definición del problema de optimización El problema de optimización se ha planteado con el objetivo de minimizar el coste monetario de la planta. La función objetivo debe incluir todas las variables controlables, que en este caso serán la potencia intercambiada con el sistema de almacenamiento, con la red eléctrica y con los vehículos (tecnología V2G). m´ ın uJ(t) = SH ∑ k=1 (JGRID(tk|t)+JBES(tk|t)+JV2G(tk|t)) 8.3 Condiciones de factibilidad de la solución 55 El término JV2G(tk|t) está formado por dos partes: la primera se relaciona con los beneficios derivados de la venta de energía a la red eléctrica que recibe la estación de carga, mientras que la segunda refleja la compensación monetaria que recibe el dueño del vehículo por permitir el uso de la batería de su vehículo como suministro de energía. JV2G(tk|t) = −ρDR(tk|t)raπsale(tk|t)(PGRID sale (tk|t) +PEV disc(tk|t))Ts+ρDR(tk|t)ruπuser(tk|t)PEV disc(tk|t)Ts La probabilidad de que ocurra un evento RD en un instante de tiempo concreto se supone un valor conocido a priori, ya que está relacionado con los picos de demanda de energía por parte de los vehículos eléctricos. 8.3 Condiciones de factibilidad de la solución Para garantizar que el problema de optimización sea resoluble en todo el espectro de estados admisibles de la instalación, deben verificarse nuevas restricciones además de las ya comentadas en el apartado 5.3. En esta ocasión las expresiones son lineales, por lo que no será necesario adaptarlas. • El flujo de potencia que suministra o absorbe el vehículo debe estar dentro de los límites establecidos: PEV min ≤PEV (tk)≤PEV max (8.3) • El estado de carga agregada de todos los vehículos conectados a la estación de carga debe estar dentro de los límites: CEV min ≤CEV (tk)≤CEV max (8.4) 8.4 Planteamiento del problema de optimización El procedimiento de optimización puede ser el mismo explicado en el apartado 5.4, la programación cuadrática de enteros mixtos (MIQP). Al programar la función objetivo será necesario añadir las nuevas variables al vector x , que contiene todas las variables cuyo valor se puede controlar, que en este caso serían PEV disc yPEV ch . Además, la expresión JV2G(tk|t)deberá ser incluida en el término proporcional c. Finalmente, las nuevas restricciones deben ser añadidas a la matriz y el vector de restricciones, siendo necesario expresar el estado de carga de los vehículos de igual forma que se hizo con el estado de carga del sistema de almacenamiento (Ec 5.27), ya que ambos son ecuaciones recursivas. 62 Bibliografía [17] J Pandya, Solving np complete timetable problem using monte carlo and simulated annealing algorithms, https://jatinp101.github.io/projects/timeTabling. [18] Enel Green Power, Sistemas de almacenamiento con baterías de litio, https://www.enelgreenpower.com/es/learning-hub/energias-renovables/almacenamiento/baterias-de- litio. [19] F Ramos, Predicción de demanda y generación renovables con deep learning: Aplicación a la optimización de estaciones de carga de vehículos eléctricos, Universidad de Sevilla (2022). [20] E Rashedi, H Nezamabadi-pour, and S Saryazdi, Gsa: A gravitational search algorithm, Information Sciences 179 (2009), no. 13, 2232–2248, Special Section on High Order Fuzzy Sets. [21] Redeia, Red eléctrica de españa, https://www.ree.es/es/datos/generacion/estructura-generacion. [22] NOAA research, Noaa solar calculator, https://gml.noaa.gov/grad/solcalc/. [23] C Robles, Optimización por colonia de hormigas: aplicaciones y tendencias, (2010). [24] Tesla, Supercargador, https://www.tesla.com/es_es/supercharger, (2022). [25] Alternative Energy Tutorials, Solar irradiance, https://www.alternative-energy-tutorials.com/solarpower/solar-irradiance.html, (2021). [26] C Vivas, F R Rubio, A López, and F Ramos, Optimized operation of an electric vehicle charging station with photovoltaic support and vehicle-to-grid implementation, Universidad de Sevilla (2022). [27] S Wang, S Bi, Y Zhang, and J Huang, Electrical vehicle charging station profit maximization: Admission, pricing, and online scheduling, IEEE Transactions on Sustainable Energy 9(2018), 1722–1731. [28] Wikipedia, Estación de carga, https://es.wikipedia.org/wiki/Estación_de_carga, (2022). [29] , Red neuronal artificial, https://es.wikipedia.org/wiki/Red_neuronal_artificial, (2022). [30] L Zachary, L Tongxin, and L Steven, ACN-Data: Analysis and Applications of an Open EV Charging Dataset, Proceedings of the Tenth International Conference on Future Energy Systems, e-Energy ’19, June 2019. [31] S Zlobec, Linear and nonlinear optimization (richard w. cottle and thapa n. mukund) springer 2017, 12 2017.