El corte transversal de cada sección está dado por lo que se conoce como perfil aerodinámico, del que dependen las fuerzas aerodinámicas de sustentación 𝐿 y arrastre 𝐷, así como sus respectivos coeficientes de sustentación y arrastre 𝐶 𝑙 y 𝐶𝑑, respectivamente. Entre secciones se produce una variación de estos perfiles buscando hacer más eficiente el trabajo del álabe, esto se logra escogiendo los perfiles según la función que cumplen en cada una de las secciones.

La raíz de la pala cumple una función estructural, por lo que para esta se necesita un perfil con un área considerable que soporte altos esfuerzos. Por otro lado, la sección interna, aunque debe soportar esfuerzos, cumple un papel de generador de energía a muy bajas velocidades, y se hace necesario un perfil con un 𝐶𝑙 alto para estas condiciones. La parte más larga del álabe es la parte central y se encarga de sacar el mayor provecho a la energía del viento. Una de las características de sus perfiles es que tienen una relación de 𝐶𝑙/𝐶𝑑 alta, esto con el fin de generar mayor energía, además de que estos perfiles tienen un ángulo elevado de entrada en pérdida, ya que, al ser la parte más larga del álabe, deben poseer una región amplia en la que puedan variar los ángulos de incidencia de sus perfiles. Por último, la parte externa que, además de generar energía, debe disminuir las pérdidas producidas por diferencias de presiones al final de la punta. Sus perfiles poseen elevados coeficientes de sustentación, pero bajos coeficientes de arrastre para disminuir las perdidas mencionadas.

Bertagnolio [13] proporciona algunas recomendaciones de perfiles NACA serie 6 utilizados generalmente en aerogeneradores, además de gráficas de coeficientes de sustentación y arrastre en función de su ángulo de incidencia. Esta serie contiene 6 números para su caracterización. Como ejemplo, tenemos el NACA 63- 215, donde el primer número, “6”, hace referencia a la serie del perfil; el segundo número, “3”, dice la localización de la presión mínima en función de la cuerda; el tercer número, “2”, es el coeficiente de sustentación de diseño, los últimos dos dígitos, “15”, representan el máximo espesor de la cuerda. Los perfiles recomendados para aerogeneradores de eje horizontal son los siguientes:

• NACA 63-215, NACA 63-218, NACA 63-221.

• NACA 63-415, NACA 63-418, NACA 63-421.

• NACA 64-415, NACA 64-421.

• NACA 65-415, NACA 65-421.

Se realizaron las polares de cada uno con ayuda del programa Javafoil [14], de donde se consiguieron los valores de 𝐶𝑙 y 𝐶𝑑 para diferentes ángulos de ataque α (entre 0° y 30°) y Reynolds (entre 50000 y 100000). Con estos valores se escogieron manualmente los perfiles teniendo en cuenta las diferentes características para cada sección. Para la sección interna se obtuvo el perfil NACA 64-421 que tiene un alto 𝐶𝑙 a bajos Reynolds y α (5° a 10°), teniendo además la mayor área con respecto a los demás perfiles. Para la sección media se eligió el NACA 63-421 con el 𝐶𝑙 más alto a diferentes Reynolds y α. Por último, en la sección externa se decidió por el NACA 63- 221 por tener un bajo 𝐶𝑑 a altos Reynolds y α (15° a 30°).

Para la optimización de perfiles se tomó como punto de partida los perfiles escogidos para cada sección del álabe. En cada uno se optimizó un parámetro diferente, dependiendo de su función ya descrita, mejorando en el perfil interior su 𝐶𝑙 a bajas velocidades (7m/s), en el medio una mejor relación entre el 𝐶𝑙 y 𝐶𝑑 y, por último, en el exterior un menor 𝐶𝑑 sin disminuir su Cl.

En la optimización de estos perfiles se usó un método de sensibilidad con respecto a la geometría mediante el módulo Adjoint Solver del software ANSYS [15]. En este caso la geometría es el perfil y la sensibilidad está dada por la relación del perfil con las fuerzas aerodinámicas. Los cambios en el perfil son producidos en los puntos de mayor sensibilidad con el fin de disminuir o aumentar estas fuerzas y que su cambio con respecto al cambio geométrico sea mínimo, de modo que se obtenga un perfil optimizado.

El método BEM es muy utilizado en la actualidad para el desarrollo de la geometría de hélices. Este se fundamenta a partir de la relación de las teorías de momento y del elemento de pala, haciendo uso de parámetros referentes a los perfiles aerodinámicos, la geometría de la turbina (número de aspas, envergadura de los álabes) y el viento (velocidad, viscosidad cinemática).

2.2.1. Teoría de momento

Una manera de explicar la teoría de momento es mediante el modelo de Rankine-Froude [7], [15]. Esta toma la hélice como un disco infinitamente delgado y divide el campo de flujo en 4 secciones, mostradas en la figura 2. En la primera sección, el viento entra con velocidad 𝑉1, presión 𝑃1 y en estado laminar. A medida que avanza hacia la sección 2 este pierde velocidad hasta alcanzar 𝑉2. En este punto atraviesa la hélice hasta la sección 3, y no se produce un cambio significativo en la velocidad, por lo que 𝑉2 = 𝑉3 = 𝑉, pero sí se produce una variación de presión debido a la energía obtenida por la hélice 𝑃2 ≠ 𝑃3. El viento continúa avanzando hasta salir por la sección 4 con una velocidad 𝑉4; aquí las presiones se igualan con la presión de entrada 𝑃1 = 𝑃4.

Figura 2.
Modelo de Rankine-Froude.
elaboración propia.

Se utiliza el teorema de Bernoulli para analizar el modelo, esto no se puede realizar a lo largo del campo de flujo debido al salto de presión producido por la hélice, así que se realiza entre las secciones 1-2 y 3-4, de lo que se obtienen las siguientes ecuaciones, donde 𝜌 es la densidad del aire:

(1)

(2)

Con las ecuaciones (1) y (2), y teniendo en cuenta las condiciones 𝑉2=𝑉3 y 𝑃1=𝑃4, se obtiene:

(3)

Teniendo las diferencias de presión en la hélice se logra hallar la fuerza de empuje ejercida sobre esta:

(4)

Además, por la ecuación de cantidad de movimiento axial y de continuidad de caudal:

𝐹𝑥=∫𝑉 𝜌 (𝑉 𝑑𝐴)0𝐴,

𝐹𝑥=𝑉4 𝜌 𝑉4 𝐴4−𝑉1𝜌 𝑉1 𝐴1,

Caudal=𝜌 𝑉4 𝐴4=𝜌 𝑉1 𝐴1=𝜌 𝑉2 𝐴,

Se obtiene la ecuación:

(5)

Con las ecuaciones (4) y (5) se obtiene:

(6)

Y al incluir un factor de inducción axial 𝑎:

(7)

Se obtiene la ecuación (8), que elimina la variable 𝑉4 de la ecuación (4):

(8)

(9)

Ahora se analiza la cantidad de momento angular para la hélice, con el fin de hallar el torque generado sobre esta.

Siendo 𝐿 el momento de inercia; 𝐼, el momento angular, y ω, la velocidad angular inducida por la hélice sobre el viento. Así que:

(10)

Al igual que en las fuerzas de empuje, se introduce un factor de inducción radial (𝑎′) a la ecuación (10) para eliminar la variable 𝜔:

(11)

Donde 𝛺 es la velocidad angular de la hélice. Recordando la ecuación (7) para eliminar 𝑉2:

(12)

Teniendo así las ecuaciones que gobiernan las fuerzas de empuje (9) y los torques (12), por medio de la teoría de momento.

Antes de seguir con la teoría de elemento de pala se debe agregar a estas ecuaciones un coeficiente de corrección de pérdidas 𝑄, debido que las fuerzas axiales y de torque son afectadas por los vórtices creados por las diferencias de presión en las puntas. Este coeficiente es estudiado por diferentes autores y, en este caso, se escoge el planteamiento de Prandtl:

(13)

Donde B es el número de álabes, R es la envergadura del álabe y 𝜙 es el ángulo de incidencia. El coeficiente se multiplica por las ecuaciones (9) y (12):

(14)

(15)

2.2.2. Teoría del elemento de pala

La teoría del elemento de pala divide la hélice en elementos diferenciales para hallar las fuerzas ejercidas sobre cada uno, sin tener en cuenta los efectos de las secciones contiguas. Por lo general, el álabe es dividido en 10 a 20 elementos, cada uno de estos asociado a una Siendo 𝐿 velocidad de flujo 𝛺𝑟, cuerda 𝑐 y ángulo de incidencia 𝜙 diferente, con el fin de aprovechar mejor la energía del viento.

El análisis comienza definiendo los ángulos de incidencia de cada sección. Para esto es necesario recordar que sobre la velocidad del aire que incide sobre el perfil 𝑉𝑤 se genera una velocidad angular 𝜔, producida por Ω y en dirección contraria a esta, como se muestra en la figura 3. Así que al analizar las velocidades tangenciales del viento que golpean el perfil se produce una velocidad relativa , mostrada en la figura 4. Junto con la velocidad axial del viento definen 𝜙 en la ecuación:3 Así que al analizar las velocidades tangenciales del viento que golpean el perfil se produce una velocidad relativa

(16)

Un factor importante en el diseño de álabes es la relación de velocidades locales de punta 𝜆𝑟 , la cual se hallará de manera gráfica más adelante, con ayuda del coeficiente de potencia y el tipo de rotor que se va a utilizar, el cual resulta de la ecuación anterior:

Las fuerzas 𝐿 𝑦 𝐷 se definen por sus coeficientes 𝐶𝑙 y 𝐶𝑑: 1

Figura 3.
Velocidades tangenciales de la hélice y el viento.
elaboración propia.

Figura 4.
Descomposición de las velocidades del viento sobre el perfil.
elaboración propia.

Ahora se procede a definir las fuerzas que actúan sobre el perfil, como se puede ver en la figura 5. Las fuerzas están asociadas a la dirección con la que el viento incide sobre el perfil (la dirección de la fuerza de sustentación L es perpendicular y de arrastre D horizontal). Descomponiéndolas tangencial y axialmente se tiene:

(17)

(18)

Las fuerzas 𝐿 𝑦 𝐷 se definen por sus coeficientes 𝐶𝑙 y 𝐶𝑑:

Figura 5.
Fuerzas sobre el perfil.
elaboración propia.

Remplazando 𝑑𝐿 𝑦 𝑑𝐷, y agregando el número de álabes 𝐵 en las ecuaciones (17) y (18) se obtiene:

(19)

(20)

Al multiplicar la fuerza tangencial por el radio se obtiene el torque generado 𝑑𝑇:

Remplazando 𝑉² 𝑒𝑛 𝑑𝑇 𝑦 𝑑𝐹 se tienen las ecuaciones que gobiernan las fuerzas de empuje y las fuerzas radiales, por medio de la teoría del elemento de pala:

(21)

(22)

𝜎′ es la solidez local

Figura 6.
Rendimiento de hélices de molinos y turbinas de viento.
[18].

2.2.3. Relación

El método BEM relaciona dichas teorías dividiendo las ecuaciones de fuerza axial y las de torques, con el fin de despejar 𝑎/(1 − 𝑎) de las axiales y 𝑎^′/(1 − 𝑎) de las tangenciales y tener las ecuaciones que definan los factores de inducción axial y tangenciales, con los cuales se hallarán la cuerda y ángulo de incidencia en cada sección:

En donde 𝐶_𝑃 es el coeficiente de potencia esperado (0.424).[17], 𝜂 es la eficiencia mecánica y eléctrica (0.9) [18] y V1 es la velocidad del viento. Así obtenemos la envergadura del álabe.

Con ayuda de las ecuaciones planteadas por el método BEM, se comienza con el desarrollo del álabe. Este es un proceso iterativo en el que se supone un valor para los factores de inducción axial 𝑎 y tangencial 𝑎′. Posteriormente, se desarrolla el método para el radio inicial de la sección interna del álabe y se halla 𝑐 𝑦 𝜙, con estos se recalcula 𝑎 𝑦 𝑎′, se comparan estos con sus valores anteriores, si sus porcentajes de error son altos se realiza nuevamente la iteración con los nuevos valores.

Por el contrario, si estos son bajos se realiza el procedimiento para el siguiente radio.

Antes de comenzar a hallar 𝑎 𝑦 𝑎′ se tendrán que definir algunos parámetros referentes a la geometría y el funcionamiento del generador (la envergadura de este 𝑅, el número de palas 𝐵 y su velocidad de rotación 𝜔). Como explica Ingram [16], la envergadura depende de la cantidad de potencia que se va a generar, y está dada por la siguiente ecuación:

En donde 𝐶𝑃 es el coeficiente de potencia esperado (0.424).[17], 𝜂 es la eficiencia mecánica y eléctrica (0.9) [18] y V1 es la velocidad del viento. Así obtenemos la envergadura del álabe.

El número de álabes es escogido de acuerdo con la figura 6, que está basada en prácticas experimentales. Teniendo 𝐶𝑝 y sabiendo que es una turbina HAWT, la gráfica nos indica que podemos utilizar dos o tres álabes (bipala o tripala). Aunque las turbinas de dos son más económicas, se escoge una de tres, debido a que esta presenta mayor estabilidad [19]. Asimismo, es posible hallar la relación de velocidades con la que podemos hallar Ω :

Para hallar un valor inicial de 𝑎 𝑦 𝑎′, Ingram [16] recomienda el uso de las siguientes ecuaciones:

Para el cálculo de la cuerda:

Desarrollando estas ecuaciones y la serie de iteracione se parametriza la geometría del álabe.

Se realiza el modelado del campo de flujo sobre los perfiles seleccionados, para esto se escoge una geometría C-mesh utilizada comúnmente en el análisis de perfiles aerodinámicos [20] y mostrada en la figura 7 para el perfil NACA 63-221.

Figura 7.
C-mesh del perfil NACA 63221.
elaboración propia.

Las condiciones de contorno utilizadas son las velocidades de entrada del aire. Estas tienen una magnitud y ángulo de ataque variante de un perfil a otro, debido a la velocidad angular de la pala y al ángulo de mayor aprovechamiento de energía, los cuales se observan en la tabla 1. Las velocidades de entrada actúan sobre los bordes inlet. Se tomó el borde airfoil como tipo Wall, y el borde outlet como Pressure-outlet. Estos se muestran la figura 8.

Tabla 1.
Velocidades y ángulos de ataque del viento sobre los perfiles de cada sección.

elaboración propia.

Figura 8.
Bordes del campo de flujo.
elaboración propia.

Para el mallado se realiza un análisis de independencia de malla a este mismo perfil, para ver la convergencia de la solución, teniendo como criterio 𝐶𝑙 y 𝐶𝑑, y un error mínimo de 2 %, los resultados se observan en la figura 9 y la tabla 2; además se puede apreciar cómo estos coeficientes se acercan a su valor real [13].

Tabla 2.
Convergencia de los coeficientes aerodinámicos para el perfil NACA 63-221.

elaboración propia.

Figura 9.
Convergencia de los coeficientes aerodinámicos para el perfil NACA 63-221.
elaboración propia.

Se usan ángulos de ataque de 5 grados, debido a que al usar ángulos mayores el 𝐶𝑑 del perfil aumenta, lo que produce una mayor fuerza de arrastre, que afecta el propósito de la turbina, además de que genera vibraciones indeseadas.

Teniendo la solución del campo de flujo, se realizó el proceso de optimización. En este se tomó 1 % como el valor de cambio entre cada una de las iteraciones, y se realizó un cálculo total cada 5 iteraciones, para tener un control sobre la deformación del perfil y verificar el cambio del valor de 𝐶𝑙 y 𝐶𝑑. Los resultados se muestran en la tabla 3.

Tabla 3.
Comparación de los coeficientes aerodinámicos con los perfiles optimizados. resultantes para los perfiles seleccionados y optimizados se muestran en la figura 10.

Se aprecia que los 𝐶𝑙 aumentaron para cada perfil según su condición de trabajo. Esta se tomó como un promedio entre secciones. Además, vemos que el 𝐶𝑑 para el NACA 63-221 disminuye como era requerido.

Los parámetros iniciales para el desarrollo del álabe fueron tomados de las turbinas Nordex N120 [17] de 1300 KW con una velocidad de viento promedio de 8 m/s [21]. Esta velocidad es la presente en el parque eólico Jepirachi [22], ubicado en la Guajira, Colombia.

Con estos valores de potencia y velocidad, se desarrolló el método BEM, con los perfiles seleccionados, teniendo como resultado la caracterización del álabe, dichas características se muestran en la tabla 4.

Para el centro de rotación de los perfiles se tomó una distancia de 0.30 veces de la cuerda a partir del borde de ataque, debido a que en dicho punto (línea de presión) inciden las fuerzas aerodinámicas [19]. También se implementó un larguero de soporte a lo largo del álabe a la misma distancia del borde de ataque. Los álabes resultantes para los perfiles seleccionados y optimizados se muestran en la figura 10.

Figura 10.
Álabes resultantes. (a) Perfiles seleccionados (morados) y optimizados (azules) y (b) estructura de los álabes.
elaboración propia.

3.2.1. Velocidades y torque sobre los álabes

Se realizó el campo de flujo sobre los álabes con un mallado de cono truncado, en el que solo se tomó un tercio de la turbina total y se analizó un solo álabe, debido a su simetría, con el fin de disminuir el gasto computacional. La malla se muestra en la figura 11.

Tabla 4.
Caracterización del álabe con perfiles seleccionados.

elaboración propia.

Figura 11.
Mallado del campo de flujo. (a) Vista isométrica y (b) vista en detalle de la superficie en contacto con el álabe.
elaboración propia.

Para las condiciones de contorno se tomó la velocidad de entrada de 8 m/s, que actúan sobre las superficies inlet e inletTop, el contorno de salida outlet es tipo outlet-pressure, la superficie en contacto con el álabe blade es tipo Wall, y por ultimo las zonas de interfaz o de simetría period1 y period2 son de tipo interface. Estos bordes se muestran en la figura 13, que outletpressure,la superficie en contacto con el álabe blade estipo Wall, y por ultimo las zonas de interfaz o de simetría period1 y period2 son de tipo interface. Estos bordes semuestran en la figura 13, que

La figura 13 muestra las velocidades tangenciales del álabe de los perfiles seleccionados, teniendo una velocidad en m/s máxima de 63 y mínima de 9.7, que concuerdan con la velocidad obtenida por el método BEM de 62 y 10, respectivamente, con un error relativo del 1.6 % y del 3 %, lo cual permite validar el mallado y las condiciones de contorno evaluadas anteriormente.

Al igual que con las velocidades tangenciales se obtuvieron y validaron los valores del torque generado para cada tipo de álabe, con el torque generado por una turbina de 1300 KW. Estos valores se pueden observar en la tabla 5. Como se observa, el torque generado por la turbina con los perfiles optimizados es mayor con respecto a la otra.

Tabla 5.
Torques generados por las turbinas.

elaboración propia.

Figura 12.
Bordes del campo de flujo.
elaboración propia.

La figura 13 muestra las velocidades tangenciales del álabe de los perfiles seleccionados, teniendo una velocidad en m/s máxima de 63 y mínima de 9.7, que concuerdan con la velocidad obtenida por el método BEM de 62 y 10, respectivamente, con un error relativo del 1.6 % y del 3 %, lo cual permite validar el mallado y las condiciones de contorno evaluadas anteriormente.

Figura 13.
Velocidades a lo largo del álabe.
elaboración propia.

Al igual que con las velocidades tangenciales se obtuvieron y validaron los valores del torque generado para cada tipo de álabe, con el torque generado por una turbina de 1300 KW. Estos valores se pueden observar en la tabla 5. Como se observa, el torque generado por la turbina con los perfiles optimizados es mayor con respecto a la otra.

Tabla 5.
Torques generados por las turbinas.

3.2.2 Esfuerzos y deformaciones sobre álabes

Para los álabes se utilizó un material ortótropo, utilizado comúnmente en la construcción de estas, y sus propiedades se muestran en la tabla 6. Los espesores seleccionados para los largueros y la superficie del álabe se encuentran en la tabla 7.

Tabla 6.
Propiedades del material del álabe.

elaboración propia.

Figura 14.
Esfuerzos sobre los álabes (a) original y (b) optimizado.

Tabla7.
Espesores para la superficie del álabe y para el larguero.

elaboración propia.