Full text
Equation Chapter 1 Section 1 Proyecto Fin de Grado Ingeniería Aeroespacial Diseño de álabes giratorios mediante simulaciones numéricas aplicado a una turbina eólica Autor: Alberto Damas Liébana Tutor: Javier Dávila Martín Dep. Ingeniería Aeroespacial y Mecánica de Fluidos Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2016
iii Proyecto Fin de Grado Ingeniería Aeroespacial Diseño de álabes giratorios mediante simulaciones numéricas aplicado a una turbina eólica Autor: Alberto Damas Liébana Tutor: Javier Dávila Martín Profesor titular Dep. Ingeniería Aeroespacial y Mecánica de Fluidos Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2016
v Proyecto Fin de Grado: Diseño de álabes giratorios mediante simulaciones numéricas aplicado a una turbina eólica Autor: Alberto Damas Liébana Tutor: Javier Dávila Martín El tribunal nombrado para juzgar el Proyecto arriba indicado, compuesto por los siguientes miembros: Presidente: Vocales: Secretario: Acuerdan otorgarle la calificación de: Sevilla, 2016
El Secretario del Tribunal
vii A mi familia A mis amigos
ix Agradecimientos Por fin el momento de los agradecimientos. Primero, agradecer a mi tutor por ofrecerme este interesante proyecto, así como a todos los profesores de la escuela, ya que de todos he aprendido algo. Como no podía ser de otra manera, agradecer a todos mis compañeros y amigos del grado, por todos esos cafés de hasta 1 hora filosofando de cualquier cosa. Agradecer también a mis compañeros de piso durante estos años, sobre todo, por aguantarme en los días malos. Y especialmente, agradecer los buenos ratos a esos tosirianos, que se enfadaban conmigo cuando tenía que estudiar. Por último, a los que quiero, aunque no se lo diga a menudo, agradecer a toda mi familia por haber hecho posible que haya llegado hasta aquí. A mi padre y mi madre por haber hecho de mi la persona que soy hoy, a mi hermano por tantísimos momentos juntos y, a mi hermana, que siempre será mi hermana pequeña. Alberto Damas Liébana Sevilla, 2016
Índice de figuras 1.1. Modelo energético europeo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.2. Instaladas y decomisadas en Europa (2015) . . . . . . . . . . . . . . . . . . . . . . . . 3 1.3. Energíaeólica......................................... 4 1.4. Potenciatotalinstalada ................................... 4 1.5. Potencia instalada durante 2015 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1. Perfil de velocidades para dos clases diferentes . . . . . . . . . . . . . . . . . . . . . . . 9 2.2. Abrigodelviento[1] ..................................... 10 2.3. Ejemplo de rosa de los vientos [1] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.4. Turbinasdeejehorizontal.................................. 12 2.5. Turbinasdeejevertical ................................... 13 2.6. Turbinasdeejevertical ................................... 14 2.7. Partes de una turbina eólica . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.1. Teoríadeldiscoactuador .................................. 20 3.2. Diferencia entre variantes realizable ystandard ...................... 24 4.1. Paso2............................................. 30 4.2. Paso4............................................. 30 4.3. Paso5............................................. 31 4.4. Paso6............................................. 32 4.5. Paso7............................................. 33 4.6. Paso8............................................. 33 4.7. Paso9............................................. 34 4.8. Paso10 ............................................ 34 4.9. Paso11 ............................................ 35 4.10.Paso12 ............................................ 36 4.11.Paso13 ............................................ 37 4.12.Paso14 ............................................ 37 4.13.Paso15 ............................................ 37 4.14.Paso16 ............................................ 38 4.15.Paso17 ............................................ 38 4.16.Paso18 ............................................ 39 III
4.17.Creacióndelcilindro..................................... 40 4.23.Mallaturbina......................................... 65 4.24.Mallalongitudinal ...................................... 65 4.25.Mallaperfil .......................................... 66 5.1. Perfilessimulados....................................... 72 5.2. Triángulodevelocidades................................... 72 5.3. Varianzadelrendimiento .................................. 74 5.4. Rendimiento ......................................... 75 5.5. Residuos............................................ 75 5.6. Númerodeceldas....................................... 77 5.7. Malla(Sección9)....................................... 78 5.8. DFBI ............................................. 81 5.9. Gradiente de velocidades adimensional (Sección 9) . . . . . . . . . . . . . . . . . . . . 82 5.10.Presión(Sección9)...................................... 83 5.11.Líneasdecorriente...................................... 85 5.12. Máximo módulo del gradiente de velocidades . . . . . . . . . . . . . . . . . . . . . . . 85 5.13.Potenciaobtenida ...................................... 86 5.14.Tiemposdesimulación.................................... 87 5.15.Memoriarequerida...................................... 88 5.18.Líneasdecorriente...................................... 95 5.19.Campodepresiones ..................................... 98 5.20.Memoriarequerida...................................... 98 5.21.Númerodeceldas....................................... 99 5.22.Tiemposdesimulación.................................... 100 5.23.................................................. 101 5.24.Líneasdecorriente...................................... 103 5.25.Campodepresiones ..................................... 104 5.26.................................................. 105 5.27.................................................. 106 5.28.................................................. 107 5.29.Líneasdecorriente...................................... 109 5.30.Campodepresiones ..................................... 110 5.31.................................................. 111 5.32.................................................. 112 5.33.Curvasderendimiento.................................... 113 5.35. Máximo gradiente de velocidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 5.36.Malla(Sección9)....................................... 117 5.37. Malla (Sección longitudinal) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 5.38.................................................. 119 5.39.Líneasdecorriente...................................... 121 IV
5.40.Campodepresiones ..................................... 123 5.41.................................................. 124 5.42.................................................. 125 5.43.................................................. 127 5.44.Malla(Sección9)....................................... 128 5.45.Valoresdelamalla...................................... 129 5.46.Líneasdecorriente...................................... 131 5.47.Campodepresiones ..................................... 132 5.48.................................................. 133 5.49.................................................. 134 5.50. Rendimiento-V elocidad angular .............................. 136 V
VI
Índice de tablas 1.1. Previsión ........................................... 2 2.1. Velocidades según rugosidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 5.1. Evolucióndelamalla .................................... 74 VII
VIII
Capítulo 1 Introducción La principal fuente de energía, de la cual proceden tanto los combustibles fósiles como las energías renovables, es el Sol. Existe una clara tendencia a aprovechar las energías renovables, frente al uso tradicionalmente mayoritario de los combustibles fósiles. Esto se debe a varios factores, como por ejemplo la contaminación que producen o el hecho de que las reservas se van agotando. Por otra lado, no suponen riesgo como la energía nuclear, ni producen ningún tipo de residuo. Por ello se van a arrojar algunos datos que ayuden a entender la importancia que tienen las energías renovables, concretamente la energía eólica ([5] y [7]). Dado que los ciclos para cambiar el modelo energético son bastante largos y que suponen una inversión inicial alta, es necesario considerar el largo plazo. El objetivo sería tener una red eléctrica globalizada, descentralizada y no contaminante, sin renunciar a la calidad, seguridad y a un continuo acceso a esta (no depender de conflictos entre países, sucesos meteorológicos, etc). Esto supone llegar a acuerdos internacionales importantes, así como a mantener un equilibrio que garantice la inversión sin caer en la falta de regulación que pueda derivar en monopolios. Cabe destacar dos niveles de alcance de cualquier proyecto de este tipo: uno local más centrado en la generación y el transporte más inmediato, y otro que asegure la conexión de diferentes puntos más distanciados. Esto supone como ya se ha comentado una cooperación importante. La organización Advisory Council on the Environment (SRU) publicó un informe en el que apostaba por un modelo energético basado única y exclusivamente en las energías renovables, argumentando que era factible para el año 2050. Por otro lado, el IPCC (Intergovernmental Panel on Climate Change) publicó un informe en el que evaluaba hasta 164 escenarios posibles. En referencia a estos objetivos, cabe mencionar que representantes de 194 países llegaron a un consenso: un 80% de la energía necesaria en 2050 podría ser renovable. En el 2000 la EWEA (European Wind Energy Asociation) publicó junto con Greenpace un informe denominado Windforce 10, el cual fue actualizando hasta que EWEA pasó a formar parte de GWEC (Global Wind Energy Council). En sus investigaciones estudiaron la energía eólica potencial para los próximos años. Manejaron tres posibles escenarios junto con Greenpace y la DLR (German Aerospace Centre): Se asume el escenario manejado por la IEA (International Energy Agency). Este escenario es el más conservador. Se asume que los países van a cumplir sus objetivos (teniendo en cuenta lo que llevan implementado hasta el momento). Se asume que los objetivos internacionales más avanzados se van a cumplir. Este escenario es el más optimista. Se muestra en la tabla 1.1 las previsiones para diferentes fechas desde su publicación partiendo de 1
CAPÍTULO 1. INTRODUCCIÓN Year Conservador (MW )Moderado (MW )Optimista(MW) 2009 158505 158505 158505 2010 185505 198717 201657 2015 295783 460364 533233 2020 415433 832251 1071415 2030 572733 1777550 2341984 Tabla 1.1: Previsión (a) 2000 (b) 2015 Figura 1.1: Modelo energético europeo lo que había entonces. Como fecha significativa hay que mencionar el año 2011, en el cual se superaron los 240 mil MW (24 ×1010 de potencia instalada a nivel mundial, de los cuales, más de 1000 MW estaban en 21 países. Por último se muestran algunos estadísticas que datan del año 2015 y a las que se puede acceder en [2]. En la figura 1.1 se puede observar el cambio del modelo energético entre el año 2000 y el 2015 en Europa. El mayor cambio es una disminución drástica de fuelóleo y un aumento de todas las energías renovables. Además, la energía nuclear ha disminuido ligeramente. Por otro lado, en la figura 1.2 se observan los cambios en el año 2015. Se pueden ver tanto las nuevas instalaciones de energías renovables como las desmanteladas más tradicionales (también nucleares). En cuanto a la energía eólica, se muestra la evolución a lo largo de los últimos años en Europa en la figura 1.3, tanto la instalada cada año como la acumulada. Se observa una clara tendencia creciente, además de un aumento de las instalaciones offshore, es decir, instalaciones en zona marítima. Por último, se va a observar las diferentes aportaciones de cada país. En la figura 1.4 se muestra la potencia instalada de cada uno de ellos (de la Unión Europea), siendo España el segundo país, por detrás de Alemania. Finalmente se observa la potencia instalada durante el año 2015 en la figura 1.5. Cabe mencionar la nula inversión de España pese a ser la segunda potencia europea. 2
Figura 1.2: Instaladas y decomisadas en Europa (2015) (a) Instalada anualmente 3
CAPÍTULO 1. INTRODUCCIÓN (b) Acumulada Figura 1.3: Energía eólica Figura 1.4: Potencia total instalada 4
2.1. VIENTOS general, en la dirección del viento dominante se dejan alrededor de 7my en la dirección perpendicular unos 4m. 2.1.1.4. Efecto túnel En este caso se aprovecha el entorno para un mayor rendimiento. Este se puede producir entre dos montañas por ejemplo, al pasar el aire por aquí se acelera. Aunque las formas de estas deben ser suaves, ya que de lo contrario se pueden producir excesivas turbulencias que anulen este efecto o incluso perjudiquen gravemente el rendimiento. 2.1.1.5. Efecto colina De igual forma se puede aprovechar una colina suave situando la torre en el punto más alto de esta. El efecto en el viento es el mismo, se acelera al disminuir la "sección"de paso. De igual forma si la colina es pronunciada la dirección del viento puede variar demasiado y perjudicar el rendimiento. 2.1.2. Rosa de los vientos Se conoce como ’Rosa de los vientos’ a una herramienta que facilita datos realmente significativos de una ubicación respecto al comportamiento del viento de forma muy compacta y visual. La figura 2.3 es un ejemplo correspondiente a Caen (Francia). Figura 2.3: Ejemplo de rosa de los vientos [1] Todos los valores dados por esta herramienta son valores relativos en porcentaje. La circunferencia se divide en sectores (por ejemplo de 30◦). Cada sector estará conformado por tres porciones. El radio de la más exterior indica el porcentaje de tiempo que el viento sopla en esa dirección. La inmediatamente posterior indica el producto de la anterior por la velocidad media del viento en esa dirección, normalizada con el total. Por último, la roja y más interior indica el producto de la primera magnitud por el cubo de la velocidad media del viento en esa dirección (también normalizado). Está 11
CAPÍTULO 2. ENERGÍA EÓLICA última es quizás la más importante como se verá en el capítulo 3, pues muestra la energía disponible del viento. El periodo de tiempo puede variar, pudiéndose referirse a un día, un mes o un año. Por tanto, en el ejemplo mostrado vemos como el viento predominante, con mayor frecuencia y más energía, se encuentra contenido en un rango de 90◦. 2.2. Turbinas eólicas Antes de entrar a describir las partes de una turbina eólica moderna, se analizarán los diferentes tipos de turbina ([6, 9]). Hay dos grupos principalmente, según el eje de giro sea vertical u horizontal. 2.2.1. Turbinas de eje horizontal (HAWT) Las figuras 2.4 muestran dos tipos de turbinas de eje horizontal. Sus principales ventajas son la capacidad para alcanzar alturas mayores, con lo que se pueden aprovechar velocidades mayores de viento, y el mayor rendimiento respecto a una de eje vertical debido a que el eje de rotación es colineal a la dirección del viento. De otra manera en ciertas fases del movimiento la fuerza del viento ejercería un efecto contraproducente, como se verá más adelante. (a) 3 palas (b) 2 palas Figura 2.4: Turbinas de eje horizontal Dentro de las desventajas se incluye la necesidad de una torre y una cimentación que soporte grandes esfuerzos, ya que todo el mecanismo de engranaje debe ir junto a las palas, es decir, en la parte alta de la torre. Además, el rotor debe situarse a ser posible de manera que sea lo primero que se encuentre el viento, ya que aguas abajo se genera una turbulencia muy negativa. Por otro lado habrá que incluir un mecanismo que varíe la posición relativa a la dirección del viento, y otro que bloqueé el rotor a una velocidad determinada para que no puedan surgir problemas estructurales. 12
2.2. TURBINAS EÓLICAS 2.2.2. Turbinas de eje vertical (VAWT) En este caso el eje principal de rotación es vertical, lo cual tiene ciertas ventajas. Todo lo referente a engranajes y generador puede estar situado bajo la torre. Esto supone que la torre solo debe soportar esfuerzos debido a su propio peso y a las palas en sí mismas. Además, no tienen que orientarse respecto a la dirección del viento, ideal para emplazamientos donde el viento es altamente fluctuante. Todo esto se traduce en un mantenimiento mucho más barato y sencillo. Normalmente este tipo de turbinas están situadas a alturas relativamente pequeñas (comparadas con las de eje horizontal), ya que es difícil situarlas en torres altas. Esto tiene un gran defecto: a alturas pequeñas la velocidad del viento es mucho menor. Sin embargo, esto se puede mitigar situándolas sobre edificios, teniendo una altura mayor (la suma de la del edificio y la torre) y aprovechando la mayor velocidad del viento debido a que tiene que superar el edificio y se deflecta hacia arriba. En general, la eficiencia de esto tipo de turbinas es menor que las de eje horizontal, principalmente porque siempre existe una resistencia de las palas en contra de la rotación, debido a que el eje de rotación es perpendicular al viento. (a) Daerrieus [6] (b) Savonius Figura 2.5: Turbinas de eje vertical Podemos distinguir dos tipos principalmente de turbinas de eje vertical, según sea la rotación generada por una fuerza sustentadora o de resistencia. La turbina Daerrieus rota por una fuerza sustentadora. Requiere una fuente externa de potencia que inicie la rotación, ya que el torque necesario para ponerlo a rotar es bastante alto. Además, sufre bastante esfuerzos de fatiga. El principio por el que gira se puede ver en la figura 2.6a. La turbina Savonius rota debido a la resistencia del viento. Su geometría hace que la resistencia en una zona sea mayor que en otra, generando así un par de rotación. Estas turbinas son eficientes en zonas donde el viento varía su dirección con asiduidad, teniendo además la capacidad de rotar a bajas velocidades de viento. Una esquema se puede ver en la figura 2.6b. 13
CAPÍTULO 2. ENERGÍA EÓLICA (a) Principio de rotación de la turbina Daerrieus (b) Principio de rotación de la turbina Savonius Figura 2.6: Turbinas de eje vertical 2.2.3. Componentes de una turbina A partir de ahora se hablará exclusivamente de turbina de eje horizontal. Como cualquier turbina eólica el objetivo de este dispositivo es generar energía eléctrica a partir del viento. Para ello son necesarios varios componentes principales (mirar figura 2.7). Palas: son las encargadas de transmitir el movimiento del viento a un eje de giro. Está formada por perfiles aerodinámicos que generaran sustentación y resistencia. En este caso, será la sustentación la responsable del giro de las palas. El número de palas vendrá determinado por la velocidad de giro principalmente. Si la velocidad es alta el número de palas será en general menor y viceversa. El objetivo es siempre aumentar el rendimiento. Deben tener un gran coeficiente resistencia/peso que garantice que soporta los elevados esfuerzos. En caso de tener dos palas los efectos vibratorios son muy importantes. Torre: es la estructura encargada de soportar el peso del resto de componentes, así como resistir los esfuerzos generados por el viento y la rotación de las palas (fatiga). Tanto su capacidad para soportar tales esfuerzos como el coste económico que tiene su instalación son factores críticos a la hora de estudiar la viabilidad de una turbina de cierto tamaño. Góndola: es la cavidad donde se alojan físicamente todos los componentes que deben estar en la parte superior de la torre. Rotor: es el conjunto de las palas y la raíz conectados al eje principal de giro. Raíz: es la pieza donde van unidas las palas, conectada al eje. Tiene una forma cuya aerodinámica es lo menos perjudicial posible. Caja de engranajes: es el mecanismo por el cual la velocidad del eje principal se transforma a través de engranajes a una velocidad de rotación mayor, en vistas a la generación de energía eléctrica. 14
2.2. TURBINAS EÓLICAS Figura 2.7: Partes de una turbina eólica Generador: es el dispositivo electrónico capaz de convertir una energía mecánica (rotación) en energía eléctrica. Transformador: es el dispositivo que varía según convenga el voltaje de la energía obtenida. Anemómetro: es el dispositivo que mide la velocidad del viento. Freno: por razones de seguridad debe haber un mecanismo de frenar la turbina si la velocidad del viento es excesiva. Se conoce como "velocidad cut-inçomo aquella velocidad del viento a partir de la cual la turbina puede generar energía eléctrica y "velocidad cut-outçomo aquella a partir de la cual la turbina podría sufrir daños estructurales si esta está rotando. Pitch control: controla el ángulo de incidencia del viento sobre las palas. Puede usarse como freno de seguridad disminuyendo la incidencia disminuyendo así la velocidad de rotación. Mecanismo de orientación: sirve para orientar el eje de la turbina según la dirección del viento. En turbinas de baja velocidad pequeñas puede consistir en una superficie que genera resistencia de tal manera que el eje y el viento se alineen. 15
CAPÍTULO 2. ENERGÍA EÓLICA 16
Capítulo 3 Análisis matemático En este capítulo se van a describir las ecuaciones que modelan el movimiento de la turbina y la atmósfera, así como las hipótesis aplicadas para su simulación ([6, 8, 4, 9]). Para la simulación solo se considerarán las palas y la raíz. En ningún momento se tendrá en cuenta otro componente, como por ejemplo la torre o la góndola. De cara a la obtención de resultados habría que tenerlas en cuenta, pero para comparar diferentes geometrías de la palas parece razonable que, la mejor geometría sin tenerlas en cuentas será la mejor cuando si se haga. Por otro lado, solo se pretende buscar el estado estacionario, por lo que las constantes de tiempo o cualquier efecto que transcurra hasta llegar a tal estado no es relevante. No obstante, se explicarán algunos métodos aplicables en la simulación para trabajar con inercias aproximadas (basadas en hipótesis) que podrían utilizarse para calcular el estado transitorio. A continuación se analizará que datos de entrada son necesarios para modelar la turbina (además de la geometría). 3.1. Ecuación del movimiento El modelo será el de un solido rígido con un eje fijo. Esto es equivalente a un sólido con dos puntos fijos (P1yP2), siendo la dirección que los une (~ν) el eje de rotación. Es decir, en las ecuaciones de conservación de la cantidad de movimiento y el momento cinético se pueden sustituir tal vínculo por dos reacciones vinculares (φ1yφ2), una en cada uno de estos dos puntos. Por tanto las ecuaciones quedan como en 3.2. ~ Ces la cantidad de movimiento, ~ Fla resultante de las fuerzas exteriores y ~ MP1 el momento resultante de tales fuerzas en el punto P1. El punto P1podría ser el centro de la turbina, por comodidad. ∂~ C ∂t =~ F+~ φ1+~ φ2(3.1) ∂~ ΓP1 ∂t =~ MP1+−−−→ P1P2×~ φ2(3.2) Proyectando ambas expresiones (son 6 ecuaciones) sobre la dirección ~ν solo queda una ecuación que no depende de las reacciones vinculares (ecuación 3.3). En la simulación el eje de rotación tendrá la misma dirección que el eje zdel sistema de referencia inercial. Por lo tanto, por comodidad usaremos que ~ν =~ k. Además, −−−→ P1P2=α~ ky~w =wz~ν. ∂~ ΓP1 ∂t ·~ k=~ MP1·~ k(3.3) 17
CAPÍTULO 3. ANÁLISIS MATEMÁTICO ~ ΓP1=I ~ ~ ·~w = IxIxy Ixz Iyx IyIyz Izx Izy Iz · 0 0 wz = Ixzwz Iyzwz Izwz (3.4) De la ecuación 3.3 queda lo siguiente. Mz=Iz ∂wz ∂t (3.5) Como observación, se aproximará el tensor de inercia de una manera muy simplificada. Dada la geometría de las palas y la raíz, se podría decir que tiene dos planos de simetría. Esto no es así pero dadas las dimensiones se podría admitir de manera aproximada. Un plano sería el plano donde ’están contenidas’ las tres palas (o cuantas haya) y el otro sería el plano perpendicular a este y que corta una de las palas en dos partes. Supuesto esto, se puede demostrar que las direcciones normales a estos planos definen ejes principales de inercia que pasan por el centro de masas. Además, también se puede afirmar que todo eje de simetría es eje principal de inercia, por tanto la intersección de ambos planos también es eje principal de inercia. Por tanto si hay tres planos perpendiculares cuyas direcciones normales son principales de inercia, el tensor está en direcciones principales de inercia (Ixy =Ixz =Iyz = 0).Y además: Iz=Ix+Iy(3.6) 3.2. Variables atmosféricas En este capítulo se describirá el modelo de atmósfera utilizado. Se considerarán las variables (temperatura, presión, densidad y viscosidad) constantes, así como un perfil de velocidades del viento uniforme. Por tanto las simulaciones se han establecido de tal manera que variando el valor de la velocidad del viento y la altura queda definidas el resto de variables. Esta simplificación puede justificarse por el diámetro de la turbina que se utilizará. Si este fuera mayor y/o se quisiera llegar a resultados más exactos sería igualmente fácil de modelar, aunque computacionalmente más caro (en términos de tiempo de computación y memoria requerida). Se tomará la dirección del viento como la del eje de revolución de la turbina, por tanto con el módulo de la velocidad esta queda definida por completo. Es decir, se compararán diferentes diseños de turbinas en el caso óptimo, en cuanto a que el viento será uniforme y constante, aunque esto no se ajuste a la realidad. Para las variables intensivas el modelo de atmósfera utilizado es el de atmósfera ISA (International Standard Atmosphere) definido por las ecuaciones 3.7. T=T0+λ(h−h0) ρ=ρ0 T T0!g Rλ−1 P=ρRT µ=1.458 ×10−6√T 1 + 110,4 T ν=µ ρ(3.7) 18
3.3. LÍMITE DE BETZ Siendo T[K]la temperatura, ρ[kg/m3],P[Pa]la presión, µ[Pa s]la viscosidad cinemática y ν[m2/s] la viscosidad cinemática. Las constantes son las siguientes: λ= 6.5×10−3[K/m] g= 9,80665[m/s2] R= 287,058[J/(KgK)] T0= 288,15[K] ρ0= 1,225[kg/m3] 3.3. Límite de Betz En esta sección se demuestra cual es el límite de energía que se puede obtener del viento. Esta teoría, denominada teoría del disco actuador hace las siguientes hipótesis: La turbina se reduce a un disco cuyo eje es el mismo que el de esta. El movimiento del aire es unidireccional. La velocidad del aire es continua a través del disco. La presión sufre una discontinuidad en el disco. En la figura 3.1 se puede apreciar la distribución de presiones y velocidades a lo largo de la dirección del eje de la turbina. El primer paso es obtener una relación entre las velocidades aguas abajo, aguas arriba y la que existe en el plano de la turbina. c=c1=c2(3.8) El gasto es: ˙m=ρAc (3.9) Fx= ˙m(cu−cd)(3.10) A continuación se hallarán dos expresiones de la fuerza axial que se relacionarán para llegar a una expresión que relacione las velocidades. Para ello se usará que las presiones de remanso deben conservarse en la región anterior al disco y en la posterior. pa+1 2ρc2 u=p1+1 2ρc2(3.11) pa+1 2ρc2 d=p2+1 2ρc2(3.12) p1−p2=1 2ρ(c2 u−c2 d)(3.13) Fx=A(p1−p2) = 1 2ρA(c2 u−c2 d)(3.14) Igualando 3.10 y 3.14 se obtiene 3.15 c=1 2(cu+cd)(3.15) Ahora hay que buscar una expresión de la potencia obtenida, haciendo uso de las entalpías de remanso. 4h0=h0u−h0D=hu+1 2c2 u−hD+1 2c2 d(3.16) 19
CAPÍTULO 3. ANÁLISIS MATEMÁTICO Figura 3.1: Teoría del disco actuador 20
Capítulo 4 Star-CCM+: Preparación de la simulación El software utilizado en el proyecto es Star-CCM+, de CD-Adapco. Este software se base en lo que se conoce como CFD (Computer Fluid Dynamics, es decir, en la resolución de problemas que incluyen fluidos haciendo uso de diferentes algoritmos. Más adelante se explicarán diferentes maneras de abordar tales algoritmos en función del tipo de problema, ya que la tipología de este hará que se puedan hacer unas u otras simplificaciones. Solo se mencionarán algunos tipos de manera cualitativa, no se entrará al detalle en las ecuaciones de fondo mostradas anteriormente, ya que son complejas y no es el objetivo del proyecto. El objetivo de la simulación además de la obtención de resultados es automatizarla todo lo posible para poder cambiar la geometría fácilmente. Para ello el software tiene la opción de grabar macros, lo que facilita esta tarea sin la necesidad de saber programar en java. En este capítulo se explicarán los pasos a seguir para llegar a simular una turbina eólica. La mayoría de opciones escogidas se indicarán mediante capturas de pantalla, mientras que las más básicas, que se manejan directamente en el árbol se señalarán mediante guiones. Por ejemplo, para crear un modelo CAD en la geometría habría que expandir el nodo Geometry, seleccionar 3D-CAD Models y hacer click derecho, y escoger nuevo. Geometry −→ 3D-CAD Models =⇒New 4.1. Esquema general En esta sección se van a presentar los diferentes pasos a seguir para llegar a resolver el problema. Tales pasos coinciden para los diferentes métodos basados en la teoría de elementos finitos o similares. El primer paso es la planificación de las regiones posteriores. Cada región tendrá fronteras en las que se aplicarán condiciones de contorno, y que además harán de interfaces con posibles regiones contiguas. En este caso habrá dos regiones, una estática y otra interior que girará y que tendrá como una de sus fronteras la turbina. A continuación hay que construir la geometría. En este paso se parametrizará todo lo posible la turbina de modo que sea fácil cambiar su geometría. La turbina se definirá mediante diez secciones a lo largo de cada pala, introduciendo un perfil aerodinámico con un cierto ángulo de ataque en cada una de ellas. En cada región hay que definir las diferentes superficies que las delimitan para poder luego generar las condiciones de contorno como ya se ha comentando. Lo siguiente es la generación de las regiones a partir de la geometría, imponiendo las condiciones de contorno oportunas en cada una de ellas y generando interfaces cuando sea necesario. 27
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Ahora se construye el mallado. Tradicionalmente se genera la malla para cada región, aunque este software también permite construirla directamente a partir de la geometría. Ya que es más sencillo hacerlo de esta manera y que parece ser la tendencia actualmente se hará así. Para continuar, hay que generar los modelos físicos necesarios. En este caso solo será necesario definir uno al que se llamará ’Aire’. Más adelante se describirá el modelo escogido. Las dos regiones tendrán este modelo. Por último queda el procesado que incluye tres pasos. Primero (pre-processing) hay que preparar los parámetros propios de la simulación y definir la tipología de esta. Este paso incluye preparar todo lo necesario para guardar los resultados que sean necesarios. Después (processing), basta con ejecutar la simulación y esperar a que la simulación termine según los parámetros establecidos en el paso anterior, o simplemente pararla manualmente pasado un tiempo. Esta opción es útil en caso de que pueda aparecer algo no previsto, ya que mientras la simulación se ejecuta se mostrarán los resultados. Por último, una vez se ha detenido, se interpretarán los resultados obtenidos (post-processing). Este paso es probablemente el más importante, ya que habrá que observar si los resultados concuerdan con nuestro problema, ya que, además de existir la posibilidad de haber introducido algún error en el modelado del problema, existe también la posibilidad de que por problemas numéricos (por la malla, por los criterios de convergencia, etc) no sea una buena solución. El objetivo de la simulación es la obtención de una curva de rendimiento en función de la velocidad de giro de la turbina con un viento y una geometría dados. Dentro de la interpretación de los resultados cabe preguntarse si estos son razonables o no, así como comprobarlos. ¿Cual es la forma adecuada de hacer esto? Supóngase que solo se quiere calcular el rendimiento para una velocidad de giro. Se pondría como velocidad inicial de giro la deseada, y se mantendría esta constante hasta que llegará a un rendimiento estacionario. Una vez llegado ese punto, se refinaría la malla en función de algún criterio, y se volvería a repetir lo anterior partiendo de la solución encontrada anteriormente. Si finalmente vuelve a la misma solución, esta sería probablemente adecuada. Este proceso se puede repetir tantas veces como se quiera, teniendo en cuenta la capacidad de computación y el tiempo requerido para ella. Dado que hacer esto para toda la curva de rendimiento es muy costoso, se ha optado por comprobar la solución solo en algún punto de la curva. Se ha supuesto que si la malla es lo suficientemente buena para que la solución converga adecuadamente en ese punto, lo es para toda la curva. En realidad, al menos cuando el fluido empieza a desprenderse (con lo que el rendimiento empieza a caer) se debería de volver a refinar. 28
4.2. GEOMETRÍA 4.2. Geometría En esta sección se va a explicar como se ha definido la geometría. El problema se va a modelar a través de dos volúmenes de control. El primero que será la geometría más exterior tendrá forma de ’bala’. El segundo será un cilindro vaciado con la geometría de la turbina que girará dentro del primero. Cada una de estas dos geometrías se llama part. 4.2.1. Bala Para escoger las dimensiones de este part se ha usado la referencia dada en el manual de StarCCM+. En el se dice que la superficie de entrada del flujo debe estar a 10-20 veces la dimensión característica del objeto (la turbina) y la superficie de salida a 20-40 veces. Para la creación de esta part se siguen los pasos que se muestran a continuación. Se omitirán los detalles del sketch ya que con la imagen queda claro lo que se ha realizado. De igual forma no se detallarán las opciones de Revolve debido a que posteriormente se volverá a utilizar esta acción para crear la turbina y ahí se profundizará más. Geometry −→ 3D-CAD Models =⇒New YZ =⇒Create Sketch Sketch 1 =⇒Revolve 4.2.2. Turbina A partir de ahora se va a generar la geometría de la turbina, por lo que para parametrizarla todo lo posible se va a grabar un macro con todos los pasos. Para ello se pulsará sobre Start Recording..., quedando en un documento de texto todas las líneas de comando necesarias para crear la geometría. 1. Geometry −→ 3D-CAD Models =⇒New 2. ZX =⇒Create Transform Sketch Plane Este paso se repetirá 15 veces, estableciendo 15 parámetros de distancia de cada perfil al centro, y 15 parámetros que indican los ángulos de ataque de cada perfil. Los primeros se nombrarán como ’L1’,’L2’,etc y los segundos como ’aoa1’,’ao2’,etc. Los parámetros que indican la distancia se elegirán de manera lineal (con un metro de distancia entre cada uno), aunque si se quisiera refinar la pala en una determinada zona se podrían poner más perfiles cerca a esa zona y aumentar la distancia entre perfiles lejanos a esa zona. Otra posibilidad sería poner más o menos planos 29
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Figura 4.1: Paso 2 Figura 4.2: Paso 4 para refinar o hacer una geometría más burda. Aparte de estos 15 planos se creará uno a 0,5m del origen para poder modelar la zona de la raíz mejor, siendo esta distancia parametrizada como ’Interseccion_raiz’. 3. TransformSketchPlane =⇒Export Coordinate System Esto se hace también para todos los planos creados anteriormente. 4. 3D-CAD Model 1 −→ Import −→ 3D Curve =⇒perfil.csv Se importa la curva una vez para cada sistema de coordenadas creado. 5. 3DCurve 1 =⇒Extrude Igualmente se hace para todas las curvas. Para que la longitud de la pala sea exactamente la que se pretende, hay que ser cuidadoso y extruir la última curva en el sentido opuesto. 6. Body 1 −→ Transform =⇒Scale De nuevo se hace esta operación para todos los cuerpos que se generaron en el punto anterior. El escalado de cada perfil hay que hacerlo en sus respectivos ejes de coordenadas. Dado que el perfil que se importó poseía una cuerda de 1m, hacer un escalado con un valor de 1 hará que la cuerda sea de 1m, y hacer un escalado de xhará que la cuerda sea de xmetros. En todas las simulaciones se ha usado la misma distribución de escalados a lo largo de la pala. Hacer un estudio detallado sobre la influencia de esta distribución en el rendimiento sería fundamental, pero dados los recursos se ha preferido comparar diferentes perfiles. La forma más cómoda de cambiar estos valores es accediendo al archivo de texto generado con el macro y cambiarlo ahí, puesto que no se puede parametrizar como se hizo con los planos distribuidos a lo largo de la pala. La distribución utilizada es la siguiente, siendo el plano 1 el más cercano a la raíz (que se 30
4.2. GEOMETRÍA Figura 4.3: Paso 5 31
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Figura 4.4: Paso 6 encuentra a 0,5mdel eje, y el plano 16 el plano donde se define la punta de la pala. Plano 1 −→ 0.4 Plano 2 −→ 1 Plano 3 −→ 1.1 Plano 4 −→ 1.2 Plano 5 −→ 1.3 Plano 6 −→ 1.4 Plano 7 −→ 1.3 Plano 8 −→ 1.2 Plano 9 −→ 1.1 Plano 10 −→ 1 Plano 11 −→ 0.9 Plano 12 −→ 0.8 Plano 13 −→ 0.7 Plano 14 −→ 0.6 Plano 15 −→ 0.5 Plano 16 −→ 0.4 7. Create Sketch From Face Edges Se selecciona una de las superficies de cada cuerpo que fueron extruidos y se aplica el comando, haciéndolo siempre sobre la superficie más cercana a la raíz. 8. YZ =⇒Create Sketch La raíz se ha modelado como un cuarto de elipse que posteriormente se revolucionará. 9. Sketch 1 =⇒Create Revolve 10. Body 17 −→ Transform =⇒Translate Para poder intersectar los perfiles con la raíz tal y como se ha diseñado esta es necesario trasladarla tal y como se indica en la figura. 11. FaceSketch (De 2 a 5) =⇒Create Loft Se seleccionan los FaceSketch de 2 a 5 en ese orden para crear una primera parte de la pala. Se podría hacer toda la pala con una sola operación pero es más proclive a que de problemas a la hora de aplicarla. En esta se ha seleccionado de tal forma que la geometría sea normal al primer 32
4.2. GEOMETRÍA Figura 4.5: Paso 7 Figura 4.6: Paso 8 33
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Figura 4.7: Paso 9 Figura 4.8: Paso 10 34
4.2. GEOMETRÍA Figura 4.9: Paso 11 plano por el cual se ha definido. Además se ha seleccionado la opción Match Vertices Manually para que el borde de salida tenga la forma deseada y no aparezcan curvas extrañas en este. Con marcar esta opción se seleccionan automáticamente los puntos que queremos, por lo que no hay que seleccionar nada más. 12. FaceSketch (De 5 a 16) =⇒Create Loft De igual forma se procede con el resto de perfiles hasta la punta de la pala. 13. FaceSketch (De 1 a 2) =⇒Create Loft En este caso se seleccionan las casillas de manera que tanto en el perfil inicial como en el final la geometría sea normal. El resto queda igual que en los dos casos anteriores. 14. FaceSketch 1 =⇒Extrude Para terminar de hacer una geometría que posteriormente se pueda unir en una sola es necesario extruir la pala hasta la raíz tal y como se puede apreciar en la figura. Hay que hacer click sobre el cuerpo que representa la raíz una vez seleccionada la opción Up to body. 15. Body (De 18 a 21) −→ Boolean =⇒Unite Con esta opción solo queda un cuerpo uniendo todas las partes. Además se le cambia el nombre al cuerpo resultante a ’Turbina’. 16. Features −→ Solid Primitive =⇒Cylinder Se crea un cilindro que ayudará a refinar la malla en el borde de ataque. Además se parametriza el radio para pode cambiarlo si se quiere considerar que el borde de ataque puede incluir una zona mayor. 17. ’Turbina’ y ’Body 22’ (cilindro) =⇒Slice Con esta operación se dividirá la geometría ’Turbina’ en dos: ’Turbina’ y ’Borde de ataque’. 18. ’Turbina’ y ’Borde de ataque’ −→ Transform =⇒Rotate Este paso habrá que repetirlo tantas veces como palas tenga la turbina. Se podría hacer de una sola vez usando la opción Circular Pattern si solo se tuviera una geometría a la que aplicarlo. Pero en ese caso, habría que generar luego más cilindros y aplicar Slice para cada pala, por tanto 35
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Figura 4.10: Paso 12 36
4.3. MALLA características y recomendaciones a tener en cuenta son las siguientes: La media de caras de cada elemento es de 14. El número máximo de elementos por procesador es de 120 millones. Por superficie requiere 5 veces menos elementos que una malla basada en tetraedros. Dado el mismo número de elementos, requiere el doble de memoria que una malla trimmed. (ii)Trimmed Mesh Esta malla se basa en hexaedros principalmente. Es una malla robusta y eficaz, capaz de modelar problemas sencillos o complejos. Para un procesador no usar más de 100.000 celdas. Usar una simulación en serie si hay menos de 100.000 celdas, pues la realimentación entre procesadores hace que la simulación en paralelo no sea óptima para esos casos. Para una simulación de tipo segregated (se explicará después), requiere 0.5GB por cada millón de celdas. 4.3.2. Opciones de mallado Como ya se ha mencionado la malla que se usara es trimmed, por lo que las opciones serán las que existen para este tipo de malla. Por lo general serán las mismas aunque en algunos aspectos puede variar algo. Algunos cambios significativos se mencionarán. Para crear la malla basta con hacer lo siguiente una vez para cada región: Geometry −→ Operations −→ Automated Mesh −→ Trimmed Cell Mesher,Surface Remesher yAutomatic Surface Repair Se nombran respectivamente: ’Malla:Rotor’ y ’Malla:Bala menos rotor’. Ahora se van a explicar las opciones escogidas en la primera, siendo extrapolables al segundo caso, pero cambiando valores. Estos se mostrarán más adelante. Meshers Se pueden mallar las partes por separado si fuera necesario, al igual que mallar las partes en paralelo o concurrentemente. Dado que solo hay dos partes, y que una de ellas tiene una malla bastante basta, se hará en serie, ya que la mejora de tiempo no es apreciable. La diferencia entre hacerlo en paralelo o de manera concurrente es que en el primer caso se divide cada part y se usan los diferentes procesadores para mallar todo, y en el segundo cada procesador malla una part. Meshers −→ Surface Remesher Esta opción permite preparar las superficies para posteriormente generar el mallado. Se dejan todas las opciones por defecto en principio. Cuando se quiere remallar en las superficies habrá que escoger una tabla en Field Function Refinement Table. 43
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Meshers −→ Automatic Surface Repair Con esta opción se reparan posibles defectos de la geometría que pueda dan problemas a la hora de efectuar el mallado. Lo más adecuado es hacerlo manualmente (requiere bastante experiencia) con una herramienta disponible en la geometría para evitar que se lleven a cabo acciones indeseadas en la geometría en el proceso. Sin embargo, en este caso se ha comprobado que esto no sucede y por lo tanto se ha usado finalmente esta opción. Meshers −→ Trimmed Cell Mesher Esta opción es la que controla el mallado del volumen. Se deja todo por defecto excepto cuando se quiera refinar la malla en el volumen. Para ello, igual que antes, se escoge una tabla en Field Function Refinement Table. Default Controls −→ Base Size Es la longitud de referencia sobre la que se indicarán otras magnitudes. Default Controls −→ Target Size Se puede dar un valor relativo a Base size o uno absoluto. En este caso daremos un valor relativo. Indica el tamaño que de una celda estándar, que no tiene ninguna limitación. Default Controls −→ Minimum Size También se dará un valor relativo e indica el mínimo tamaño de una celda. En caso de querer refinar, es posible que haya que disminuir su valor ya que puede ser limitante. Default Controls −→ Surface Growth Rate Indica la velocidad de crecimiento entre celdas contiguas que se encuentran en una superficie, concretamente es el cociente entre el tamaño característico de los lados de ambas. Default Controls −→ Volume Growth Rate Indica la velocidad de crecimiento de las celdas en un volumen. Si se elige la opción custom, un número menor hará que crezca más rápido y viceversa. Custom Controls =⇒New −→ Surface Control 44
4.3. MALLA Esta acción se realiza tres veces, una para los bordes de ataque, otra para la turbina y otra para la intersección del rotor con la bala. Se nombran por tanto ’Borde de ataque’, ’Borde de salida’ e ’Interseccion con bala’. Al crear ’Malla:Bala menos rotor’ se hará lo mismo con la intersección, escogiendo los mismos valores (absolutos) para esta. De esta forma las dos mallas tienen un tamaño similar en la intersección. Una vez creados estos controles para refinar la malla en ciertos lugares, se pueden escoger los diferentes valores para las magnitudes descritas anteriormente (Base Size, Target Size,etc). Se pueden escoger iguales a estos (Use Parent Value) o escogiendo valores nuevos. Además se puede concretar un nuevo parámetro: Custom Controls −→ Borde de ataque −→ Values −→ Trimmer Surface Growth Rate Esta magnitud indica el número de celdas que siguen los valores del control de superficie. Si se especifica un valor de 100 habrá más celdas del tamaño especificado para la superficie que con un valor de 10. (a) Trimmer Surface Growth Control= 10 (b) Trimmer Surface Growth Control= 100 A modo de ejemplo se muestra la imagen anterior, en la que se aprecia una sección longitudinal de un cilindro, y cuya base (a la izquierda) es una superficie de control. Cabe señalar dos parámetros que aparecen si la malla es de tipo polyhedral que sirven para modificar el tamaño y el crecimiento de la malla: Default Controls −→ Mesh Density −→ Density El valor por defecto de este parámetro es 1. Si se dobla el número de elementos se doblará aproximadamente y si vale 0.5 será la mitad. Default Controls −→ Mesh Density −→ Growth Factor Sigue la misma idea que el anterior pero para modelar el crecimiento desde las superficies. El comportamiento de ambos se muestran a continuación: (a) Density= 10,Growth Factor= 1 (b) Density= 0,1,Growth Factor= 1 (c) Growth Factor= 10,Density= 1 (d) Growth Factor= 0,1,Density= 1 45
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN 4.4. Continua En esta sección se puede modelar la malla como una característica que se asignará a las regiones que se han creado posteriormente. Sin embargo, como ya se explicó la tendencia es hacerlo en la geometría por lo que se ha obviado explicar este procedimiento. Por lo general, el procedimiento es similar al ya expuesto para la geometría. Por otro lado, también se modelan las características de los modelos físicos y de resolución. Continua =⇒New −→ Physics Continuum Se renombra como ’air’. Continua −→ Aire −→ Models =⇒Select Models... Como se puede apreciar en la imagen hay diferentes opciones en cada modelo, y a medida que se van eligiendo van apareciendo nuevos. Además, por defecto se autoseleccionan algunos durante este proceso, aunque estos pueden revocarse. En los modelos asociados al espacio (Space), se selecciona Three Dimensional. En los modelos de tiempo (Time) se selecciona Implicit Unsteady. Lo que se pretende en la simulación es modelar el movimiento real de la turbina, por lo que no puede ser un modelo estacionario. Esto implica que la simulación utilizará un paso de tiempo. Además, se ha escogido el modelo implícito porque el explícito solo funciona para fluidos en régimen viscoso y laminar (tal y como menciona la ayuda del software), no pudiendo modelar régimen turbulento. Cabe destacar también que surge una incompatibilidad con el explícito debido a la existencia de un cuerpo en movimiento. Para profundizar algo más sobre el método explícito se podría destacar como, en una gran parte 46
4.4. CONTINUA de los problemas, es una cuestión de escalas de tiempo: si el problema tiene un paso de tiempo suficientemente pequeño se podría considerar laminar y teóricamente podría resolverse. Para juzgar que paso de tiempo sería adecuado se podría tener en cuenta el número de Courant, que se define como sigue: C=4t 4x u . La idea es que el cociente del paso de tiempo y el tiempo de residencia del fluido en una celda sea pequeño. Esto es, C≤1. La principal diferencia a la hora de resolver es que ,con el explícito, para resolver en un nodo en t+∂t se usa la solución conocida en ese nodo y en los nodos cercanos en t. Mientras tanto, en el implícito solo se usa la solución en se nodo en tjunto con la relación que debe guardar con los nodos cercanos en t+∂t. En caso de que esta relación entre nodos sea muy determinante en problema podría ser prácticamente irresoluble con el método explícito. Una vez escogidos estos, queda escoger el respectivo al material. Aquí se irán abriendo opciones a medida que se vayan escogiendo otras. Para empezar se selecciona gas. A continuación hay que escoger entre Segregated Flow oCoupled Flow. Para tomar una decisión se muestran algunas reseñas dadas en la ayuda: 1. El primero requiere menos memoria. 2. El segundo es más robusto y preciso para problemas con fluido compresible (imprescindible si hay ondas de choque). 3. El segundo es más robusto si existen una alta convección (alto número de Rayleigh). Se dan también las siguientes recomendaciones: 1. Utilizar Coupled Flow yCoupled Energy si hay flujo compresible, es un problema con convección importante o hay grandes energías. 2. Si no hay problemas en cuanto a recursos para la computación, utilizar Coupled Flow para fluidos incompresibles o isotérmicos. 3. Utilizar Segregated Flow para fluidos incompresibles o con baja compresibilidad. Teniendo en cuenta todo esto, se ha escogido Segregated Flow. El siguiente paso es escoger Constant Density. Esto es razonable ya que las velocidades son bajas en comparación a la del sonido, por tanto se puede aproximar como incompresible. Más adelante se comentará como se implementaría un cambio de la densidad con la altura que podría ser relevante si el diámetro de la turbina es lo suficientemente grande. A continuación se escoge Turbulent. En principio se habrá autoseleccionado Reynolds-Averaged Navier-Stokes, aunque esto se podría cambiar si fuera necesario. Se dejará por defecto esta opción. Para escoger a continuación el modelo adecuado entre las diferentes opciones (k−,k−ω, Spalart-Allmaras o RSEM), se van a mostrar algunos consejos dados en la ayuda de Star-CCM+. Spalart-Allmaras es adecuado para casos en los que la capa límite esta fuertemente adherida o pueda existir un desprendimiento muy leve. k−ofrece una solución de compromiso entre precisión, coste computacional y robustez. Es capaz de modelar problemas con recirculación o transmisión de calor. k−ωes similar a k−pero con la diferencia de que una de las dos variables en las ecuaciones difiere. Requiere de valores más pequeños de y+que el modelo k−y sería útil en casos en los que el número de Reynolds fuera pequeño (teniendo la capa límite mayor importancia) o en casos en los que la malla es bastante fina. Por lo general, está recomendado para aplicaciones similares que el modelo Spalart-Allmaras. RSEM está recomendado para aplicaciones en las que existe una turbulencia muy pronunciada, ya que la modela muy bien pero a un coste computacional muy alto. Por ejemplo, sería útil para 47
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN modelar un separador ciclónico. Por tanto, teniendo en cuenta todo esto se va a escoger el modelo k−. Por todo lo explicado en el capítulo 3 se mantiene la elección por defecto Realizable Two-Layer All-y+wall treatment. Además se autoseleccionan Turbulence Supression yTransition Boundary Distance. Con estos seleccionados el software define regiones en las que la turbulencia es despreciable, ahorrando tiempo de resolución. A continuación hay que seleccionar el material, en este caso saldrá por defecto. Si se quisiera cambiar los pasos serían los siguientes: Continua −→ Air −→ Models −→ Gas −→ Air =⇒Replace with... Ya solo queda escoger las propiedades. Las que no se mencionen se dejarán por defecto: Continua −→ Air −→ Models −→ Gas −→ Air −→ Material Properties −→ Density −→ Constant −→ $Densidad aire Continua −→ Air −→ Models −→ Gas −→ Air −→ Material Properties −→ Dynamic Viscosity −→ Constant −→ $Dynamic Viscosity Continua −→ Air −→ Reference Values −→ Reference Pressure −→ Value −→ $Presion aire La utilidad de este valor es reducir errores, ya que una diferencia de presiones baja puede tener grandes efectos. Usando el modelo de densidad constante no tiene ninguna relevancia. Continua −→ Air −→ Initial Conditions −→ Pressure −→ Constant −→ $Presion aire Continua −→ Air −→ Initial Conditions −→ Velocity −→ Constant −→ [0.0, 0.0, -$Velocidad ref] m/s 4.5. Regiones En esta sección se van a modelar finalmente las regiones del problema que fueron creadas a partir de la geometría. El primer paso es asignar a cada región su malla y el tipo de material tal y como indican las imágenes. Regions −→ Bala menos rotor −→ Propiedades Regions −→ Rotor −→ Propiedades 48
4.5. REGIONES (a) Bala menos rotor (b) Rotor Ahora toca modelar las superficies. Todas excepto las que siguen a continuación se quedarán con el tipo Wall que se les asigna por defecto. Regions −→ Bala menos rotor −→ Boudaries −→ Entrada −→ Propiedades −→ Type −→ Velocity Inlet En esta superficie se escogerá la velocidad del viento. En este caso se ha elegido constante e igual a la velocidad inicial que se escogió para la región, aunque podría ser variable usando una field function o una tabla (se explicarán más adelante). •Physics Conditions −→ Velocity Specifications −→ Components •Physics Values −→ Velocity −→ Method −→ Constant −→ [0.0, 0.0, -$Velocidad ref] m/s Regions −→ Bala menos rotor −→ Boudaries −→ Salida −→ Propiedades −→ Type −→ Pressure Outlet Esta superficie debe tener la presión del aire como condición de contorno del problema, en este caso constante. •Physics Values −→ Pressure −→ Method −→ Constant −→ $Presion aire$ Para completar adecuadamente esta sección se debe crear hacer uso de una herramienta que hay que crear previamente. Tools −→ Motions =⇒New −→ DFBI Embedded Rotation Habiendo hecho esto se creará un nodo en el árbol llamado DFBI (Dynamic Fluid Body Interaction). Esta herramienta permite calcular las fuerzas y momentos a los que se ve sometido un cuerpo en el seno de un fluido. Este nuevo nodo se modelará después. Solo creándolo ya se puede terminar el nodo Continua. Ahora hay que asignarle a la región ’Rotor’ el movimiento, mientras que ’Bala menos rotor’ debe estar fija. Regions −→ Bala menos rotor −→ Physics Value −→ Motion Specification −→ Motion −→ Stationary Regions −→ Rotor −→ Physics Value •Axis −→ Coordinate System −→ Laboratory •Axis −→ Origin −→ [0,0,0] 49
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN •Axis −→ Direction −→ [0.0, 0.0, 1.0] El eje se ha escogido de manera que la velocidad de la turbina sea negativa cuando esta este funcionando. Esto se ha hecho así dado que hay que escoger varios ejes durante la simulación y es más sencillo cogerlos todos positivos para que no puedan cometerse ciertos fallos. •Motion Specification −→ Motion −→ DFBI Embedded Rotation •Motion Specification −→ Reference Frame −→ Lab Reference Frame Por último queda crear interfaces entre las regiones. Seleccionar Boundaries de ambas regiones con el mismo nombre =⇒Create Interface Se creará un nodo denominado Interfaces que quedará por defecto, siendo estas de tipo interno y topología In-place. Hasta aquí se ha seguido el árbol de forma ordenada. Ahora se va a explicar de la forma más sencilla para llegar a la simulación que se busca. 4.6. Field Functions Se conoce como Field Function a cualquier expresión (escalar, vector o tensor) que pueda estar asociado a una región o superficie. Es una herramienta que se encuentra dentro del nodo Tools. Aunque también pueden usarse como una forma de hacer cálculos entre magnitudes que pueda servir para alguna representación de resultados. A continuación se van a escribir las definiciones de todas las funciones que han sido necesarias o podrían serlo en alguno de los casos considerados. Algunas definiciones no tendrán sentido hasta que se expliquen otros nodos. Cabe destacar que, por comodidad, para encontrar las funciones que en principio pueden ser modificadas se han nombrado como sigue: 0_Nombre. También se debe mencionar que se pueden poner unidades a cada una, aunque por defecto se trabaja con unidades del Sistema Internacional, por lo que no se van a detallar estas. 0_Altura: 50 Cambiando este valor se obtienen las demás magnitudes atmosféricas. 0_Base Size: 8 Ya que no se pueden utilizar field functions en la geometría (incluido el mallado), tener un valor igual al que se ponga en el mallado puede ser útil a la hora de remallar. 0_Compensador_barrido: (${Time}==0) ? 0 : ((${Time}<${maxtiempoaceleradorMonitor}) ? -1e13 : -${6-DOFBodyMomentzReport}) 0_Compensador_discreto: (${Time}==0) ? 0 : -${6-DOFBodyMomentzReport} 0_L_turbina: 15 Por la misma razón que antes, este valor debe tener el mismo valor que la longitud de la turbina. Esto es, debe ser igual a la distancia que se tiene el último plano al hacer la geometría. Haciendo esto no es necesario entrar a cambiar ninguna expresión, evitando posibles errores. 0_Masa_turbina: 3000 50
4.6. FIELD FUNCTIONS Teniendo un valor para esta magnitud se puede calcular la inercia por si se quisiera calcular algún transitorio. 0_R_rotor: 18 Al igual que antes, debe valer lo mismo que lo que se haya puesto en la geometría. (trimmed) 0_Refinar malla gradiente velocidades: (${Gradiente velocidades adim}<=4) ? 0.1*${Base size} : ((${Gradiente velocidades adim}>6) ? 2*${Base size} : 0) Esta función es un ejemplo de como se podría refinar la malla en caso de ser trimmed. (polyhedral) 0_Refinar malla gradiente velocidades: (${Gradiente velocidades adim}<=4) ? 0.1*1.2*pow(${Volume},1/3) : ((${Gradiente velocidades adim}>6) ? 2*1.2*pow(${Volume},1/3) : 0) Esta función sería un ejemplo si la malla fuera de tipo poliédrico. 0_Refinar malla wally+: (${WallYplus}>100) ? 0.75*1.2*pow(${Volume},1/3) : 0 Esta función podría servir para refinar la malla poliédrica en las superficies en función de como sea y+. 0_Tiempo acelerador: (${Time}<=5) ? 0 : ((${varianzarendimientoReport}>=0.001) ? 0 : (${Time}+1*${Time step})) Esta función servirá valdrá cero mientras consideremos que la turbina está en transitorio, excepto en todos los primeros pasos de transitorio (en este caso 1 pasos de tiempo) que indicarán hasta cuando debe estar acelerándose en cada transitorio. 0_Time step: 1 Dado que la simulación es de tipo Implicit Unsteady este valor no es muy relevante siempre y cuando no sea demasiado alto. Dadas las velocidades del problema y el número de iteraciones que se tendrá en cada paso este valor es razonable. 0_Velocidad ref: 8 Densidad aire: 1.225*pow(${Temperatura}/288.15,9.80665/(287.058*6.5e-3)-1) Densidad media turbina: ${Masa turbina}/(3.14159*4*pow(${R_rotor},2)-${VolumenrotorReport}) Dynamic Viscosity: (1.458e-6*sqrt(${Temperatura}))/(1+110.4/${Temperatura}) Gradiente velocdiades: mag(grad(mag($${Velocity}))) Gradiente velocidades adim: log10(${MáximogradienteReport}/(abs(${Gradiente velocidades})+1e-8)) De esta forma se obtiene un número adimensional que compara órdenes de magnitud del gradiente de velocidades en cada punto con el del máximo gradiente. Si el número resultante está próximo a cero, quiere decir que en ese punto el gradientes está cerca del máximo. Si en cambio vale 2, por ejemplo, el gradiente en ese punto vale unas 100 veces menos que el máximo. El factor de 51
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN 10−8está para asegurar que el denominador no se acerca demasiado a cero, ya que esto daría error. Izcilindro completo aire: 1/12*${Densidad aire}*3.14159*pow(${R_rotor},2)*(3*pow(${R_rotor},2)+pow(4,2)) Izturbina: ${Densidad media turbina}*${Densidad aire}*(${I_z cilindro completo aire}- ${I_zrotor(aire)Report}) Izturbina para todo tiempo: [5e13,5e13,1e14] Hay que tener en cuenta cual es el paso de tiempo, ya que junto con este vector determinará cuanto se acelerará la turbina entre estacionario y estacionario. Si se quiere aproximar la inercia en lugar de poner una aleatoria, se debería calcular en el primer paso de tiempo, es decir: (${Time}==0) ? [5e13,5e13,1e14] : [${I_z turbina}/2,${I_z turbina}/2,${I_z turbina}] Integrando Izrotor (aire): ${Densidad aire}*pow($${Position}[2],2) Kynematic Viscosity: ${dynamic viscosity}/${Densidad aire} Limite de Betz: 1/2*16/27*pow(${velocidad ref},3)*3.14159*pow(${L_turbina},2)*${Densidad aire} Presion aire: ${Densidad aire}*${Temperatura}*287.058 Rendimiento: ${6-DOFBodyAngularVelocity1Report}*${6-DOFBodyMomentzReport}/${Límite de Betz}*100 Reynolds 9: ${Densidad aire}*$${Position}("3D-CAD Model 1 9")[0]*$${Velocity}/${dynamic viscosity}/(1*2) Reynolds 13: ${Densidad aire}*$${Position}("3D-CAD Model 1 13")[0]*$${Velocity}/${dynamic viscosity}/(0.6*2) Reynolds 14: ${Densidad aire}*$${Position}("3D-CAD Model 1 14")[0]*$${Velocity}/${dynamic viscosity}/(0.5*2) En las tres anteriores se ha tenido en cuenta cuanto mide la cuerda en cada una de las secciones respectivamente. Temperatura: 288.15-6.5e-3*(${Altura}) Uno: 1 Velocidad viento (y): -${Velocidad ref}*pow(($${Position}[1]+${Altura})/${Altura},1/7) Esta expresión es típica para tener en cuenta el gradiente de velocidades cercano al suelo, similar a la expresión 2.1. 52
4.9. REPORTS,MONITORS YPLOTS Frecuency: 15 Start: 0 Se ha elegido el número de iteraciones como la mitad de las totales para cada paso de tiempo. De esta forma se actualiza dos veces para cada paso. Por tanto serán 15 iteraciones si la simulación es un barrido, o 25 si solo se simula una velocidad de giro. Aquí es donde se elige la función que define el compensador según la simulación también. Por tanto ya se han visto todas las variaciones entre los dos tipos de simulación. Resumiendo, hay que cambiar la definición del compensador en este report, las iteraciones totales de cada paso y las correspondientes a este report y la velocidad de rotación inicial. •Report: Tiempo acelerador ◦Tipo: Volume average report ◦Scalar Field Function: 0_Tiempo acelerador ◦Parts: Rotor ◦Monitor: Trigger:Time Step Frecuency: 1 Start: 0 Usar un report de tipo Expression report da ciertos problemas, por lo que dado que se necesita un número se ha escogido la media de una field function cuyo valor es igual en todas las celdas. La part escogida podría ser cualquiera. Ahora hay que hacer algunas operaciones algo tediosas y aparentemente innecesarias, pero se ha tratado de simplificar todo lo posible, pero simplificándolo más no funciona como se espera. Primero es necesario crear un par de monitor que no estan asociados a ningún report, pero que ofrecen otras posibilidades. •Monitors =⇒New Monitor −→ Field Variance Se nombra ’Rendimiento acelerador’ y sus características son las siguientes: •Monitor: Rendimiento acelerador ◦Definición: Field variance monitor ◦Part: Rotor ◦Sample Count: 4 ◦Enable Sliding Sample Window:X ◦Field Function: Rendimiento ◦Trigger:Time Step ◦Frecuency: 1 ◦Start: 0 ◦Sliding Window −→ Sliding Sample Window Size: 4 Este monitor calcula la varianza de un escalar (el rendimiento) en un paso de tiempo en cada celda teniendo en cuenta todas las celdas correspondiente a la part escogida. Dado que el valor del rendimiento es el mismo la varianza de un paso de tiempo es cero. Se podría haber escogido cualquier part. Cabe señalar que el resultado no es un número, sino un número en cada celda (una field function. Si se escoge un número mayor que uno para Sliding Sample Window Size se calcula la varianza en los últimos cuatros pasos de tiempo en cada celda. De esta forma se sabrá si se esta llegando a un estacionario o no. Si el valor resultante es cero querrá decir que el rendimiento se ha mantenido constante en los últimos cuatro pasos. Además, se puede publicar la media por si hubiera que realizar alguna operación con ella (Publish Mean). 59
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Para tener un número solo queda calcular la media, ya que todas las celdas tienen el mismo valor. •Report: Varianza rendimiento ◦Tipo: Volume average report ◦Scalar Field Function: Variance of Rendimiento ◦Parts: Rotor ◦Monitor: Trigger:Time Step Frecuency: 1 Start: 0 La misma idea se va a aplicar en el siguiente monitor. Hasta ahora se tiene una función del tiempo (’Tiempo acelerador’) que indica los tiempos a los cuales se debe acelerar y en cuales no mediante picos en instantes determinados. Sin embargo, lo que se necesita es una función cuyo valor sea el máximo de la función anterior desde el inicio hasta ese instante de la simulación, para posteriormente comparar ese valor con el tiempo actual. •Monitors =⇒New Monitor −→ Field Max Se renombra como ’max tiempo acelerador’. ◦Definición: Field maximum monitor ◦Part: Rotor ◦Sample Count: 100 ◦Enable Sliding Sample Window:X ◦Field Function: Report: tiempo acelerador ◦Trigger:Time Step ◦Frecuency: 1 ◦Start: 0 ◦Sliding Window −→ Sliding Sample Window Size: 100 Se ha escogido 100 para Sliding Sample Window Size porque se ha supuesto que entre aceleración y aceleración de la turbina no va a haber más de 100 pasos de tiempo, siendo por tanto más que suficientes. •Report: Tiempo maximo ◦Tipo: Volume average report ◦Scalar Field Function: Max of Report: tiempo acelerador ◦Parts: Rotor •Monitor: ◦Trigger:Time Step ◦Frecuency: 1 ◦Start: 0 Otros Se van a mostrar alguno más que puede tener información extra. Por ejemplo, se han creado reports en los que se calcula la desviación típica de la componente perpendicular a la cuerda en varios perfiles. Este valor puede dar una idea del desprendimiento. EL siguiente ejemplo muestra la seccioón 14. •Report: Desviación 14 ◦Tipo: Surface standard deviation report ◦Scalar Field Function: Velocity in 3D-CAD Model 1 14[j] ◦Parts:[Copy of perfil 14] 60
4.10. CONVERGENCIA-RESIDUALS ◦Monitor: Trigger:Time Step Frecuency: 1 Start: 0 También se puede aproximar el número de Reynolds en una sección. Por ejemplo, también en la sección 14 sería como sigue: •Report: Desviación 14 ◦Tipo: Line integral report ◦Scalar Field Function: Reynolds 14 ◦Parts:[Copy of perfil 14] ◦Monitor: Trigger:Time Step Frecuency: 1 Start: 0 Por último, un report para mostrar el máximo módulo del gradiente de velocidades para poder observar como varía al refinar la malla. Igual se podría hacer con el valor de y+. •Report: Máximo gradiente velocidades ◦Tipo: Maximum value report ◦Scalar Field Function: Gradiente velocidades ◦Parts: [Bala menos rotor;Rotor] ◦Monitor: Trigger:Time Step Frecuency: 1 Start: 0 Otra representación interesante se consigue guardando el rendimiento cada vez que se llega a estacionario. Para ello, a partir del report del rendimiento se escoge Create Monitor and Plot from Report y en lugar de escoger Time Step en el Monitor como Trigger, se escoge Update Event yGuardar scenes. De esta forma solo guardará el rendimiento cuanto el criterio del Update Event se cumpla (se verá a continuación). Se puede cambiar el eje de abscisas para mostrar Rendimiento-Velocidad angular. Para terminar, se pueden crear plots con varias curvas para representar por ejemplo Reynolds en varias secciones o la desviación de igual forma. Para ello: Plots =⇒New Plot −→ Monitor Plot Monitor Plot 1 −→ Data Series =⇒Add Data Se escogen los monitors que interesen y se elige ’Guardar scenes’ como trigger si se quieren únicamente los valores en lo que se considera estacionario. 4.10. Convergencia-Residuals En esta sección se trata la existencia por defecto de los residuos (residuals), tanto en monitors como en plots, así como la convergencia de la solución. Los residuos consisten en varias cantidades cuya finalidad es tratar de aclarar sobre una posible convergencia o no. En general, la tendencia que seguirían los residuos de una solución que converga sería decreciente. Los residuos dependen bastante de los primeros pasos en la simulación, ya que están normalizados, por lo que más que valores absolutos de estos, conviene observar cuantos órdenes de magnitud dis61
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN minuyen en un mismo paso de tiempo. La tendencia entre dos pasos de tiempo no es relevante, debe disminuir en un mismo paso de tiempo. En general, tres órdenes de magnitud es suficiente para poder asegurar la convergencia, aunque como luego se verá este criterio es muy relativo. Los residuos que aparecen por defecto son: Continuity X-momentum Y-momentum Z-momentum Tke (Turbulent kinetic energy) Tdr (Turbulent dissipation rate) Existen dos tipos de errores: disipativos, de primer orden, y dispersivos, de segundo. Los primeros tienden a estabilizar la solución haciendo que los residuos disminuyan con cada iteración, mientras que los segundos pueden hacer que a medida que aumentan las iteraciones la tendencia sea aumentarlos. La mejor solución para juzgar la convergencia es por tanto seguir una o varias magnitudes del problema. En la simulación se seguirá principalmente el rendimiento, aunque también se observarán los campos de presiones y las líneas de corriente entre otras cosas. Si no se aprecian cambios bruscos y los valores tienden a un valor concreto se podrá considerar que la solución converge, no sin antes razonar si la solución tiene sentido físico. Por ejemplo, no puede salir un rendimiento mayor que 1, ya que se obtendría más energía del viento que la dada por el límite de Betz. Como ya se ha explicando anteriormente, los criterios que se establecerán para la convergencia van a variar dependiendo de si la simulación calculará un barrido o no. 4.11. Update events Esta herramienta se encuentra dentro de Tools. El objetivo es tener una forma de guardar los resultados obtenidos (plots yscenes). Como ya se advirtió antes, hay ciertos problemas en guardar más de un plot haciendo uso de esto, aunque si funciona correctamente para las scenes. Tools −→ Update Events −→ New Event −→ Logic Event Se renombra como ’Guardar scenes’ y se escoge And como Logic Operator. Se le dará un margen de 5 segundos (5 pasos) antes de poder considerar que converge la solución. Guardar scenes −→ Update Events =⇒New Event −→ Monitor Range Se repite esta operación dos veces: •Monitor Range 1 ◦Monitor: varianza rendimiento Monitor ◦Range Operation: <= ◦Range Value: 0.001 El valor de Range Value será 0.001 si la simulación es de un barrido, será menor (5e-4 por lo general) si es de una sola velocidad de giro. Estos valores se han elegido tras haber observado como varía el rendimiento. Cuanto mayor sean los valores menos precisa será la solución. •Monitor Range 2 ◦Monitor:Physical Time ◦Range Operation: >= ◦Range Value: 5.0 s Lo más natural sería elegir Monitor Asymptote. Esta opción se activa si un número de muestras (a elegir) no difiere más de un valor (a escoger) entre su máximo y su mínimo. El problema es que para trabajar con las field functions no se puede incluir esta herramienta, y en su lugar se modelo 62
4.12. STOPPING CRITERIA haciendo uso de la varianza del rendimiento. Por tanto, para que los resultados se guarden acorde con la definición de tales funciones se ha modelado de esta forma. 4.12. Stopping Criteria En esta sección se explica los criterios por los cuales la simulación para, tanto para pasar de un paso de tiempo a otro, como por completo. Como ya se ha explicado en la sección 4.11, se va a usar la varianza en lugar de utilizar la herramienta asintótica. Se van a diferenciar dos tipos de criterios: uno que afecta a las inner iteration y otro que afecta a las outer iteration. Esto es, uno que para las iteraciones dentro de un paso de tiempo, y otro que para la simulación por completo. El primer tipo se usará en ambas simulaciones (sea un barrido o no), mientras que el segundo solo se utilizará si la simulación es para remallar. El remallado se efectuará (según el macro) una vez se pare la simulación por completo con este criterio. En la simulación del barrido no se ha establecido ningún criterio que la pare por completo, se ha parado manualmente una vez se ha visto que el rendimiento empieza a disminuir (desprendimiento) claramente. Para parar las iteraciones en un paso de tiempo se ha optado por el camino más simple. Se podrían haber creado criterios en función de los residuos, pero finalmente se ha establecido un número máximo de iteraciones suficientemente alto (mayor si la simulación es de una sola velocidad de giro). Como ya se adelantó, este número es de 50 si solo se simula una velocidad y 30 si es un barrido. Se desactivan los criterios de tiempo máximo y número de pasos. Maximum Inner Iterations •Maximum Inner Iterations: 10 •Logical Rule: Or Para parar la simulación se han usado los mismos criterios que en Update Event para que sea congruente. Stopping Criteria =⇒Create New Criterion −→ From Monitor... Se crean dos, uno se renombra como ’Physical Time Criterion’ y otro como ’varianza rendimiento Monitor Criterion’, que solo se activarán en caso de que la simulaión no sea un barrido. Physical Time Criterion •Monitor:Physical Monitor •Criterion Option:Maximum •Logical Rule:And •Stop Inner Iteration:X •Stop Outer Iteration:X •Maximum Value: 5.0 s varianza rendimiento Monitor Criterion •Monitor: varianza rendimiento Monitor •Criterion Option:Minimum •Logical Rule:And •Stop Inner Iteration:X •Stop Outer Iteration:X •Maximum Value: 5e-4 63
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN 4.13. Scenes En esta sección se explican los pasos para crear scenes, que consisten en representaciones sobre las regiones. Estas pueden basarse en cada celda y mostrar un valor escalar, un vector o una representación de la malla. También pueden mostrar la geometría simplemente. Para crear por ejemplo una representación de un escalar: Scenes =⇒New −→ Scalar Cada scene tiene dos nodos principales: Displayers yAttributes. Para simplificar la explicación se van a destacar las dos únicas opciones relevantes de Attributes, que serán prácticamente las mismas para todas las representaciones. Scenes −→ Scalar 1 −→ Attributes −→ View •Focal Point: Vector que especifica la posición del foco. •Position: Vector que especifica la posición de la cámara. •View Up: Vector que especifica la dirección en la que apunta la cámara. Las tres magnitudes anteriores son vectores que definen la vista de la scene. Lo más comodo es ajustarla con el ratón y si se quiere usar siempre la misma copiar y pegar esos vectores en las diferentes scenes. •Coordinate System: Específica el sistema al que están referidos los vectores anteriores. Por ejemplo, si se quiere tener una sección de la turbina, habrá que escoger unos ejes que roten junto a la turbina. •Projection Mode: ofrece la posibilidad de una proyección paralela o con ’perspectiva’, esto es, no manteniendo los paralelismos para permitir visualizar más detalles. Scenes −→ Scalar 1 −→ Attributes −→ Update Al igual que en los plots permite escoger cuando se deben guardar. Se pueden escoger varios formatos, entre los que se encuentra .sce. Este formato esta asociado al propio software y permite guardar una misma scene en diferentes momentos en un mismo archivo. Para ello hay que marcar la opción Append. El resto de opciones son equivalentes a las que ya se explicaron. En las propiedades de cada displayer hay que seleccionar Volume Mesh en la característica Representation para que se aplique los resultados a las celdas, excepto si se quiere mostrar la geometría, en cuyo caso habrá que seleccionar Geometry. Primero se van a crear scenes para ver el mallado, por tanto se selecciona Mesh. Se van a renombrar según lo que vayan a mostrar. Se van a mostrar solo las características que no están por defecto. Malla turbina −→ Displayers −→ Displayer 1 −→ Parts −→ ’Turbina’ y los tres ’Borde de ataque’ 64
4.13. SCENES Figura 4.23: Malla turbina Malla longitudinal −→ Displayers −→ Displayer 1 −→ Parts −→ Longitudinal Figura 4.24: Malla longitudinal Malla longitudinal −→ Displayers −→ Plano 9 −→ Parts −→ Longitudinal En este caso se ha puesto como ejemplo el perfil del plano 9. En View hay que seleccionar los ejes asociados a la turbina. 65
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Figura 4.25: Malla perfil Ahora se van crear aquellas scenes asociadas a una magnitud escalar. En general las opciones son las mismas, excepto un conjunto de opciones que antes no aparecía y que está asociado a la magnitud que se quiere representar. Scalar Scene 1 −→ Displayers −→ Displayer 1 −→ Scalar Field •Function: Se selecciona la magnitud escalar a representar. •Auto Range: Permite seleccionar varias opciones según se quiera mostrar los valores máximos y mínimos, o por el contrario, se quieran escoger manualmente. •Clip: Permite mostrar solo valores por encima del mínimo, por debajo del máximo o entre ambos. Junto con la opción anterior permite representar solo aquellas celdas que estén entre los valores que interesen. Dado que se mostrarán los resultados más adelante no se van a describir todas las scenes creadas. A continuación se van a explicar los detalles básicos para crear una animación con las líneas de corriente. Para ello se puede crear cualquier tipo de scene. Para ver la turbina de forma nítida, se creará una de tipo Geometry, para ver la superficie sin necesidad de ver la malla de esta. Geometry 1 −→ Displayers =⇒New Displayer −→ Streamline Streamlines 1 Se selecciona en Mode la opción Tubes, para que las líneas de corriente se muestren como tubos y sean fáciles de visualizar. Streamlines 1 −→ Parts Se seleccionar la derived part correspondiente, en este caso solo hay una: Streamline. Streamline 1 −→ Animation −→ Tubes Streamline 1 −→ Markers •Display Start Point Markers:X •Display End Point Markers:X Por último queda explicar lo más relevante de las Vector Scenes. Vector Scene 1 −→ Displayers −→ Displayer 1 −→ Vector Field Este nodo tiene las mismas opciones que en las Scalar Scenes, con la única diferencia que cuando antes se hacía referencia a un valor, ahora se hace referencia al módulo del vector. Vector Scene 1 −→ Displayers −→ Displayer 1 −→ Line Integral Convolution Con esta opción se pueden representar las líneas de corriente a partir de los vectores de velocidad. Para personalizar como deben mostrarse estas existen en este nodo opciones con las que hacerlo. Para los resultados que se mostrarán se han escogido las siguientes: 66
4.14. TABLES •Enhanced:X •Number Steps: 50 •Blending Factor: 0.5 •Step Size (in px): 0.6 4.14. Tables Esta herramienta se encuentra en Tools. Se usará para remallar en función de alguna magnitud. Tools −→ Tables =⇒New Table −→ XYZ Internal Table Se van a crear dos tablas, una para cada región, que refinará la malla en base al módulo del gradiente de velocidades. Para ello se creo una función en Field Functions. De esta forma a cada celda, que se guarda con su respectiva posición (x-y-z), se le asocia un nuevo valor del tamaño de la malla en esa posición. Es posible que haya que disminuir el valor Minimum Surface Size manualmente, ya que si este fue el limitante en la malla original impedirá que se disminuya el tamaño. Otra opción es escoger un valor suficientemente pequeño desde el principio. Se podría crear una tabla para refinar según el valor de wally+, pero dado los recursos computacionales se ha visto que disminuir este valor significativamente es excesivamente caro. No obstante, esto se puede ignorar dado que se va a refinar la malla y se va a ver que la solución converge. Como ejemplo se muestra una de ellas, renombrada como ’Bala menos rotor gradiente velocidades’. Tools −→ Tables −→ Bala menos rotor gradiente velocidades •Scalars: [0_Refinar malla gradiente velocidades] •Parts: Bala menos rotor •Representation:Volume Mesh •Update ◦Enabled:X ◦Autoextract:X ◦Trigger:Time Step 4.15. Proceso de remallado Se van a resumir los pasos para remallar de manera adecuada. 1. Geometry −→ Parts −→ Rotor =⇒Transform −→ Coordinate System Dado que la turbina habrá rotado antes de remallar por primera vez, y que la malla está referida a cada una de las part, habrá que hacer alguna modificación porque cada part esta referida a los ejes globales. 67
CAPÍTULO 4. STAR-CCM+: PREPARACIÓN DE LA SIMULACIÓN Con esto se evita que la turbina quede rotada respecto al sistema de coordenadas Body 1-Csys, lo cual haría entre otras cosas que los cortes con los planos (derived parts) no dieran lugar al perfil aerodinámico escogido, sino a un corte transversal o directamente no cortara a la turbina. Debe quedar como se observa en la siguiente figura: 2. Geometry −→ Operations −→ Malla:Rotor −→ Meshers −→ Surface Remesher 3. Geometry −→ Operations −→ Malla:Rotor −→ Meshers −→ Trimmed Cell Mesher Se seleccionan las tablas preparadas anteriormente. Adicionalmente, se puede seleccionar el sistema de coordenadas Body 1-Csys, para que la malla se alineé con este en lugar de con los ejes globales (esto solo se haría en ’Malla:Rotor’). Igual se hace para la otra región si es necesario. (a) Surface Remesher (b) Trimmed Cell Mesher 68
5.1. AH 93-W-145 Rendimiento Monitor Plot Rendimiento Monitor 0 5 10 15 20 25 30 35 40 Physical Time (s) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 Rendimiento Monitor Figura 5.4: Rendimiento Residuals Residual 1e-04 0.001 0.01 0.1 1 Iteration 100 200 300 400 500 600 700 800 900 1000 1100 1200 Continuity X-momentum Y-momentum Z-momentum Tke Tdr Figura 5.5: Residuos 75
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS Por otro lado se muestra la evolución de la malla en la figura 5.7. La diferencia en el parámetro Trimmer Surface Growth Rate en las dos superficies (’Borde de ataque’ y ’Borde de salida’) es lo que más se aprecia. Se muestran en las figuras 5.8 algunas magnitudes como la aceleración angular, la velocidad angular, la componente zdel momento que ejerce el fluido sobre la turbina (obviamente tendrá una igual y contrario) y la suma de esta componente de todos los momentos. Los resultados son los esperados, una aceleración angular prácticamente nula y un momento respecto al eje de giro negativo (pues este momento se define como el que aporta el fluido). La representación de la velocidad parece no mostrar nada, pero es un problema que presenta el software cuando el valor de la magnitud apenas varía. De hecho, el eje de ordenadas muestra el mínimo y el máximo valor de esta magnitud, y como se puede apreciar, son prácticamente iguales. Esto es, la velocidad es en la práctica constante, como se esperaba. Por último, la suma de momentos respecto al eje zdebería ser nula aproximadamente. En algunos casos esto podría ponerse en duda a la vista de los resultados.En cualquier caso, puede considerar que esto es así dado que los órdenes de magnitud de los momentos por separados son de 1×105, y al sumarlos hasta un resultado de orden 1×102se puede considerar despreciable. El hecho de que todos ellos, excepto la resistencia en la dirección del eje, tiendan a un valor con las tres mallas confirma la convergencia de la solución. Debe explicarse la razón por la que la resistencia no converge. Las fuerzas que originan la sustentación de un perfil se deben al gradiente de presiones alrededor del perfil. Sin embargo, las fuerzas de resistencia se deben principalmente (en dos dimensiones) a las fuerzas viscosas en la capa límite. Por tanto, para modelar correctamente este tipo de efectos la capa límite debe estar correctamente modelada (malla más fina). En el caso tridimensional además se superponen velocidades transversales al perfil debidas también al gradiente de presiones a lo largo de la pala, así como torbellinos en las puntas de esta principalmente. También se muestra el máximo módulo del gradiente de velocidades en la figura 5.12. En principio esta magnitud no tiene por que tender a la misma solución con diferentes mallas, ya que el gradiente depende de las velocidades en celdas contiguas, y según sean estas la velocidad se aproximará mejor (de manera más continua) o no. Junto a esto, se muestra la field function: ’Gradiente de velocidades adim’. Se muestran tres imágenes (figura 5.9), una para el estacionario de cada malla, en la misma sección que antes. Según su definición, el valor es más pequeño cuanto mayor es el gradiente. Es por ello que cerca de la turbina tiene menores valores. Es interesante a su vez observar el campo de presiones (figura 5.10, que permite ver el pico de succión. Se distinguen perfectamente extradós e intradós. Se muestran también las líneas de corriente en varias secciones (figura 5.11, así como las velocidades respecto a la turbina para poder observar las posibles zonas de desprendimiento. En este caso dado que el ángulo de ataque que ve la corriente es φ= 10◦, no se observa demasiado este fenómeno, aunque en el borde de salida del extradós se ve una zona en las que las velocidades disminuyen y parece intuirse cierta recirculación. Esto se apreciara mejor cuando se muestren con ángulos de ataque mayores. En este caso solo se van a mostrar para la malla más fina. Además se muestra la potencia real obtenida de la turbina, esto es, el producto del límite de Betz por la potencia teórica que se puede obtener (figura 5.13). 76
5.1. AH 93-W-145 Elementos bala menos rotor Monitor Plot Elementos bala menos rotor Monitor 320000 320500 321000 321500 322000 322500 323000 Physical Time (s) 5 10 15 20 25 Elementos bala menos rotor Monitor (a) ’Bala menos rotor’ Elementos rotor Monitor Plot Elementos rotor Monitor 4200000 4400000 4600000 4800000 5e+06 5200000 5400000 Physical Time (s) 5 10 15 20 25 Elementos rotor Monitor (b) ’Rotor’ Figura 5.6: Número de celdas 77
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Malla original (b) Primer remallado (c) Segundo remallado Figura 5.7: Malla (Sección 9) 78
5.1. AH 93-W-145 6-DOF Body Angular Acceleration Monitor Plot 6-DOF Body Angular Acceleration Monitor (radian/s^2) -5.5e-12 -5e-12 -4.5e-12 -4e-12 -3.5e-12 -3e-12 -2.5e-12 -2e-12 -1.5e-12 -1e-12 -5e-13 0 Physical Time (s) 5 10 15 20 25 6-DOF Body Angular Acceleration Monitor (a) Aceleración angular 6-DOF Body Angular Velocity 1 Monitor Plot 6-DOF Body Angular Velocity 1 Monitor (rpm) -15.0000000005 -15.00000000045 -15.0000000004 -15.00000000035 -15.0000000003 -15.00000000025 -15.0000000002 -15.00000000015 -15.0000000001 Physical Time (s) 5 10 15 20 25 6-DOF Body Angular Velocity 1 Monitor (b) Velocidad angular 79
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS 6-DOF Body Moment z Monitor Plot 6-DOF Body Moment z Monitor (N-m) -40000 -35000 -30000 -25000 -20000 -15000 -10000 -5000 0 Physical Time (s) 5 10 15 20 25 6-DOF Body Moment z Monitor (c) Momento de eje z 6-DOF Body Sumatorio Momentos Monitor Plot 6-DOF Body Sumatorio Momentos Monitor (N-m) -100 -80 -60 -40 -20 0 20 40 60 80 100 Physical Time (s) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 6-DOF Body Sumatorio Momentos Monitor (d) Suma de momentos de eje z 80
5.1. AH 93-W-145 6-DOF Body Force 1 Monitor Plot 6-DOF Body Force 1 Monitor (N) -20 -15 -10 -5 0 5 10 15 Physical Time (s) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 6-DOF Body Force 1 Monitor (e) Resistencia en dirección z Figura 5.8: DFBI 81
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Malla original (b) Primer remallado (c) Segundo remallado Figura 5.9: Gradiente de velocidades adimensional (Sección 9) 82
5.1. AH 93-W-145 (a) Malla original (b) Primer remallado (c) Segundo remallado Figura 5.10: Presión (Sección 9) 83
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Sección 4 (b) Sección 9 84
5.1. AH 93-W-145 tiempo acelerador Monitor Plot tiempo acelerador Monitor (s) 0 20 40 60 80 100 120 140 160 180 200 220 Physical Time (s) 20 40 60 80 100 120 140 160 180 200 220 tiempo acelerador Monitor (c) Tiempo de aceleración 6-DOF Body Angular Velocity 1 Monitor Plot 6-DOF Body Angular Velocity 1 Monitor (rpm) -22 -20 -18 -16 -14 -12 -10 -8 -6 -4 -2 Physical Time (s) 20 40 60 80 100 120 140 160 180 200 220 6-DOF Body Angular Velocity 1 Monitor (d) Velocidad angular 91
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Velocidad 1 5.18), desde ω= 0,5rpm hasta el mayor valor de la velocidad pasando por puntos intermedios. Se puede apreciar como la adherencia de la capa límite evoluciona, no pudiendo distinguirse y observando grandes torbellinos a velocidades bajas , y empezando a desprenderse una vez se ha alcanzado el punto óptimo en cuanto a rendimiento. Como ya se explicó, se ve en los resultados que una malla más fina podría aproximar mejor la realidad, ya que la velocidad mínima aumenta cuando debería mantenerse cercana a cero cerca de la pared. Para las mismas velocidades se van a mostrar el campo de presiones (figura 5.19). Se puede observar como el punto de remanso avanza hacia el borde de ataque, así como la evolución de la zona de succión, que se va haciendo mayor a medida que la capa límite se encuentra más adherida. Es por ello que en la última imagen vuelve a hacerse más pequeña. Se van a mostrar también magnitudes relativas a la computación como se hizo antes. Se puede ver los tiempos de simulación en la figura 5.22. En ambos la pendiente se mantiene aproximadamente constante, ya que la malla lo es. También se muestran el número de celdas de cada región (figuras 5.21. Como se explicó antes, la única variación se debe a la interfase entre ambas regiones. Por último se muestra la memoria utilizada (5.20), que como cabía esperar, es mucho menor dado que la malla es menos fina. 92
5.1. AH 93-W-145 (b) Velocidad 2 (c) Velocidad 3 93
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (d) Velocidad 4 (e) Velocidad 5 94
5.1. AH 93-W-145 (f) Velocidad 6 Figura 5.18: Líneas de corriente (a) Velocidad 1 95
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (b) Velocidad 2 (c) Velocidad 3 96
5.1. AH 93-W-145 (d) Velocidad 4 (e) Velocidad 5 97
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (f) Velocidad 6 Figura 5.19: Campo de presiones Maximum Memory Report Monitor Plot Maximum Memory Report Monitor (GiB) 5 5.5 6 6.5 7 7.5 8 8.5 9 9.5 10 10.5 Physical Time (s) 20 40 60 80 100 120 140 160 180 200 220 Maximum Memory Report Monitor Figura 5.20: Memoria requerida 98
5.1. AH 93-W-145 Elementos bala menos rotor Monitor Plot Elementos bala menos rotor Monitor 320000 320500 321000 321500 322000 322500 323000 323500 Physical Time (s) 20 40 60 80 100 120 140 160 180 200 220 Elementos bala menos rotor Monitor (a) ’Bala menos rotor’ Elementos rotor Monitor Plot Elementos rotor Monitor 4092500 4093000 4093500 4094000 4094500 4095000 4095500 4096000 Physical Time (s) 20 40 60 80 100 120 140 160 180 200 220 Elementos rotor Monitor (b) ’Rotor’ Figura 5.21: Número de celdas 99
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS Total Solver Elapsed Time Monitor Plot Total Solver Elapsed Time Monitor (hr) 0 5 10 15 20 25 30 35 Physical Time (s) 20 40 60 80 100 120 140 160 180 200 220 Total Solver Elapsed Time Monitor (a) Tiempo real de simulación Total Solver CPU Time Monitor Plot Total Solver CPU Time Monitor (hr) 20 40 60 80 100 120 140 Physical Time (s) 20 40 60 80 100 120 140 160 180 200 220 Total Solver CPU Time Monitor (b) Suma de los tiempos de todos los procesadores Figura 5.22: Tiempos de simulación 100
5.1. AH 93-W-145 5.1.3. φ= 15◦,ω= 15 rpm Por último, se va ha simulado este perfil para otro ángulo de ataque efectivo. Rendimiento Monitor 2 Plot Rendimiento estacionario Monitor 0 10 20 30 40 50 60 70 80 90 100 6-DOF Body Angular Velocity 1 Monitor (rpm) -22 -20 -18 -16 -14 -12 -10 -8 -6 -4 -2 0 Rendimiento estacionario Monitor (a) Rendimiento −V elocidad angular Rendimiento Monitor Plot Rendimiento Monitor 0 5 10 15 20 25 30 35 40 Physical Time (s) 5 10 15 20 25 30 Rendimiento Monitor (b) Rendimiento ω= 5 rpm Figura 5.28 Se vuelven a mostrar las líneas de corriente (figura 5.29) y las presiones (figura 5.30) para la misma sección (sección 9) para tres velocidades ordenadas de menor a mayor. 107
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Velocidad 1 (b) Velocidad 2 Por último se muestra el número de elementos de cada región (figura 5.31 para ’Rotor’ y 5.32 para ’Bala menos rotor’). Las funciones usadas para refinar la malla son: (${Gradiente velocidades adim}<=4) ? 0.1*${Base size} : ((${Gradiente velocidades adim}>6) ? 2*${Base size} : 0) (${Gradiente velocidades adim}<=4) ? 0.02*${Base size} : ((${Gradiente velocidades adim}>6) ? 2*${Base size} : 0) 108
5.1. AH 93-W-145 (c) Velocidad 3 Figura 5.29: Líneas de corriente (a) Velocidad 1 109
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (b) Velocidad 2 (c) Velocidad 3 Figura 5.30: Campo de presiones 110
5.1. AH 93-W-145 Elementos rotor Monitor Plot Elementos rotor Monitor 4012500 4013000 4013500 4014000 4014500 4015000 Physical Time (s) 50 100 150 200 250 Elementos rotor Monitor (a) ’Rotor’ sin remallado Elementos rotor Monitor Plot Elementos rotor Monitor 4420000 4430000 4440000 4450000 4460000 4470000 Physical Time (s) 5 10 15 20 25 30 Elementos rotor Monitor (b) ’Rotor’ con remallado Figura 5.31 111
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS Elementos bala menos rotor Monitor Plot Elementos bala menos rotor Monitor 321000 321500 322000 322500 323000 323500 Physical Time (s) 50 100 150 200 250 Elementos bala menos rotor Monitor (a) ’Bala menos rotor’ sin remallado Elementos bala menos rotor Monitor Plot Elementos bala menos rotor Monitor 320000 321000 322000 323000 324000 325000 326000 327000 328000 Physical Time (s) 5 10 15 20 25 30 Elementos bala menos rotor Monitor (b) ’Bala menos rotor’ con remallado Figura 5.32 112
5.1. AH 93-W-145 5.1.4. Comparación de los tres casos con perfil AH 93-W-145 Se muestran en la figura 5.33 las curvas de rendimiento de los casos anteriores. En función de la velocidad del viento más probable que daría la rosa de los vientos habría que estudiar que curva ofrece un mayor rendimiento. -25 -20 -15 -10 -5 0 -5 0 5 10 15 20 25 30 35 φ=10º,ω=15 rpm φ=10º,ω=5 rpm φ=15º,ω=15 rpm Figura 5.33: Curvas de rendimiento Como ya se ha visto en cada caso, la malla parece no ser lo suficientemente fina para modelar el desprendimiento. Además, solo se ha mostrado la función definida para refinar a alguna velocidad de referencia (no en casos cercanos al desprendimiento), por lo que puede parecer que su definición no serviría tampoco para refinar la malla a otras velocidades. Para ello se muestra a modo de ejemplo la figura 5.34a, en la que se puede observar que a otras velocidades se refinaría en zonas en las que se desprende la corriente (φ= 15◦,ω= 15 rpm), aunque quizás sería mejor no haber introducido el logaritmo en la definición y hacerlo adimensionalizarlo de manera lineal. Por otro lado, se podría usar alguna otra magnitud como por ejemplo la energía cinética turbulenta. Se muestra un ejemplo de esta magnitud en la figura 5.34b para el mismo caso en el tercer caso (φ= 15◦,ω= 15 rpm). Para ello se crea un report que guarde el mayor valor de la energía cinética turbulenta para poder adimensionalizar esta magnitud. 113
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Gradiente velocidades adimensional (b) Energía cinética turbulenta 114
5.2. FX 84-W-218 5.2. FX 84-W-218 Se muestran ahora los resultados del segundo perfil aerodinámico para el caso en que φ= 10◦y ω= 15 rpm. Se muestran los rendimientos en la figura 5.38. En este caso se ha buscado no solo aumentar el número de elementos entre mallado y mallado, sino que la malla se vea refinada en una zona mayor, no de manera tan localizada. Como se puede observar, el rendimiento parece converger a un valor ligeramente superior (2 décimas) con la malla más fina. Esto se podría esperar ya que se puede modelar mejor el gradiente de presiones. Parece razonable suponer que el valor encontrado estará cercano al real y que por lo general, si converge será a un límite inferior (obviando otro tipo de efectos como el desprendimiento). En la figura 5.36 se pueden ver los diferentes mallados, siendo notable el aumento de elementos. Las funciones utilizadas para remallar se muestran a continuación: (${Gradiente velocidades adim}<=4) ? 0.1*${Base size} : ((${Gradiente velocidades adim}>6) ? 2*${Base size} : 0) (${Gradiente velocidades adim}<=4) ? 0.005*${Base size} : ((${Gradiente velocidades adim}>6) ? 2*${Base size} : 0) En la figura 5.37 se puede ver la sección longitudinal que contiene las dos regiones. Se puede observar como aumenta el número de celdas, incluyendo en la zona posterior a la raíz donde existirá recirculación pronunciada. También se muestra el máximo gradiente de presiones con las diferentes mallas en la figura 5.35, donde se puede observar como aumenta debido a que los elementos se hacen más pequeños, por lo que el cociente de la variación de la velocidad entre su tamaño característico aumenta. Se vuelven a mostrar las líneas de corriente (figura 5.39) y las presiones (figura 5.40) para la sección 9 y tres velocidades de rotación ordenadas de menor a mayor. Además se muestran ambas para la velocidad de referencia (caso más cercano a la segunda velocidad) con la malla más fina. En tal caso, se puede ver como el campo de presiones es mucho más suave y mejor modelado en el extradós. Por último también se muestra el número de elementos de cada región (figura 5.41 para ’Rotor’ y 5.42 para ’Bala menos rotor’). 115
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS Mˆ¡ximo gradiente Monitor Plot Mˆ¡ximo gradiente Monitor (/s) 50000 55000 60000 65000 70000 75000 80000 85000 90000 95000 1e+05 105000 Physical Time (s) 5 10 15 20 25 30 35 Mˆ¡ximo gradiente Monitor Figura 5.35: Máximo gradiente de velocidad 116
5.2. FX 84-W-218 (c) Velocidad 3 (d) Velocidad de referencia Figura 5.40: Campo de presiones 123
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS Elementos rotor Monitor Plot Elementos rotor Monitor 1302000 1302500 1303000 1303500 1304000 1304500 1305000 1305500 Physical Time (s) 20 40 60 80 100 120 140 160 180 Elementos rotor Monitor (a) ’Rotor’ sin remallado Elementos rotor Monitor Plot Elementos rotor Monitor 3800000 4e+06 4200000 4400000 4600000 4800000 5e+06 Physical Time (s) 5 10 15 20 25 30 35 Elementos rotor Monitor (b) ’Rotor’ con remallado Figura 5.41 124
5.2. FX 84-W-218 Elementos bala menos rotor Monitor Plot Elementos bala menos rotor Monitor 320000 320500 321000 321500 322000 322500 323000 Physical Time (s) 20 40 60 80 100 120 140 160 180 Elementos bala menos rotor Monitor (a) ’Bala menos rotor’ sin remallado Elementos bala menos rotor Monitor Plot Elementos bala menos rotor Monitor 320000 320500 321000 321500 322000 322500 323000 323500 324000 324500 325000 325500 Physical Time (s) 5 10 15 20 25 30 35 Elementos bala menos rotor Monitor (b) ’Bala menos rotor’ con remallado Figura 5.42 125
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS 5.3. MH 102 Se muestran ahora los resultados del tercer perfil aerodinámico para el caso en que φ= 10◦yω= 15 rpm. Se muestra el rendimiento de la misma manera en la figura 5.43. En este caso se observa una disminución ligera del valor al que converge con el segundo remallado. Comparando el valor obtenido remallando y el valor de la curva para ω= 15 rpm se observa una diferencia apreciable. Se puede concluir que la malla original en este caso no es adecuada, y que por tanto, habría que buscar una malla mejor intentando minimizar los recursos necesarios para simularla. De hecho, usando las mismas funciones apenas se aprecia la mejora en la malla en la sección 9 (figura 5.44), aunque si se pueden apreciar en la figura 5.45 las diferentes dimensiones que asigna la malla a cada celda. En este caso si se observan zonas donde la malla se ha refinado. (${Gradiente velocidades adim}<=4) ? 0.1*${Base size} : ((${Gradiente velocidades adim}>6) ? 2*${Base size} : 0) (${Gradiente velocidades adim}<=4) ? 0.05*${Base size} : ((${Gradiente velocidades adim}>6) ? 2*${Base size} : 0) Se muestran las líneas de corriente (figura 5.46) y las presiones (figura 5.47) para la sección 9 y tres velocidades de rotación ordenadas de menor a mayor. Por último se muestra el número de elementos de cada región (figura 5.48 para ’Rotor’ y 5.49 para ’Bala menos rotor’). 126
5.3. MH 102 Rendimiento Monitor 2 Plot Rendimiento estacionario Monitor 0 10 20 30 40 50 60 70 80 90 100 6-DOF Body Angular Velocity 1 Monitor (rpm) -25 -20 -15 -10 -5 0 Rendimiento estacionario Monitor (a) Rendimiento −V elocidad angular Rendimiento Monitor Plot Rendimiento Monitor 0 5 10 15 20 25 30 35 40 Physical Time (s) 5 10 15 20 25 30 35 Rendimiento Monitor (b) Rendimiento ω= 15 rpm Figura 5.43 127
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Malla original (b) Primer remallado (c) Segundo remallado Figura 5.44: Malla (Sección 9) 128
5.3. MH 102 (a) Malla original (b) Segundo remallado Figura 5.45: Valores de la malla 129
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (a) Velocidad 1 (b) Velocidad 2 130
5.3. MH 102 (c) Velocidad 3 Figura 5.46: Líneas de corriente (a) Velocidad 1 131
CAPÍTULO 5. STAR-CCM+: ANÁLISIS DE RESULTADOS (b) Velocidad 2 (c) Velocidad 3 Figura 5.47: Campo de presiones 132