Repositorio Institucional de Documentos

Proyecto Fin de Carrera MODELADO DEL FLUJO DE AIRE A TRAVÉS DE GRUPOS DE ÁRBOLES Autor José Miguel Maldonado Jiménez Director Mathias Cehlin Ponente Guillermo Hauke Departamento de Mecánica de Fluidos Escuela de Ingeniería y Arquitectura 2012 RESUMEN I Resumen Los emplazamientos más apropiados para instalar un parque eólico, están ocupados actualmente o por motivos medioambientales no se puede hacer uso de los mismos. Por este motivo, es necesario comprobar la posibilidad de si otras localizaciones disponibles, donde el terreno es más complejo, son viables. Para este fin, se utilizan herramientas de Fluido dinámica Computacional o CFD (Computational Fluid Dynamics). Para la realización de este trabajo se ha seleccionado el software llamado OpenFOAM, un software libre que, a diferencia de otros comerciales, su código (programado en C++) es suficientemente flexible y el programa es gratuito. El proyecto básicamente consiste en desarrollar una aplicación que permita simular el efecto de los árboles y bosques sobre las corrientes de aire. Hasta el momento esta tarea se venía haciendo considerando los arboles como parte de la rugosidad del terreno, pero es recomendable añadir la distorsión producida por estos de manera separada a la rugosidad. Se puede separar el proyecto en dos partes, coincidiendo cada una con una aplicación a desarrollar: • El primer objetivo es desarrollar una aplicación que cargue los parámetros que definen el tipo de árboles y el tamaño de éstos. Esta herramienta de pre-procesado situará los diferentes grupos de árboles a lo largo de la malla del terreno en estudio. • El segundo objetivo consiste en programar un “solver” el cual resuelve una configuración teniendo en cuenta el efecto producido por el follaje sobre el flujo de aire. Gracias a estas aplicaciones, se podrá estudiar de manera más detallada si un terreno es adecuado para instalar aerogeneradores. La finalidad de este tipo de simulaciones es prever la velocidad con la que el viento llega al molino, para así estimar la producción energética anual del mismo. ÍNDICE III Índice. 1. Introducción. ......................................................................................................................... 1 1.1. Energía eólica. ............................................................................................................... 1 1.1. Motivación. ................................................................................................................... 3 1.2. Objetivos. ...................................................................................................................... 4 1.3. Método. ......................................................................................................................... 5 2. Fundamentos Teóricos. ......................................................................................................... 7 2.1. Turbulencia. ................................................................................................................... 7 2.1.1. Efecto de la turbulencia en las Ecuaciones del Transporte promediadas............. 7 2.2. Modelado Numérico. .................................................................................................. 10 2.2.1. Modelos Turbulentos .......................................................................................... 10 2.2.2. Aproximación de Boussinesq .............................................................................. 11 2.2.3. Modelo turbulento k-epsilon. ............................................................................. 12 2.3. Modelado de la fuerza de arrastre debida al follaje. .................................................. 13 2.3.1. Parámetros característicos. ................................................................................. 13 2.4.2 Implementación de la fuerza de arrastre. ........................................................... 14 3. OPENFOAM ......................................................................................................................... 15 3.1. Introducción. ............................................................................................................... 15 3.2. Aplicaciones en OpenFOAM. ....................................................................................... 16 3.2.1. Código principal de la aplicación, archivo “.C” .................................................... 16 3.2.2. Codificación de las ecuaciones en C++. ............................................................... 18 3.3. Casos en OpenFOAM. .................................................................................................. 19 3.3.1. Fichero “System” ................................................................................................. 20 3.3.2. Fichero “Constant” .............................................................................................. 21 3.3.3. Fichero “0”........................................................................................................... 22 4. Proceso y resultados ........................................................................................................... 23 4.1. Caso de estudio, la colina Bolund. .............................................................................. 23 4.2. Adición del follaje como variable en OpenFOAM. ...................................................... 25 4.3. Aplicación de pre-procesado. ...................................................................................... 27 4.3.1. Archivos “.map”................................................................................................... 27 4.3.2. Aplicación de pre procesado “leafIndexToFoam” ............................................... 28 4.3.3. Simulación de zonas arboladas. .......................................................................... 31 4.4. TreeFoam “solver” test. .............................................................................................. 33 ÍNDICE IV 4.4.1. Sin bosques o zonas arboladas definidas en el dominio. .................................... 34 4.4.2. Simulación con “treeFoam” añadiendo bosques al dominio. ............................. 36 5. Conclusión. .......................................................................................................................... 39 6. Líneas futuras. ..................................................................................................................... 41 7. Bibliografía. ......................................................................................................................... 43 APPENDIX I .................................................................................................................................. 45 A.I.1 Power of wind. ................................................................................................................ 45 A.I.2 Betz Law .......................................................................................................................... 48 APPENDIX II ................................................................................................................................. 49 Appendix II.1 ........................................................................................................................... 49 Appendix II.2 Experimental researches for obtaining the drag coefficient (CD) ..................... 52 APPENDIX III ................................................................................................................................ 55 Appendix III.1 leafIndexToFoam.C .......................................................................................... 55 Appendix III.2 CreateCdAlphaL.H ............................................................................................ 56 Appendix III.3 Readmapfile_CdAlphaL.H ................................................................................ 57 Appendix III.4 Check_inside_cell.H ................................................................................... 61 Appendix III.5 tree.map (Map file example) ........................................................................... 62 Appendix III.6 Cd and alphaL initial files. ................................................................................ 63 APPENDIX IV ................................................................................................................................ 67 Appendix IV.1 ControlDict. ...................................................................................................... 67 Appendix IV.2 fvSchemes ........................................................................................................ 68 Appendix IV.3 fvSolution ......................................................................................................... 70 Appendix IV.4 transportProperties ......................................................................................... 72 Appendix IV.5 turbulenceProperties ....................................................................................... 72 Appendix IV.6 points ............................................................................................................... 73 Appendix IV.7 faces. ................................................................................................................ 74 Appendix IV.8 boundary .......................................................................................................... 75 Appendix IV.9 owner ............................................................................................................... 76 APPENDIX V ................................................................................................................................. 79 APPENDIX VI ................................................................................................................................ 83 APPENDIX VI.1 U file ................................................................................................................ 83 APPENDIX VI.2 p file ................................................................................................................ 84 APPENDIX VI.3 Cd and alphaL files. ......................................................................................... 84 APPENDIX VII Memoria original en lengua inglesa ..................................................................... 91 ÍNDICE V Lista de Figuras Fig 1.1 Producción eléctrica en Suecia en 2009………………………………………………..…………………..… 2 Fig 1.2 Producción eléctrica en España en 2009………………………………………………………….….......... 2 Fig 1.3 Producción eléctrica en Dinamarca en 2009……………………………………………..…..……………. 3 Fig 1.4 La colina Bolund representada en OpenFOAM…….……..…………………………………….….……. 5 Fig 3.1 Estructura de OpenFOAM …………………………………………………………………………………………. 16 Fig 3.2 Estructura de una aplicación ……………………………………………………………………………….….… 17 Fig 3.3 Estructura de un caso en OpenFOAM …………………………………………………………………..…… 20 Fig 4.1 Vista panorámica de la colina Bolund …………………………………………………………………….…. 23 Fig 4.2 La colina Bolund ………………………………………………………………………………………………………….24 Fig 4.3 Malla correspondiente a la colina Bolund ……………………………………………………….…….….. 24 Fig 4.4 Proceso seguido para leer un archivo .map.……………………………………..…………………….…. 29 Fig 4.5 Representación de “alphal” en “Paraview” vista normal al eje ………………….……………... 32 Fig 4.6 Representación de “alphal” en “Paraview” vista normal al eje Y……………………..……….. 33 Fig 4.7 Representación de “alphal” en “Paraview” vista normal al eje X a -30 m…………………... 33 Fig 4.8 Representación de “alphal” en “Paraview” vista normal al eje X a 45 m …….………………33 Fig 4.9 Perfil de velocidades (m/s) normal al eje Z a 20 m ………………………………..……..…….……. 34 Fig 4.10 Perfil de velocidades (m/s) normal al eje Z a 40 m …………………….……………….………..… 35 Fig 4.11 Perfil de velocidades (m/s) normal al eje Z a 10 m ………………………………………...………. 35 Fig 4.12 Perfil de velocidades (m/s) a lo largo de la dirección del viento ……………………..…..…. 35 Fig 4.13 Perfil de velocidades (m/s) normal al eje Z a 20 m.……………………….…………..……..……. 37 Fig 4.14 Perfil de velocidades (m/s) normal al eje Z a 40 m …………………………………………………. 37 Fig 4.15 Perfil de velocidades (m/s) normal al eje Z a 10 m ……………………..………………..…..…… 37 ÍNDICE VI Fig 4.16 Perfil de velocidades (m/s) a lo largo de la dirección del viento.……………………..…….. 38 Fig A.I.1 Wind rose diagram…………………………………………………………………………………………………… 46 Fig A.I.2 Weibull diagram f(v)……………………………………………………………………………………...………… 46 Fig A.I.3 Power curve P(v)…………………………………………………….………………………………………..……… 47 Fig A.V.1 Mesh generated with blockMesh utility (cavity case)……………………………………………… 81 ÍNDICE VII Lista de Tablas Tabla 2.1 Modelos RANS suministrados por OpenFOAM ……………………………………………………… 11 Tabla 2.2 Coeficientes en las ecuaciones del transporte k-epsilon …………………………………….…. 13 Table 4.1 Características de la picea mariana, el pinus banksiana y aspen forest …………...……31 Tabla 4.2 Características del bosque utilizado en la simulación …………………………………..……….37 Table A.V.1 BlockMesh keywords ………………………………………………………………………………..…..….. 81 1. INTRODUCCIÓN 6 Simulación 1: Se utiliza un solver estándar creado por OpenFOAM (“simpleFoam”) de esta manera se obtienen datos validos sobre el flujo de aire en el caso de estudio. Simulación 2: La simulación se realiza con el nuevo solver (“treeFoam”), pero sin cargar las zonas arboladas previamente. Así se puede comparar la simulación 2 frente a la simulación 1, para comprobar que el nuevo solver al menos funciona correctamente. Simulación 3: Se realizan varias pruebas cargando previamente los bosques con la aplicación “leafIndexToFoam”, de esta manera se comprueba si el solver “treeFoam” tiene en cuenta el efecto de los árboles sobre el flujo de aire. Se varían las características del bosque para comprobar que el viento se ve distorsionado de manera correcta. 2. FUNDAMENTOS TEÓRICOS 7 2. Fundamentos Teóricos. A la hora de trabajar con herramientas de fluido dinámica computacional, se debe tener en cuenta que existen varios modelos numéricos que pueden otorgar una solución a un problema. La diferencia de estos modelos determinará su precisión, es decir su similitud con la situación real y también su complejidad a la hora de resolver el problema. Asimismo, las condiciones del flujo escogidas o asumidas también juegan un papel importante. Por ello se debe aceptar o elegir varias hipótesis:  Flujo incompresible. La densidad del fluido permanece contante, por lo que las variaciones de la densidad debidas a la altura y a la temperatura no se tienen en cuenta.  Comportamiento NEUTRAL del flujo en la “Capa límite atmosférica”. Se considera neutral debido a que la densidad del fluido es constante (ver APÉNDICE VII apartado “2.1 Atmospheric Boundary Layer”).  Efecto de fuerzas externas como la fuerza de Coriolis, o la rugosidad del terreno, se desprecian. Estas fuerzas externas pueden dificultar el análisis del efecto del follaje sobre la corriente de aire, por lo que es mejor añadirlas tras comprobar el buen funcionamiento del nuevo solver.  Modelo turbulento a utilizar. Se elige el modelo turbulento k-epsilon, por ser el mas común en las aplicaciones ingenieriles dentro de los métodos de resolución RANS (ver Apartado 2.3). 2.1. Turbulencia. 2.1.1. Efecto de la turbulencia en las Ecuaciones del Transporte promediadas. Primero de todo, se puede encontrar una breve definición de la turbulencia en el apéndice VII apartado 2.1.1 (“Turbulence” “Definition”). Calcular el comportamiento de los vórtices en una corriente turbulenta, no es una tarea sencilla; además no existe una capacidad computacional tal para poder resolver las ecuaciones 2. FUNDAMENTOS TEÓRICOS 8 de Navier-Stokes para flujos turbulentos, de manera completamente detallada, con dependencia del tiempo y con altos números de Reynolds. Normalmente el promedio de las propiedades del flujo, como la velocidad media o las tensiones medias, suministran suficiente información para los usuarios de CFD. No es necesario predecir el comportamiento de cada vórtice. Las propiedades como la velocidad (u) o la presión (p), se pueden establecer como la suma de una componente media estacionaria (U o P) y una componente fluctuante que varía con el tiempo (u´ o p´): La media de las componentes fluctuantes son cero por definición.  ∫ ( ) De la media cuadrática “root-mean-square” (rms) de las fluctuaciones, se puede obtener la componente de la velocidad fluctuante. √( )       [ ∫ ( ) ] ⁄ La media cuadrática de la velocidad se puede medir con una sonda (por ejemplo un “hot-wire” anemómetro) y un circuito eléctrico. Con urms la velocidad fluctuante (u´) se puede extraer. Introduciendo las componentes fluctuantes en las ecuaciones de Navier-Stokes, se puede estudiar el efecto de la turbulencia sobre el flujo (ver APÉNDICE II.1). ( ) ( ) ( ( )) 2. FUNDAMENTOS TEÓRICOS 9 Donde ρ es la densidad del fluid y SM representa las fuerzas externas. Utilizando coordenadas cartesianas: ( ) ( ) ( ) Las ecuaciones de Navier-Stokes promediadas son: ( ) ( ) ( ) [                ] ( ) ( ) [                ] ( ) ( ) [                 ] Estas son las llamadas ecuaciones de Reynolds. Aparecen seis términos nuevos después del proceso de promediado, estos términos son conocidos como “las tensiones de Reynolds”. Tres de ellos son tensiones normales y las otras tres son tensiones cortantes. Se necesita un modelo turbulento para calcular estas nuevas tensiones de manera suficientemente precisa, pero sin olvidar que una solución extremadamente precisa requiere una capacidad computacional no disponible. 2. FUNDAMENTOS TEÓRICOS 10 2.2. Modelado Numérico. 2.2.1. Modelos Turbulentos Resolver las fluctuaciones turbulentas de una manera muy detallada, es innecesario. Por eso, normalmente se buscan los efectos turbulentos que distorsionan el flujo medio. Un código CFD que sea útil, el modelo turbulento en el que se basa debe ser: suficientemente preciso, simple, informáticamente económico y sobre todo debe ser aplicable a un amplio número de casos. Existen varios métodos para estudiar el comportamiento de flujos turbulentos, a continuación se mencionan algunos de ellos.  Direct Numerical Simulations (DNS): Este modelo resuelve directamente las ecuaciones de Navier-Stokes, sin promedios ni aproximaciones. Es la manera mas precisa de resolver el comportamiento de un flujo turbulento, sin embargo conlleva un coste computacional tal que no se puede llevar a la práctica.  Large Eddy Simulation (LES): En los flujos turbulentos se pueden encontrar varias escalas de tiempo y espacio. Las escalas mas grandes tienen una capacidad de transporte mayor y albergan casi toda la energía. Este método se basa en resolver las grandes escalas del flujo y modelar las menores. Exige gran capacidad computacional y tiempo de cáculo.  Reynolds-Average Simulation (RAS) o Reynolds-Averaged Navier-Stokes (RANS): Es el modelo más utilizado para aplicaciones ingenieriles. Este modelo se basa en el promedio de las ecuaciones de Navier-Stokes, que da lugar a las ecuaciones de Reynolds.  Detached Eddy Simulation (DES): Este método es un híbrido entre RANS y LES. Se utiliza principalmente el método RANS pero para algunas regiones del flujo se usa la metodología LES. 2. FUNDAMENTOS TEÓRICOS 11 Como se ha dicho, el método de resolución RANS es el más común a la hora de resolver flujos turbulentos en el mundo de la ingeniería. Antes de que la turbulencia apareciera en el flujo existían cuatro incógnitas (u,v,w y p), sin embargo si se tiene en cuenta la turbulencia a través del promedio (time-averaging) de Reynolds, aparecen seis incógnitas nuevas. Para poder resolver el nuevo sistema de ecuaciones, formado por las ecuaciones de Reynolds (Eqn 2.7, Eqn 2.8.a, Eqn 2.8.b, Eqn 2.8.c), se requiere un modelo turbulento para cerrar el sistema. Dentro del método de resolución RANS existen varios modelos turbulentos, a continuación se mencionan aquellos que provee OpenFOAM para flujos incompresibles (Tabla 2.1). De todos estos, dos de ellos se utilizan comúnmente en aplicaciones ingenieriles, el modelo kepsilon standard y el modelo k-epsilon RNG. En este proyecto se utiliza el primero de ellos por lo que posteriormente, se desarrolla con mayor detalle. 2.2.2. Aproximación de Boussinesq Se puede suponer que las tensiones viscosas y las tensiones de Reynolds, están conectadas entre ellas y que ambas tienen casi el mismo efecto sobre el flujo medio. Por lo que las tensiones debidas a la turbulencia se pueden escribir como: ( ) Laminar Dummy turbulence model for laminar flow kEpsilon Standard k-ε model kOmegaSST k-ε-SST model RNGkEpsilon RNG k-ε model LaunderSharmaKE Launder-Sharma low-Re k-ε model LRR Launder-Reece-Rodi RSTM LaunderGibsonRSTM Launder-Gibson RSTM realizableKE Realizable k-ε model SpalartAllmaras Spalart-Allmaras 1-eqn mixing-length model Tabla 2.1 Modelos RANS suministrados por OpenFOAM 2. FUNDAMENTOS TEÓRICOS 12 De acuerdo con Boussinesq, las tensiones de Reynolds se pueden relacionar con la velocidad media de deformación, de tal manera que las tensiones turbulentas aumentan con el la velocidad media de deformación.       ( ) Donde μt es la viscosidad turbulenta (dimensiones Pa.s), se utiliza para modelar las fuerzas turbulentas como fuerzas viscosas. Si se obtiene μt, el sistema se puede resolver. La viscosidad total (μTotal) es la suma de la viscosidad del fluido (μ) y de la viscosidad turbulenta (μt) [2 y 3]. 2.2.3. Modelo turbulento k-epsilon. El modelo turbulento k-epsilon es un modelo de dos ecuaciones. Se introducen dos nuevos términos con sus respectivas ecuaciones de transporte, las cuales representan las propiedades turbulentas que aparecen en el flujo medio [2 y 3]. Estas dos nuevas variables son la energía cinética (k) y la disipación de energía cinética turbulenta (ε). La última de ellas determina la escala de la turbulencia mientras que la primera representa la energía en la turbulencia. La energía cinética turbulenta por unidad de masa se puede escribir como: (              ) Según el modelo de k-epsilon estándar (Launder y Spalding 1974), estas dos variables k y ε, se utilizan para definir las escalas de longitud (l) y velocidad (ν) características: ⁄ ⁄ La escala de velocidad y de longitud, definen la viscosidad turbulenta: Donde Cμ es una constante adimensional. 2. FUNDAMENTOS TEÓRICOS 13 Finalmente las ecuaciones del transporte que introduce el modelo k-epsilon son: Ecuación de la energía cinética: ( ) ( ) [( ⁄) ] Disipación de la energía cinética: ( ) ( ) [( ⁄) ] Las constantes que aparecen en las ecuaciones anteriores (Eqn 2.15 and 2.16) toman el siguiente valor: 2.3. Modelado de la fuerza de arrastre debida al follaje. 2.3.1. Parámetros característicos. Coeficiente de arrastre (CD). El comportamiento del flujo de aire a través de un bosque y la arquitectura del follaje de los árboles, están conectadas con el coeficiente de arrastre (CD). Para poder definir el mapa de árboles que se deben insertar en el domino es necesario calcular este parámetro. El problema radica en que no existe una tabla o lista donde se clasifique el coeficiente CD para cada tipo de árbol. Cμ Cε1 Cε2 σk σε 0,09 1,44 1,92 1,0 1,3 Tabla 2.2 Coeficientes en las ecuaciones del transporte k-epsilon. 2. FUNDAMENTOS TEÓRICOS 14 Hay varios estudios que intentan obtener ese parámetro mediante métodos experimentales, pero no hay ninguna ecuación que establezca una manera fiable para calcularlo (ver APÉNDICE II.2). Densidad de hojas por área (α). Este parámetro significa la superficie de follaje (hojas, ramas o tronco) por metro cúbico de bosque. Al contrario que CD se puede medir y no presenta problemas para su obtención. 2.4.2 Implementación de la fuerza de arrastre. A pesar de la dificultad para calcular el CD, la fuerza de arrastre que produce el follaje sobre el flujo de aire se puede añadir a la ecuación del momento como una fuerza externa de arrastre (Dalpe and Masson 2008 [7; 8]): | | El modelo turbulento elegido para realizar las simulaciones durante este estudio, es el modelo k – ε, el cual debe ser modificado para insertar la turbulencia producida por el follaje en la corriente de aire. Aparecen dos nuevos términos en las ecuaciones del transporte (SK y Sε). ( | | | |) ( | | | |) 3. OPENFOAM 15 3. OPENFOAM 3.1. Introducción. OpenFOAM es un software de CFD el cual es gratis y cuyo código esta abierto, fue desarrollado a finales de los años ochenta por la Imperial College de Londres; sin embargo no sería un software libre hasta el 2004. Al ser un software libre, se debe utilizar un sistema operativo también libre, como Linux, para poder utilizarlo. Como es un programa abierto, los usuarios pueden modificar, ampliar, modificar o simplemente personalizar las aplicaciones actuales de OpenFOAM; por este motivo es una herramienta CFD común entre la comunidad ingenieril. Se pueden resolver una amplia gama de problemas fluido dinámicos, flujos complejos con turbulencia, transferencia de calor e incluso reacciones químicas. También se puede utilizar para analizar la dinámica del sólido y sistemas electromagnéticos. Además de los “solvers”, OpenFOAM facilita herramientas de pre-procesado como por ejemplo, creadores de mallas (e.g. BlockMesh o snappyHexMesh) o conversores para adaptar mallas de otros softwares (e.g. FLUENT) a OpenFOAM. También incluye herramientas de posproceso, como paraView, que permiten trabajar con la solución del solver de manera gráfica. Existe una gran comunidad en la red, la cual comparte sus aplicaciones y ofrece ayuda entre los usuarios de OpenFOAM [14]. Las capacidades de los “solvers” estándares, las herramientas de pre y post procesado, así como las librerías de las que dispone OpenFOAM, están listadas en el apéndice VII apartado 3.1 “OPENFOAM Introduction”. Se puede decir que OpenFOAM es solo una librería, la cual permite crear archivos ejecutables (aplicaciones) para resolver problemas fluido dinámicos. Existen dos tipos de aplicaciones:  “Solvers”, los cuales resuelven los sistemas de ecuaciones que se aplican en cada caso.  Utilidades, que se utilizan para tratar los datos antes y después de ser resuelto un problema. 3. OPENFOAM 22 owner: En este archivo asocia una celda a cada cara (ver APÉNDICE IV.9). neighbour: Es una lista de las celdas vecinas por su etiqueta. blockMeshDict: Este archivo aparece si la malla se define utilizando “blockMesh”. Este archivo determina el procedimiento que la herramienta sigue para crear la malla (ver APÉNDICE V). La malla se puede definir a mano, lo cual es una tarea costosa si se quiere estudiar una situación levemente compleja; utilizando la utilidad “blockMesh” (ver APÉNDICE V) o la aplicación “snappyHexMesh” (ver [14] apartado 5.4 “Mesh generation with the snappyHexMesh utility”) también suministrada con OpenFOAM; o importando la malla desde otra herramienta de mallado como por ejemplo “GAMBIT”. 3.3.3. Fichero “0” Las variables que se utilizan en el caso de estudio se deben declarar en unos archivos, los cuales se encuentran en el interior del fichero “0”. La velocidad, presión, energía cinética o las nuevas variables necesarias para insertar el efecto de los bosques sobre el flujo de aire como el coeficiente de arrastre (CD) y la densidad de hojas por área (α). En estos archivos se fija el valor inicial para cada variable. El perfil de velocidades inicial se define en el archivo “U”; la presión en el archivo “p”, y así sucesivamente (ver APÉNDICE VI.1 y VI.2). Para las nuevas variables se crean dos archivos “Cd” y “alphaL” (ver APÉNDICE III.6) en los cuales los valores de las variables se escriben tras ejecutar la herramienta de pre-procesado “leafIndexToFoam” que carga los bosques al dominio (ver APÉNDICE VI.3). 4. PROCESO y RESULTADOS 23 4. Proceso y resultados 4.1. Caso de estudio, la colina Bolund. El caso Bolund es una colina real, la cual se ha estudiado en múltiples ocasiones mediante procesos experimentales y programas CFD. Ha sido utilizada para validar modelos del flujo sobre terreno complejo. Los efectos de cizalla y la turbulencia son muy importantes en el comportamiento del flujo de aire, y estos se acrecientan en terrenos complejos. Es por esta razón que las herramientas de CFD se están haciendo cada vez más imprescindibles cuando se quiere encontrar la mejor localización para instalar aerogeneradores. Los modelos que los programas de CFD suelen utilizar solo se comprueban con túneles de viento o sobre terreno simple. La malla que simula este caso, será la utilizada para testear las aplicaciones desarrolladas. La colina tiene 130 m de largo, 12 m de alto y 75 m de ancho, está situada en la costa, al norte de of Risø DTU (Dinamarca). Fig 4.1 Vista panorámica de la colina Bolund. 4. PROCESO y RESULTADOS 24 Esta colina es un buen desafío para cualquier solver de flujos en 3-D, ya que cambia radicalmente de rugosidad, del agua a la hierba, y además posee una gran pendiente entre el nivel del agua y la superficie de la colina que proporciona una cara vertical totalmente expuesta al flujo. Sin embargo como se puede apreciar en las imágenes (Fig 4.1 y Fig 4.2) no hay árboles sobre la colina, por lo que se probarán las aplicaciones mediante la creación de bosques virtuales de forma rectangular (para facilitar la prueba) sobre la colina. Por este motivo el proyecto necesitará de una comprobación empírica para asegurar que funciona correctamente al 100 %. Fig 4.2 La colina Bolund. Fig 4.3 Malla correspondiente a la colina Bolund. 4. PROCESO y RESULTADOS 25 4.2. Adición del follaje como variable en OpenFOAM. Como se ha dicho con anterioridad, el efecto del follaje de los árboles sobre el flujo de aire, se introducirá mediante el coeficiente de arrastre (CD) y la densidad de hojas por área (α). Ambos parámetros modifican las ecuaciones de Navier-Stokes que se utilizan para resolver el comportamiento del flujo. Este efecto aparece en forma de fuerza exterior (Su) en la ecuación del momento. El modelo turbulento k-epsilon también se verá afectado por estos parámetros, y los términos Sk y Sε (Eqn 2.19 y 2.20) aparecerán en las ecuaciones de transporte, en ambas, en la ecuación de la energía cinética y en la ecuación de disipación de la energía cinética (Eqn 2.15 y 2.16 respectivamente). Para que el simulador pueda tener en cuenta estos parámetros, primero deben ser declarados e introducidos al caso de estudio, para ello se ha desarrollado una herramienta de preprocesado que convierte una zona arbolada dada (en un archivo .map) en archivos que OpenFOAM puede leer y trabajar con ellos (archivos “Cd.H” y “alphaL.H” en el directorio “0”). Para que el solver pueda tener en cuenta las variables CD y α, primero de todo se deben modificar algunas librerías que el solver utilizará. Las librerías que incorpora OpenFOAM se encuentran en la carpeta “src” dentro del directorio “OpenFOAM”, por ello para evitar modificar el núcleo del programa, se crea un nuevo directorio “src” en el directorio de usuarios ($WM_PROJECT_USER_DIR), donde se crean las nuevas librerías. Las nuevas librerías también se deben compilar antes de correr el simulador, pero se usa el comando “wmake libso” [14] en lugar del utilizado para compilar las aplicaciones “wmake”. Los términos Sk y Sε se añaden modificando los archivos kEpsilon.C y kEpsilon.H los cuales definen el modelo turbulento kepsilon. CD y α son definidos en estos archivos también, al igual que todas las constantes necesarias que utilizan Sk y Sε. A continuación se muestra parte del código modificado en estos archivos para introducir las fuerzas externas debidas al efecto del paso del flujo de aire a través del follaje: 4. PROCESO y RESULTADOS 26 // Dissipation equation tmp<fvScalarMatrix> epsEqn ( fvm::ddt(epsilon_) + fvm::div(phi_, epsilon_) - fvm::Sp(fvc::div(phi_), epsilon_) - fvm::laplacian(DepsilonEff(), epsilon_) == C1_*G*epsilon_/k_ - fvm::Sp(C2_*epsilon_/k_, epsilon_) + fvm::Sp(Cd_*alphaL_/k_*(CEps4_*betap_*pow(mag(U_),3)- CEps5_*betad_*k_*mag(U_)), epsilon_); //Source Term) // Turbulent kinetic energy equation tmp<fvScalarMatrix> kEqn ( fvm::ddt(k_) + fvm::div(phi_, k_) - fvm::Sp(fvc::div(phi_), k_) - fvm::laplacian(DkEff(), k_) == G - fvm::Sp(epsilon_/k_, k_) + fvm::Sp(Cd_*alphaL_/k_*(betap_*pow(mag(U_),3)- betad_*k_*mag(U_)), k_); //Source term) // Momentum equation tmp<fvVectorMatrix> UEqn ( fvm::div(phi, U) + turbulence->divDevReff(U) + 2*(Omega ^ U) //coriolis force 2*(omega x U) + gravity //gravity force + fvm::Sp(Cdalpha*mag(U)*U,U) //canopy drag force == sources(U) 4. PROCESO y RESULTADOS 27 ); El nuevo “solver” desarrollado (“treeFoam”) requerirá las nuevas librerías cuando empiece la simulación. El “solver” incluye el nuevo término (Su) en la ecuación del momento, la cual se define en el archivo UEqn.H. Una vez se compile el código del “solver”, la aplicación será creada por el compilador en el directorio $FOAM_USER_APPBIN. El problema que surge al crear una nueva librería para la nueva aplicación, es que esas librerías están conectadas con una gran cantidad de otros archivos y librerías, por eso es necesario un buen programador para poder encajar las nuevas con todas las demás que también son necesarias. Finalmente tras muchos intentos infructuosos para crear estas, se optó por desbloquear los archivos de OpenFOAM y modificar kEpsilon.C y kEpsilon.H directamente en el núcleo del programa y añadir ahí las nuevas fuerzas de arrastre debidas al follaje y sus constantes. 4.3. Aplicación de pre-procesado. 4.3.1. Archivos “.map” Las características de las distintas zonas arboladas que se encuentran en la malla de estudio, se suministran a la simulación a través de un archivo con extensión “.map”. En este archivo se hayan: el coeficiente de arrastre (CD), la densidad de hojas por área (α) y la altura del follaje de los árboles (zo and ∆z). El área sobre el cual se extienden los distintos bosques también esta definida en este archivo, mediante una serie de coordenadas en el plano X-Y que definen polígonos que se corresponden con la situación de los árboles en el mapa. Las primeras cuatro líneas del archivo “.map” el código no las tiene en cuenta, en ese encabezado aparece el nombre del archivo, las dimensiones de CD y α, por lo que no es información relevante para la aplicación que convierte el archivo “.map” a un formato de trabajo para OpenFOAM. Tras esas líneas, se definen los polígonos y características del follaje que determinan un bosque. Se separan en grupos, cada grupo representa una zona boscosa y estos están formados por (ver APÉNDICE III.5):  CD, α, ∆z, zo y el número de vértices del polígono que define el área del bosque sobre el plano X-Y. Estos se hayan en la primera línea de cada grupo de datos. 4. PROCESO y RESULTADOS 28  Después del encabezado del grupo de datos, se declaran las coordenadas en 2-D de cada vértice, estas están ordenadas en sentido horario. Con todos los datos mencionados anteriormente, se define un volumen en el dominio, cuyas celdas contienen los valores de CD y α. En este caso se definirá un bosque con forma de prisma rectangular, cuya proyección en el plano X-Y será un simple rectángulo. 4.3.2. Aplicación de pre procesado “leafIndexToFoam” Esta aplicación convierte el archive “.map” a otros archivos con un formato que puede ser leído y cargado al solver de OpenFOAM (Cd.H and alphaL.H). Para conseguir este objetivo, el código implementado comprobará cada celda de la malla, y si el punto central de esta pertenece a un volumen de los definidos en el archivo “.map”, los valores de CD y α serán copiados a los campos internos de la celda. Cuando toda la malla se haya comprobado, los valores CD y α se escribirán en los archivos Cd.H y alphaL.H, los cuales se encuentran en el directorio “0” (ver Fig 3.3). A continuación se puede ver un esquema que explica el proceso seguido por la aplicación “leafIndexToFoam” para leer el archivo “.map” y escribir los archivos Cd.H y alphaL.H. 4. PROCESO y RESULTADOS 29 Fig 4.4 Proceso seguido para leer un archivo .map. 4. PROCESO y RESULTADOS 30 Para llevar a cabo el programa mostrado en la imagen anterior (Fig 4.4) se necesitan codificar en C++ los siguientes archivos:  leafIndexToFoam.C : Al compilar este archivo se creará el ejecutable que se usa para lanzar la aplicación, este archivo conecta todos los demás y los invoca según son necesarios. Es el código principal del programa.  createCdAlphaL.H : Como su nombre indica, crea los archivos donde se almacena los valores de CD y α de cada celda de la malla.  Check_inside_cell.H : En este archivo se codifica la parte del proceso donde se comprueba si el punto central de las celdas pertenecen al volumen definido en el archivo “.map”. Repite el proceso para todas las celdas de la malla.  Readmapfile_CdAlphaL.H : Abre el archivo “.map”, lee el polígono y los parámetros (CD,α,zo,∆z), tras ello invoca al código Check_inside.H. Repite el proceso tantas veces como polígonos son definidos. El código complete de cada archivo puede verse en el APÉNDICE III. Se pueden utilizar varios métodos para comprobar si un punto pertenece a un polígono. Se pueden ver algunos de ellos en la referencia [19], en internet. Sin embargo se exponen solo el método seleccionado a continuación. El método consiste en comparara las coordenadas X-Y de cada pareja de vértices consecutivos que forman el polígono, con las coordenadas X-Y del punto central de la celda. Existen dos grupos condiciones, cada grupo se compone a su vez de dos requisitos; el punto debe satisfacer alguna de para ser considerado perteneciente al polígono, en otro caso estará fuera. Condiciones A: ( ) ( ) ( ) ( ) Condiciones B: 4. PROCESO y RESULTADOS 31 ( ) ( ) ( ) ( ) Este método se puede utilizar para comprobar si un punto pertenece a un polígono en 2-D, por eso en principio es insuficiente para localizar un punto en un volumen. Sin embargo, una vez se comprueba que la proyección del punto central de la celda sobre el plano X-Y pertenece al polígono, únicamente se debe comprobar si la coordenada Z del punto se encuentra entre zo y (zo + ∆z). Si se cumplen ambas condiciones el punto estará alojado en el interior del volumen. El código completo que comprueba la existencia de un punto en un volumen puede verse en el APÉNICE III.4 (Check_inside_cell.H). 4.3.3. Simulación de zonas arboladas. Un archivo “.map” se ha creado para poder llevar a cabo la simulación de la aplicación “leafIndexToFoam”. Se desea comprobar si el programa carga correctamente los datos y características de los bosques definidos en el archivo “.map” al dominio en estudio. A continuación se muestran las características de algunas especies de árboles más comunes en Europa (Tabla 4.1). ∫ Donde LAI es el índice de hojas por área (LAI = Leaf Area Index) y α es la densidad de hojas por área (m2/m3). La altura a la cual el follaje comienza es (zo), la cual se asume de 3 m para estas especies árboles. Como la superficie de la colina Bolund está a 12 m sobre el nivel del mar, zo = 15 m. Se Picea Mariana Pinus Banksiana Aspen forests ∆z(m) 10 15 10 LAI 9,19 2,22 3,57 CD 0,15 0,45 0,4 α 0,919 0,148 0,357 Tabla 4.1 Características de la picea mariana, el pinus banksiana y aspen forest. [8] 4. PROCESO y RESULTADOS 38 Primero de todo, se puede apreciar una gran diferencia entre la solución obtenida por “simpleFoam” o por “treeFoam” sin zonas boscosas, frente a cualquiera de los dos perfiles de velocidades que ofrecen las simulaciones con “treeFoam” añadiendo bosques. Por lo que el nuevo solver tiene en cuenta el efecto que produce el follaje cuando el viento lo atraviesa. La segunda conclusión que se puede obtener de las dos últimas simulaciones, es que si se observan las figuras anteriores (Fig 4.15, 4.16, 4.17, 4.18), se puede afirmar que para valores mayores de α, el flujo de aire se ve distorsionado en mayor medida cerca del bosque. Por eso se declarar que se ha testeado “treeFoam” y funciona correctamente. No obstante para que pueda ser 100 % comprobado, se deben enfrentar simulaciones informáticas con CFD frente a los datos empíricos obtenidos en un experimento. α =0,919 α =1,838 Fig 4.16 Perfil de velocidades (m/s) a lo largo de la dirección del viento. 5. CONCLUSIÓN 39 5. Conclusión. Para empezar, este proyecto tenía una dificultad añadida, el software empleado, OpenFOAM, no es uno de los programas habituales que se estudian en las aulas, principalmente porque su manejo no es tan sencillo como el de otros programas CFD y trabaja sobre el sistema operativo Linux. Por lo que este proyecto y los que de él deriven, requieren de conocimientos básicos con Linux, con el lenguaje de programación C++ y especialmente con OpenFOAM. Sin embargo, como OpenFOAM es un software libre, es un programa CFD perfecto para desarrollar nuevos “solvers”, sin incurrir en los gastos de adquirir un software comercial. Por otro lado ya que el código fuente está abierto, los usuarios pueden modificar las aplicaciones existentes para añadir mejoras o simplemente adaptarlos a sus casos particulares más específicos. Por último, el trabajo de cada usuario se puede compartir con la comunidad online para que lo prueben e incluso aconsejen o sugieran algunas mejoras o correcciones. Se han desarrollado dos aplicaciones en este proyecto. La primera de ellas (“leafIndexToFoam”) es una aplicación de pre procesado la cual carga las zonas arboladas al domino de estudio. Los valores de α y CD se leen del archivo “.map” y se escriben satisfactoriamente en los archivos alphaL.H y Cd.H; gracias a la herramienta paraView se pueden observar estos valores de manera gráfica para cada celda de la malla, como se muestra en las figuras (Fig 4.7, 4.8, 4.9 y 4.10). Con este visualizador se puede ver como el bosque cargado al dominio es fiel al descrito en el archivo “.map”. La segunda aplicación es “treeFoam”, este es un “solver” el cual resuelve el sistema de ecuaciones de Reynolds, aplicando un modelo turbulento. En este caso se añadió al “solver” las ecuaciones pertinentes para simular el efecto producido en el flujo de aire cuando atraviesa el follaje de los árboles. Las imágenes (Fig 4.15, 4.16, 4.17 y 4.18) muestran este efecto, el perfil de velocidades varía de manera apreciable cuando el viento se topa con un conjunto de árboles. Se puede observar también que cuanto mas frondoso es el bosque (el coeficiente α es mayor) el perfil de velocidades se distorsiona en mayor medida. Para poder presentar “treeFoam” como un “solver” 100% fiable, se debería realizar un ensayo experimental y comparar los resultados obtenidos con la simulación del ensayo. Si los resultados son semejantes se podría confirmar como un “solver” CFD totalmente válido. Sin embargo esta tarea no es sencilla porque incluso la obtención del coeficiente de arrastre (CD) 5. CONCLUSIÓN 40 no esta claramente definida por el momento para cada especie de árbol, por lo que suele ser estimada. 6. LÍNEAS FUTURAS 41 6. Líneas futuras. El objetivo mas importante que se propone para seguir mejorando este “solver”, es automatizar la ejecución de una simulación. De esta manera se permitirá al usuario realizar varias simulaciones con diferentes perfiles de flujo de aire, cambiando su velocidad y dirección, sin tener que preparar todos los archivos necesarios para simular un caso en OpenFOAM. El software incorpora una aplicación “python” facilita esta automatización, este sería un buen punto de partida. En segundo lugar se sugiere añadir también el efecto que produce un aerogenerador sobre el flujo de aire, ya que normalmente en una granja eólica se hallan varios de estos instalados y cada uno distorsiona el flujo de aire. Por último se debería cambiar el proceso del estado estacionario al modelo transitorio, para poder realizar simulaciones más fieles a las situaciones reales. 7. BIBLIOGRAFÍA 43 7. Bibliografía. [1] International Energy Agency (IEA). [2] H.K. Versteeg and W. Malalasekera. An introduction to Computational Fluid Dynamics. The Finite Volume Method. Longman Group Ltd, 1995. [3] Lars Davidson. An introduction to turbulence Models Göteborg, Sweden: Chalmers University of Technology, 2011. [4] D. Poggi, A. Porporato, C. Ridolfi, J.D. Albertson and G.G. Katul. The effect of vegetation density on canopy sub-layer turbulence. Kluwer Academic Publishers, 2004. [5] B. Marcolla, A.Pitacco and A.Cescatti. Canopy architecture and turbulence structure in a coniferous forest. Kluwer Academic Publishers, 2003. [6] A.Cescatti and B.Marcolla. “Drag coefficient and turbulence intensity in conifer canopies.” Agricultural and Forest Meteorology, vol. 121, pp. 197-206, 2004. [7] B.Dalpé and C.Masson. “Numerical simulation of wind flow near a forest edge.” Journal of Wind Engineering and Industrial Aerodynamics, vol. 97, pp. 228-241, 2009. [8] B.Dalpé and C.Masson. “Numerical Study of Fully Developed Turbulent Flow Within and Above a Dense Forest.” Wind Energy, vol.11, pp. 503-515, 2008. [9] Xabier Pedruelo. Modeling of Wind Flow over Complex Terrain Using OpenFOAM. Vattenfal Research and Development AB, 2009. [10] Eric Furbo. Evaluation of RANS turbulence models for flow problems with significant impact of boundary layers. Uppsala University, 2010. [11] Benjamin Martinez. Wind resource in complex terrain with OpenFOAM. Risø DTU, National Laboratory for Sustainable Energy, 2011. [12] Håkan Nilsson. Basic of C++ and object orientation in OpenFOAM. Chalmers University. [13] Jan Skansholm C++ From the Beguinning. Addison Wesley Longman, 1977. [14] OpenFOAM User´s Guide. OpenFOAM Foundation, 2011. 7. BIBLIOGRAFÍA 44 [15] OpenFOAM Programmer´s Guide. OpenFOAM Foundation, 2011 [16] El Ideal. “Reserva ‘ricotí ‘, un territorio protegido para la Alondra de Dupont”. Internet: http://www.ideal.es/granada/20120708/local/granada/reserva-ricoti-territorio-protegido- 201207081726.html [17] Godfrey Boyle. Renewable Energy: Power for a Sustainable Future.Oxford. [18] Wikipedia. “Wind rose”. Internet: http://en.wikipedia.org/wiki/Wind_rose [19] Paul Bourke. “Determining if a point lies on the interior of a polygon”. Internet: http://local.wasp.uwa.edu.au/~pbourke/geometry/insidepoly/ [20] OpenFOAM programming tutorial. Politecnico de Milano, Chalmers. [21] Wind Energy THE FACTS. “The Annual Variability of Wind Speed”. Internet: http://www.wind-energy-the-facts.org/en/part-i-technology/chapter-2-wind-resource- estimation/local-wind-resource-assessment-and-energy-analysis/the-annual-variability-of- wind-speed.html [22] WIND POWER PROGRAM. “Wind Turbine power output variation with steady wind speed.” Internet: http://www.wind-power-program.com/turbine_characteristics.htm