Aplicación del método de los volúmenes finitos al cálculo del parámetro de masa añadida de un cilindro vibrante en el seno de un fluido

Trabajo Fin de M´aster Aplicaci´ on del M´ etodo de los Vol´ umenes Finitos al c´ alculo del par´ ametro de masa a˜ nadida de un cilindro vibrante en el seno de un fluido Las Palmas de Gran Canaria, Julio de 2015 Alumno D. H´ ector Rub´ en D´ ıaz Ojeda Tutor Dr. D. Luis Alberto Padr´ on Hern´ andez Universidad de Las Palmas de Gran Canaria Espa˜na SIAN I 2014-2015 M´ aster Universitario en Sistemas Inteligentes y Aplicaciones Num´ ericas en Ingenier´ ıa Agradecimientos No quería dejar pasar esta oportunidad para recordar y agradecer a todas aquellas personas que han hecho posible que haya llegado hasta aquí y que me instarán a continuar. Entre ellas, cabe nombrar a mis padres y hermana por el apoyo y sacrificio que realizan diariamente para posibilitarme esta valiosa formación que he estado adquiriendo durante estos años y por ser fuente de apoyo tanto en las buenas como en las malas. A mi tutor, D. Luis Alberto Padrón Hernández, primero por a haberme permitido trabajar junto a él y segundo por la formación, consejos y ayuda que desde el sosiego me ha transmitido siendo de gran valor para mí. Finalmente me gustaría agradecer a los compañeros del Máster SIANI por su comportamiento, amabilidad y unión durante este año. Aunque, especialmente señalar a D. Jorge Ignacio Pirez Mora por el gran apoyo y por la motivación que me transmite diariamente llegando a ser una de las personas que más apoyo me han dado para acabar este trabajo. Y no dejando escapar la oportunidad de agradecer enormemente también a los compañeros Ingenieros Navales Pedro Gallego Martínez y Tanausú Almeida Medina por su gran interés, ayuda y ganas de progresar que hacen que mi interés por continuar trabajando duro aumente. Gracias. I Resumen Enmarcado dentro del campo de la Mecánica de Fluidos y desarrollado en conjunto con la división de Mecánica de los Medios Continuos y Estructuras del Instituto Universitario SIANI, se desarrolla un trabajo fin de título correspondiente al máster en Sistemas Inteligentes y Aplicaciones Numéricas en Ingeniería. En este trabajo académico se propone estudiar la interacción fluido-estructura para el caso de estructuras infinitamente rígidas (2D), vibrando en contacto con un flujo viscoso e incompresible. Para la discretización y resolución de las ecuaciones Navier–Stokes se hará uso del Método de Volúmenes Finitos que se encuentra implementado en la herramienta de código abierto OpenFoam. Más concretamente se realizará la comparativa del coeficiente de masa añadida y del coeficiente de arrastre para una estructura cilíndrica que vibra según un movimiento armónico en un fluido inicialmente en reposo con los resultados de problemas clásicos de los que existen diversas publicaciones. No obstante, para poder desarrollar y obtener unos resultados válidos, se realizará un tratamiento adecuado de la malla de estudio. En conclusión, con los resultados que se obtienen en este trabajo, se estima conseguir la base metodológica para posteriormente pasar a calcular no solo geometrías más complicadas sino también geometrías en 3D. III Abstract Belonging to the Fluid Mechanic field and developed together with the Continuum Mechanics and Structures group at the University Institute SIANI, a Master Thesis is developed corresponding to the Intelligent Systems and Numerical Applications in Engineering master. In this academic work, we do intend to study the fluid – structure interaction for indefinitely rigid structures (2D), vibrating in contact with a non-compressible and viscous flow. For the discretization and the solution of Navier–Stokes′equations we will make use of the Finite Volume Method which is implemented in the open source tool OpenFoam. Concretely, a comparative will be done between the added mass coefficient and drag coefficient for a cylinder which will be oscillating in vacuum under a harmonic movement with classicproblemsaboutwhichseveral research papers havebeen published. Nonetheless, in order to develop and get reliable results, we had better make an accuracy process to the mesh we are going to deal with. To sum up, with the results that are supposed to be getting, it is deemed that we have already got the methodological basis to subsequently go in not only with more complicated geometries but also with 3D geometries. V Índice general Lista de figuras 5 Lista de tablas 7 1. Introducción 9 1.1. Motivación .................................. 9 1.2. Antecedentes ................................. 10 2. Conceptos básicos de Dinámica de Fluidos 13 2.1. Introducción.................................. 13 2.2. Las ecuaciones de Navier - Stokes . . . . . . . . . . . . . . . . . . . . . . 14 2.2.1. Ecuación de conservación de masa . . . . . . . . . . . . . . . . . . 14 2.2.2. Ecuación de la conservación de momento . . . . . . . . . . . . . . 15 2.3. Modelos de flujo simplificados . . . . . . . . . . . . . . . . . . . . . . . . 15 2.3.1. Flujo incompresible . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.3.2. Flujonoviscoso ........................... 17 2.3.3. Flujo potencial y flujo ideal . . . . . . . . . . . . . . . . . . . . . 17 2.4. Turbulencia.................................. 18 2.4.1. Aproximaciones numéricas para el tratamiento de la turbulencia . . 19 2.4.2. Modelok–épsilon........................... 20 1 2.5. Resolución y discretización de las ecuaciones de la dinámica de fluidos (CFD) 20 3. Metodología 23 3.1. Introducción.................................. 23 3.2. El Método de los Volúmenes Finitos . . . . . . . . . . . . . . . . . . . . . 23 3.2.1. Conceptos generales . . . . . . . . . . . . . . . . . . . . . . . . . 23 3.2.2. Discretización espacial . . . . . . . . . . . . . . . . . . . . . . . . 25 3.2.3. Definiciones generales de la metodología numérica . . . . . . . . . 26 3.3. Herramienta numérica: OpenFoam . . . . . . . . . . . . . . . . . . . . . . 27 3.3.1. Definición de OpenFoam . . . . . . . . . . . . . . . . . . . . . . . 27 3.3.2. Capacidades de OpenFoam . . . . . . . . . . . . . . . . . . . . . . 28 3.3.3. Casos................................. 29 3.4. Simulación en OpenFoam . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 3.4.1. Generación de mallas mediante BlockMesh . . . . . . . . . . . . . 30 3.4.2. Elección del solver . . . . . . . . . . . . . . . . . . . . . . . . . . 41 3.4.3. Post-proceso ............................. 42 3.5. Concepto de masa añadida . . . . . . . . . . . . . . . . . . . . . . . . . . 43 3.6. EcuacióndeMorison............................. 44 3.6.1. MétododeMorison.......................... 45 3.6.2. Cálculo del ángulo 𝛿mediante la transformada rápida de Fourier . . 46 4. Resultados 49 4.1. Introducción.................................. 49 4.2. Solución numérica del problema del cilindro oscilante . . . . . . . . . . . . 49 4.2.1. Estructura de la malla utilizada para los estúdios numéricos . . . . 50 4.2.2. Condiciones de contorno . . . . . . . . . . . . . . . . . . . . . . . 51 4.2.3. Mallados ............................... 52 2 Capítulo 1 Introducción 1.1. Motivación Dentro del área de la mecánica de fluidos, surgen numerosos problemas como la tarea de estimar las cargas tanto de corrientes como de oleaje en el diseño de las plataformas petrolíferas (Figura 1.1), molinos eólicos Offshore y otras estructuras relacionadas. Esto, en otras palabras, tiene como partida la determinación del coeficiente de masa añadida y del coeficiente de arrastre siendo esto en lo que se centrará el presente trabajo. Dicho trabajo se ha dividido en seis capítulos que se exponen a continuación: Introducción. Conceptos básicos de Dinámica de Fluidos. Metodología. Resultados. Conclusiones. Referencias. En los capítulos mencionados se comentará lo relativo al problema que se ha estudiado. Dicho problema consiste en la interacción fluido-estructura en el caso de estructuras infinitamente rígidas vibrando en contacto con un flujo viscoso e incompresible, donde validar una metodología basada en el Método de los Volúmenes Finitos para el cálculo de los parámetros 9 10 1. INTRODUCCIÓN masa añadida y coeficiente de arrastre para el caso bidimensional será el principal objetivo. Esta validación se realizará a través de la comparación de resultados con un problema clásico para el que existen diversos resultados publicados y haciendo uso de la herramienta de código abierto OpenFoam. Esto, es de suma importancia ya que se deberá ser capaz de diseñar estructuras que puedan operar en zonas con condiciones adversas y que han sido de acceso imposible en las últimas décadas, además de tener la intención de aumentar la vida útil de los materiales no solo para disminuir el aspecto económico sino también para aumentar la seguridad de operación de dichas estructuras. Siendo aquí donde los parámetros comentados entran en juego ya que son entre otros los que permiten modelar la influencia del fluido en contacto con una estructura en vibración. La consecución del objetivo arriba expuesto permitirá, en un paso posterior, utilizar la misma metodología para el cálculo de estos parámetros, primero en problemas de carácter tridimensional con geometrías sencillas (por ejemplo de cilindros de longitud reducida), y posteriormente de geometrías complicadas (como pueda ser el caso de estructuras tipo Jacket para estructuras Offshore). 1.2. Antecedentes Indudablemente el uso y la importancia del CFD viene dado principalmente por el gran desarrollo que la computación ha experimentado desde 1950. Esto es así ya que, hasta que este desarrollo no se produjo, era inviable el uso de estas técnicas que requieren el cálculo de gran número de operaciones por unidad de tiempo. Estas técnicas permiten que se lleve a cabo en este trabajo final de máster necesario para la obtención del título de máster en Sistemas Inteligentes y Aplicaciones Numéricas en la Ingeniería una serie de simulaciones que indicarán los valores de los coeficientes comentados en el apartado anterior. Por otro lado, en este trabajo se ha seguido una extensa y variada bibliografía que va desde textos elementales de Mecánica de Fluidos hasta textos más específicos relacionados conla materiaquenos involucra (Massieand Journée[8], White[17].)Todo estopasando por libros como los de Ferziger and Perić [3], Versteeg and Malalasekera [11] o Fernandez Oro [13] que hacen un recorrido desde los conceptos de la mecánica de fluidos hasta los métodos numéricos aplicados para resolver diferentes problemas como el que se ha propuesto en este trabajo fin de máster. Sin lugar a dudas debemos destacar tres textos base en lo relacionado con la aplicación que será desarrollada en este trabajo: Chen [2], Wilson [4] y Sarpkaya [15]. Estos son de alguna manera los textos más conocidos y más referenciados por diversas publicaciones científicas, ya que básicamente contienen todo lo elemental sobre aplicaciones y problemas 1.2. ANTECEDENTES 11 Figura 1.1: Ejemplo real de estructuras cilindricas en plataforma Offshore. específicos del sector Offshore. Finalmente, es de destacar que existe una amplia literatura e investigaciones sobre el estudio de la interacción fluido–estructura (Chaplin [1], D𝑢tsch et al. [6], Gadelho et.al. [7], Bonfiglio et.al. [10] y Keulegan and Carpenter [9]), artículos que se centran en el cálculo de los coeficientes hidródinamicos para diversas consideraciones de la estructura como por ejemplo para una estructura oscilatoria propuesta por Bonfiglio et.al. [10]. Aunque nuestro trabajo tendrá una estrecha relación con Sarpkaya [14] y sobre todo con Troesch and Kim [16]. Esta última publicación muestra unos resultados experimentales obtenidos para 𝐶𝑚y 𝐶𝑑los cuales serán de gran interés a la hora de validar la metodología que aquí se expone. Capítulo 2 Conceptos básicos de Dinámica de Fluidos 2.1. Introducción A la hora de analizar los problemas de mecánica, existen dos posibilidades: 1. El enfoque lagrangiano donde la solución buscada es la evolución temporal de las magnitudes de todas y cada una de las partículas de sistema que intervienen en el problema. Matemáticamente las magnitudes las consideraremos funciones de la posición inicial de la partícula y del tiempo (usado en mecánica de sólidos). 2. El enfoque euleriano en el que se calculan los valores las magnitudes de todas las partículas del sistema que en cada instante 𝑡están ocupando un volumen de control fijado previamente. Matemáticamente las magnitudes que se analizan son funciones de la posición en el volumen de control y en el tiempo (usado en mecánica de fluidos). En este trabajo, se hace uso de la segunda posibilidad. Es por ello, que en este capítulo se introducirán de una forma somera las ecuaciones pertinentes de los fluidos con la única finalidad de poner en contexto y dar la base elemental al lector para que se pueda lograr un mayor entendimiento del trabajo que será desarrollado posteriormente. 13 14 2. CONCEPTOS BÁSICOS DE DINÁMICA DE FLUIDOS 2.2. Las ecuaciones de Navier - Stokes Si consideramos una variable específica 𝜙definida sobre un volumen de control (VC) y calculamos la variación temporal de la variable en dicho volumen de control, teniendo en cuenta que existen dos mecanismos fundamentales responsables de la generación de un flujo, como son la difusión, originado a nivel molecular, microscópico, y la convección, asociado al movimiento del fluido a nivel macroscópico, es posible llegar a la ecuación 2.1: 𝜕(𝜌𝜙) 𝜕𝑡 ⏟⏟⏟ 𝑇 𝑒𝑚𝑝𝑜𝑟𝑎𝑙 + ∇ ⋅(𝜌𝑣𝜙) ⏟⏞⏞⏟⏞⏞⏟ 𝐶𝑜𝑛𝑣𝑒𝑐𝑡𝑖𝑣𝑜 = ∇ ⋅(Γ∇𝜙) ⏟⏞⏞⏟⏞⏞⏟ 𝐷𝑖𝑓𝑢𝑠𝑖𝑣𝑜 +𝑆 ⏟⏟⏟ 𝐹𝑢𝑒𝑛𝑡𝑒 (2.1) Por tanto y como se define en la ecuación 2.1 aparecen cuatro términos: Término temporal, representa la variación local con el tiempo en el interior del volumen de control; es decir, la acumulación o disminución de 𝜙. Término convectivo, que representa el transporte de la variable de un punto a otro del dominio por medio de la velocidad de flujo. Término difusivo, se corresponde con alguno de los fenómenos de transporte que ocurren a nivel molecular: la ley de Fourier para la difusión de calor; la ley de Fick para la difusión de masa o la ley de Newton para la difusión de cantidad de movimiento por efectos viscosos. Término fuente, para tener en cuenta fuentes de generación o destrucción de la variable transportada. Finalmente cabe señalar que las ecuaciones de gobierno para el flujo y la transferencia de calor y de masa, así como aquellas que rigen el transporte de otras variables se pueden enunciar de forma conservativa como en la ecuación 2.1. En cada caso, los coeficientes de transporte, así como la variable específica que sustituye a 𝜙, adoptan formas diferentes. 2.2.1. Ecuación de conservación de masa La ecuación de conservación de masa, o ecuación de continuidad, establece de forma general que el incremento de masa en el interior de un elemento fluido es consecuencia del flujo neto de masa hacia dicho elemento. Puesto que en general no puede crearse ni 2.3. MODELOS DE FLUJO SIMPLIFICADOS 15 destruirse masa (excepto en casos donde intervengan reacciones nucleares), la expresión no estacionaria, tridimensional, en un punto para un fluido general es: 𝜕𝜌 𝜕𝑡 + ∇ ⋅(𝜌𝑣)=0 (2.2) En el caso particular de flujo incompresible la densidad del fluido no varía temporal ni espacialmente en el dominio, por lo que la ecuación se reduce a que la divergencia de la velocidad debe ser nula. ∇⋅𝑣 = 0 (2.3) 2.2.2. Ecuación de la conservación de momento La expresión general de la ecuación de movimiento para un fluido newtoniano (por tanto, que presenta una relación lineal entre las tensiones cortantes y deformaciones resultantes según 𝜏=𝜇𝜕𝑢 𝜕𝑦 ) e incompresible es (White [17]): 𝜌𝜕𝑣 𝜕𝑡 +𝜌(𝑣 ⋅∇)𝑣 = −∇𝑝+𝜌𝑔 + ∇(𝜇∇𝑣)(2.4) Aprovechándose de las propiedades de los campos vectoriales, es relativamente fácil reformular el término convectivo según: (𝑣 ⋅∇)𝑣 ≡𝑣 ⋅∇𝑣 = ∇ ⋅(𝑣 ⋅𝑣) − 𝑣 ⋅(∇ ⋅𝑣)(2.5) Y teniendo en cuenta que al ser incompresible el flujo es adivergente , se obtiene directamente sustituyendo 2.5 en 2.4: que es la ecuación en su forma conservativa. 𝜕(𝜌𝑣) 𝜕𝑡 + ∇(𝜌𝑣𝑣) = ∇(𝜇∇𝑣)−∇𝑝+𝜌𝑔 (2.6) Se obtendría una expresión similar para el caso compresible pero para el caso que se estudia aquí no procede enunciarla. 2.3. Modelos de flujo simplificados Aunque no sea objeto de este trabajo analizar todas las posibles simplificaciones de las ecuaciones de gobierno, sí es conveniente repasarlas brevemente con el objeto de presentar las técnicas numéricas más apropiadas para resolver cada una de ellas. En primer lugar, se puede decir que el comportamiento dinámico de fluidos está gobernado por las ecuaciones de conservación de la masa, la cantidad de movimiento y la energía. 16 2. CONCEPTOS BÁSICOS DE DINÁMICA DE FLUIDOS En el caso más general, las ecuaciones completas de Navier-Stokes presentan las siguientes características: Son dependientes del tiempo. Tridimensionales. Incluye esfuerzos viscosos. Considera la compresibilidad del fluido. Son ecuaciones en derivadas parciales elípticas y altamente no lineales. Se encuentran acopladas entre sí. Sin embargo, en la mayoría de los casos es posible reducir el análisis, simplificando las ecuaciones mediante: •Imponiendo unas condiciones de contorno que no varíen con el tiempo (estacionarios). •La reducción de dimensiones (problemas estacionarios, bidimensionales o axisimétricos). •Aproximando el comportamiento del fluido (incompresible, no viscoso). •Simplificando las fuerzas que actúan sobre el fluido (hidrostático, potencial). •Aproximaciones basadas en valores medios (Reynolds medio, profundidad media (depth-averaged). 2.3.1. Flujo incompresible Los fluidos, en general, son todos compresibles en mayor o menor grado. Sin embargo, las variaciones de densidad motivadas por cambios de presión pueden ser despreciados si: El número de Mach, 𝑀𝑎=𝑢 𝑐<1. Donde: •u = velocidad del fluido. •c = velocidad del sonido. Las variaciones de temperatura son pequeñas. Cabe resaltar que el término incompresible, no significa que la densidad sea constante; de hecho, existen un gran número de casos donde la densidad es la responsable del movimiento de los fluidos; como por ejemplo: las variaciones de salinidad en los océanos o las variaciones de temperatura en la atmósfera. Las consecuencias numéricas de los flujos incompresibles son: 2.3. MODELOS DE FLUJO SIMPLIFICADOS 17 La energía interna es irrelevante (las ecuaciones no están acopladas). La conservación de la masa lleva implícita una ecuación para la presión. Solución mediante un método basado en la presión. Con lo dicho con anterioridad, podremos concluir diciendo que en condiciones de flujo isotermo, la hipótesis de incompresibilidad implica que la solución de las variables primitivas de flujo (presión y velocidad) se haga independiente del campo de temperaturas. El sistema de ecuaciones requerido queda reducido a la ecuación de continuidad adivergente y a la ecuación de cantidad de movimiento, que expresadas vectorialmente establecen: ∇⋅𝑣= 0 𝜌𝜕𝑣 𝜕𝑡 +𝜌(𝑣 ⋅∇)𝑣 =𝜌 𝑓− ∇𝑝+𝜇∇2𝑣 (2.7) Nota: la ecuación 2.7 corresponde a un flujo incompresible viscoso. 2.3.2. Flujo no viscoso Todos los fluidos reales poseen viscosidad. Sin embargo, muchas veces se comportan como si no la tuvieran. Por esta razón, resulta útil conocer la dinámica de un fluido ideal que se supone incompresible y carente de viscosidad. El análisis del movimiento de un fluido ideal es más simple que el de un flujo viscoso debido a que no existen esfuerzos cortantes. Los esfuerzos normales son los únicos que se deben considerar en dicho análisis. Esto deja la ecuación como la ecuación 2.8 como la correspondiente para este caso y que se conoce como la ecuación de Euler: 𝜌𝑔 − ∇𝑝=𝜌𝑑 𝑣 𝑑𝑡 (2.8) 2.3.3. Flujo potencial y flujo ideal En el caso de flujos con números de Reynolds suficientemente altos, despreciar los efectos viscosos y de conducción resulta una aproximación bastante cómoda, pues elimina los términos difusivos de segundo orden en las ecuaciones, convirtiéndolas en ecuaciones de primer orden. 18 2. CONCEPTOS BÁSICOS DE DINÁMICA DE FLUIDOS Manteniendo la condición de flujo subsónico, al igualar la viscosidad a cero, se obtienen las ecuaciones de Euler, en las que todo el transporte de cantidad de movimiento se debe a fenómenos macroscópicos. Esto es: ∇⋅𝑣= 0 𝜌(𝜕𝑣 𝜕𝑡 + (𝑣 ⋅∇)𝑣) = 𝜌 𝑓− ∇𝑝(2.9) En el flujo estacionario, se cumple que: ∇(𝑣2 2+𝑝 𝜌+𝑔𝑧)=0 (2.10) por lo que H no varía a lo largo de una línea de corriente. Si el flujo es irrotacional se cumple además que ∇x𝑣 = 0, por lo que H se hace constante en todo el dominio del flujo. En estas condiciones, denominadas flujo potencial, se puede definir un potencial de velocidad, 𝑣 =∇𝜙, que introducido en la ecuación de continuidad proporciona la ecuación de Laplace. Esta ecuación se resuelve con la condición de componente normal nula en la superficie de un cuerpo sumergido y con la especificación de la velocidad en los puntos lejanos al cuerpo. Una vez obtenido el campo de velocidades, la ecuación de Euler proporciona el campo de presiones. La ecuación de Laplace es lineal y posee soluciones exactas simples que pueden superponerse para crear nuevas soluciones. Cuando la velocidad no se puede expresar a partir de un potencial, el flujo ideal se suele resolver introduciendo la vorticidad (∇x𝑣), que ya no es nula, en las ecuaciones de cantidad de movimiento y energía. 2.4. Turbulencia La turbulencia es un estado caótico e irregular del movimiento de un fluido que se establece a partir de la aparición de irregularidades en las condiciones iniciales o de contorno de la corriente fluida. Estas inestabilidades se amplifican y se retroalimentan de forma cíclica, creando vórtices (eddies) turbulentos que se crean y se destruyen. La turbulencia es una característica de los flujos, no de los fluidos como tales, sino que su aparición exige de la existencia de un fluido en movimiento, en el que los fenómenos de convección (inerciales) asociados a la velocidad sean varios órdenes de magnitud superiores a los efectos difusivos (disipativos)relacionados con la viscosidad del fluido. Esta relación es 3.2. ELMÉTODO DE LOS VOLÚMENES FINITOS 25 Figura 3.2: Ejemplo de malla bidimensional.[5] Este esquema se puede aplicar antes o después de integrar. Si la integral es factible de hacer entonces primero se hace la integral, y luego se aplica el esquema de diferencias finitas. Si la integral no es resoluble, se hace el esquema de diferencias fintas dentro del integrando y luego se integra. 4. Construcción de una matriz global. El procedimiento anterior conduce a un sistema de ecuaciones, que una vez resuelto proporciona la solución buscada. 3.2.2. Discretización espacial Básicamente, existen dos tipos de mallado: 1. Mallados estructurados: cada punto de la malla está inequívocamente identificado por 26 3. METODOLOGÍA los índices i, j, k, en coordenadas cartesianas. Las celdas de la malla son cuadriláteros en 2–D, y hexaedros en 3–D. 2. Mallados no estructurados: las celdas y los nodos de la malla no tienen un orden particular, es decir, los celdas o nodos cercanos a uno dado, no pueden identificarse directamente por sus índices. Los elementos de la malla, en este caso, son una mezcla de cuadriláteros y triángulos en 2–D y tetraedros y hexaedros en 3–D. El tipo de discretización espacial a emplear, depende del tipo de discretización de las ecuaciones empleado, así como de la estructura interna de datos empleados para resolver el flujo. El método de las diferencias finitas, que discretiza las ecuaciones de Navier-Stokes en forma diferencial, requiere una malla estructurada de puntos en los que se guardan las variables de flujo. Por otro lado, en el método de los volúmenes finitos es necesario definir puntos de control en cada volumen generado. Precisamente es, en esos vértices, donde se guardan las variables, que dependiendo del método que se emplee pueden ser definidos en el centro del elemento (cell-centred) o en los vértices (cell-vertex) de las celdas. Figura 3.3: Ejemplo de mallas con punto de control en el centro o en los vértices de los elementos.[18] 3.2.3. Definiciones generales de la metodología numérica La forma general de una ley de conservación para una cantidad escalar U, con fuentes volumétricas Q, sobre un volumen finito que incorpora flujos por las caras del volumen de control viene dada por: 𝜕 𝜕𝑡 ∫Ω 𝑈 𝑑Ω + ∮𝑆  𝐽⋅𝑑 𝑆=∫Ω 𝑄 𝑑Ω(3.1) Usando la ecuación general de transporte en su formulación integral que ya había sido definida 𝜕 𝜕𝑡 ∫𝑉 𝜌𝜙 𝑑𝑉 +∮𝐴 (𝜌𝑣𝜙 − Γ∇𝜙)⋅𝑑 𝐴=∫𝑉 𝑆𝜙𝑑𝑉 (3.2) 3.3. HERRAMIENTA NUMÉRICA: OPENFOAM 27 Basta con identificar que U = 𝜌 𝜙, 𝐽=𝜌 𝑣 𝜙 -Γ ∇ 𝜙y Q = 𝑆𝜙, así como que el volumen se representa por Ωen lugar de Vy el área por Sen vez de A, para entender que estamos ante la misma expresión. Esta nueva forma de expresar la ecuación de transporte es muy compacta e idéntica de forma matemática al sentido de los términos involucrados: término temporal (unsteady), término de flujo por las superficies de control (fluxes) y término fuente. 3.3. Herramienta numérica: OpenFoam 3.3.1. Definición de OpenFoam OpenFoam, es un paquete de software de libre divulgación que utiliza diferentes métodos de discretización para resolver ecuaciones en derivadas parciales. Este es principalmente usado en dinámica de fluidos aunque abarca otros campos (se especificará mejor en el capítulo 3.3.2). Se trata de una librería C++ a partir de la cual se crean ejecutables, conocidos como aplicaciones. Existen dos tipos de aplicaciones: solvers y utilidades. Los solvers se ocupan de resolver las ecuaciones que gobiernan un fenómeno físico. Las utilidades son aplicaciones que se utilizan antes o después de haber empleado el solvers y que básicamente son aplicaciones que nos permiten manipular datos y ajustarlos al problema que se está tratando y que será tratado por el solver. Un esquema de la organización global de OpenFoam se muestra en la figura. Figura 3.4: Esquema general de OpenFoam. [19] Queda claro que los pasos para usar OpenFoan vienen dados por la definición del problema, pre-procesar haciendo uso de las utilidades, resolver en función del problema establecido y finalmente post procesar con las diversas utilidades que se puedan necesitar. No obstante, 28 3. METODOLOGÍA no siempre se resolverán casos que estén implementados en Openfoam por lo que será necesario la modificación o creación mediante el lenguaje de programación C++ de solvers que puedan resolver dichos casos. 3.3.2. Capacidades de OpenFoam Una de las principales características que desde el diseño inicial se buscado es que sea un sistema abierto. Esta característica junto con sus características de software libre ha posibilitado el desarrollo de numerosos solvers así como herramientas de pre y post procesado para adaptarse a un gran abanico de casos. Entre los distintos solvers pueden destacarse: Fluidos incompresibles. •boundaryFoam: Solver estacionario para flujos turbulentos de 1D. •icoDyMFoam: Flujo laminar no estacionario con movimientos de malla asociados. •icoFoam: Flujo laminar no estacionario. •nonNewtonianIcoFoam. Flujo laminar no estacionario para flujos no estacionarios. •turbFoam: Flujos turbulentos incompresibles no estacionarios. Fluidos compresible. •coodles: LES para flujos compresibles. •rhoPimpleFoam: Flujos turbulentos para ventilación y transferencia de calor. Flujos multifase. •bubbleFoam: Dos fluidos incompresibles con una fase dispersa. •lesCavitatingFoam: Procesos no estacionarios de cavitación con turbulencia mediante LES. Combustiones. •dieselEngineFoam: Análisis de procesos de combustión con inyección en diésel spray. Transmisión de calor. 3.3. HERRAMIENTA NUMÉRICA: OPENFOAM 29 •buoyantFoam: Flujos turbulento para flujos compresibles para fenómenos de ventilación y transmisión. de calor, los efectos gravitacionales son considerados. •buoyantSimpleRadiationFoam: Flujos turbulentos para fluidos compresibles con fenómenos de radiación. Análisis de fuerzas. •solidDisplacementFoam: Solver no estacionario para pequeñas deformaciones elásticas, con efectos térmicos. 3.3.3. Casos En un problema que se quiera resolver es necesario una serie de carpetas y archivos que definan el problema y también que definan como se va a resolver dicho problema. Dichas carpetas son tres carpetas: constant,system y una carpeta temporal inicial. Esta última no tiene porque corresponder con el tiempo 0, es decir, puede partir en un tiempo posterior al tiempo 0 ya que puede darse el caso que se tengan simulaciones previas y se necesite continuar la simulación. Figura 3.5: Estructura carpetas de un caso en OpenFoam.[19] La carpeta constant contiene una carpeta llamada polyMesh, donde se encuentra toda la información acerca de como es la malla que ha creado, además de archivos que son necesarios para el solver ya que le proporcionan información sobre propiedades características del problema que se esté tratando. 30 3. METODOLOGÍA En system están los archivos correspondientes a la resolución. Debe haber como mínimo tres: controlDict,fvSchemes, y fvSolution. El controlDict es donde se imponen parámetros que tienen que ver con la ejecución del problema que se esté tratando, esto es por ejemplo el tiempo de comienzo y de finalización, el paso de tiempo o cada cuántas iteraciones se desea guardar el resultado. El fvSolution permite determinar como afrontar y como solucionar las ecuaciones de nuestro trabajo. Este tiene algoritmos de resolución de esquemas numéricos como el Gauss-Seidel. Por otro lado , el archivo fvSchemes es el que se contiene los esquemas numéricos que se usarán para resolver las ecuaciones de nuestro problema. Finalmente, para empezar la simulación se debe tener una carpeta temporal inicial en la que se establezcan las condiciones iniciales de las diferentes variables implicadas en el problema que se vaya a tratar. (a) Sin simular (b) Simulada Figura 3.6: Ejemplos de directorio antes y después de llevar a cabo una simulación 3.4. Simulación en OpenFoam Una vez introducidas las características principales de OpenFoam, se hace necesario destacar los diferentes pasos que se necesitan para hacer una simulación de forma más específica. 3.4.1. Generación de mallas mediante BlockMesh Lo primero que se hace es crear la geometría y proceder a mallarla. Dado que tratar de crear una malla discreta “a mano” es un trabajo largo y algo complicado, OpenFoam nos 3.4. SIMULACIÓN EN OPENFOAM 31 presenta la utilidad blockMeshDict que permite realizar dicha tarea. 3.4.1.1. Cabecera Lo primero que se encuentra al abrir el fichero blockMesh es la cabecera. En esta simplemente se observa que el fichero es un fichero OpenFoam y las unidades en las que se trabajaran. /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.1.1 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile {version 2.0; format ascii; class volScalarField; object p; } //*************************************// convertToMeters 1; 3.4.1.2. Definición de vértices A la hora de definir la geometría se necesitarán las coordenadas espaciales (x , y, z) de los diferentes puntos de la geometría que tienen importancia significativa, es decir, puntos que representan un extremo, un inicio de otra geometría, etc. vertices (//back up 32 3. METODOLOGÍA (0.1015 0 -0.5) //Punto 0 (0.203 0 -0.5) //Punto 1 (0.812 0 -0.5) //Punto 2 (0.812 0.143543 -0.5) //Punto 3 (0.143543 0.143543 -0.5) //Punto 4 (0.071771 0.071771 -0.5) //Punto 5 (0.812 0.812 -0.5) //Punto 6 (0.143543 0.812 -0.5) //Punto 7 (0 0.812 -0.5) //Punto 8 (0 0.203 -0.5) //Punto 9 (0 0.1015 -0.5) //Punto 10 (-0.1015 0 -0.5) //Punto 11 (-0.203 0 -0.5) //Punto 12 (-0.812 0 -0.5) //Punto 13 (-0.812 0.143543 -0.5) //Punto 14 (-0.143543 0.143543 -0.5) //Punto 15 (-0.071771 0.071771 -0.5) //Punto 16 (-0.812 0.812 -0.5) //Punto 17 (-0.143543 0.812 -0.5) //Punto 18 . . . . . //front down (0.812 -0.143543 0.5) //Punto 51 (0.143543 -0.143543 0.5) //Punto 52 (0.071771 -0.071771 0.5) //Punto 53 (0.812 -0.812 0.5) //Punto 54 (0.143543 -0.812 0.5) //Punto 55 (0 -0.812 0.5) //Punto 56 (0 -0.203 0.5) //Punto 57 (0 -0.1015 0.5) //Punto 58 (-0.812 -0.143543 0.5) //Punto 59 3.4. SIMULACIÓN EN OPENFOAM 33 (-0.143543 -0.143543 0.5) //Punto 60 (-0.071771 -0.071771 0.5) //Punto 61 (-0.812 -0.812 0.5) //Punto 62 (-0.143543 -0.812 0.5) //Punto 63 ); 3.4.1.3. Definición de bordes Para geometrías curvas se necesita definir los puntos que intervienen en ella, ya que de no definirsedichospuntos se uniráncomo si estuvieranenen linea entresí. Paraello definiremos que por esos puntos pasa una línea circular mediante la palabra arc, que básicamente lo que hace es indicarle al programa que esos puntos están unidos por una curva. Para poder llevar a cabo esto necesitaremos tres puntos; el inicial, el final y el intermedio (son los que definen por donde va geometría). El inicial y el final ya se tienen definidos de cuando se definieron los vértices, por lo tanto queda definir el central, haciendo uso otra vez de las coordenadas espaciales. edges (//up arc 0 5 (0.093774 0.038842 -0.5) arc 5 10 (0.038842 0.093774 -0.5) arc 19 24 (0.093774 0.038842 0.5) arc 24 29 (0.038842 0.093774 0.5) arc 11 16 (-0.093774 0.038842 -0.5) arc 16 10 (-0.038842 0.093774 -0.5) arc 30 35 (-0.093774 0.038842 0.5) arc 35 29 (-0.038842 0.093774 0.5) . . . . . 34 3. METODOLOGÍA arc 1 39 (0.187548 -0.077685 -0.5) arc 39 44 (0.077685 -0.187548 -0.5) arc 20 52 (0.187548 -0.077685 0.5) arc 52 57 (0.077685 -0.187548 0.5) arc 47 44 (-0.077685 -0.187548 -0.5) arc 12 47 (-0.187548 -0.077685 -0.5) arc 31 60 (-0.187548 -0.077685 0.5) arc 60 57 (-0.077685 -0.187548 0.5) ); 3.4.1.4. Definición de bloques La definición de los bloques consiste en definir una lista de vértices que definen una geometría y untipo específico demallado. Paraempezaradefinir dichos bloques esnecesario indicárselo al programa, por ello se usa la palabra blocks. Una vez hecho esto se definirán tres cosas: 1. Los vértices que componen el bloque (aunque estemos trabajando geometrías 2D siempre se definen en 3D los bloques). 2. Numero de celdas con la que se procederá a mallar según las coordenadas espaciales. 3. El ratio de expansión de las celdas que básicamente es su orientación en el bloque según las coordenadas espaciales. blocks // hex (5 4 9 10 24 23 28 29) (38 24 1) simpleGrading (1 1 1) //0 hex (0 1 4 5 19 20 23 24) (38 24 1) simpleGrading (1 1 1) //1 hex (1 2 3 4 20 21 22 23) (144 24 1) simpleGrading (1 1 1) //2 hex (4 3 6 7 23 22 25 26) (144 144 1) simpleGrading (1 1 1) //3 hex (9 4 7 8 28 23 26 27) (24 144 1) simpleGrading (1 1 1) //4 hex (15 16 10 9 34 35 29 28) (38 24 1) simpleGrading (1 1 1) //5 hex (12 11 16 15 31 30 35 34) (38 24 1) simpleGrading (1 1 1) //6 hex (13 12 15 14 32 31 34 33) (144 24 1) simpleGrading (1 1 1) //7 hex (14 15 18 17 33 34 37 36) (144 144 1) simpleGrading (1 1 1) //8 hex (15 9 8 18 34 28 27 37) (24 144 1) simpleGrading (1 1 1) //9 3.4. SIMULACIÓN EN OPENFOAM 41 out 336 674 ok (non-closed singly connected) Top 336 674 ok (non-closed singly connected) Down 336 674 ok (non-closed singly connected) in 336 674 ok (non-closed singly connected) cylinder 192 384 ok (non-closed singly connected) back 117888 118656 ok (non-closed singly connected) front 117888 118656 ok (non-closed singly connected) Checking geometry... Overall domain bounding box (-0.812 -0.812 -0.5) (0.812 0.812 0.5) Mesh (non-empty, non-wedge) directions (1 1 0) Mesh (non-empty) directions (1 1 0) All edges aligned with or perpendicular to non-empty directions. Boundary openness (9.49811e-17 -7.52451e-17 4.38594e-15) OK. Max cell openness = 2.8836e-16 OK. Max aspect ratio = 2.39071 OK. Minimum face area = 8.98725e-06. Maximum face area = 0.00664288. Face area magnitudes OK. Min volume = 8.98725e-06. Max volume = 2.80826e-05. Total volume = 2.60502. Cell volumes OK. Mesh non-orthogonality Max: 43.8649 average: 7.51961 Non-orthogonality check OK. Face pyramids OK. Max skewness = 0.488529 OK. Coupled point location match (average 0) OK. Mesh OK. End Finalmente, observamos el resultado en ParaView usando el comando paraFoam. 3.4.2. Elección del solver Para poder elegir el solver es muy importante tener claro las características del problema en cuestión y como no, que se pretende calcular. Como se ha comentado en el apartado 42 3. METODOLOGÍA Figura 3.7: Malla creada. 3.3.2 existen numerosos solvers. Nosotros hemos elegido el solver pimpleDyMFoam ya que es un solver que permite analizar fluidos newtonianos, incompresibles para una malla en movimiento. Simplemente se debe indicar que a la hora de utilizar este solver, se presenta la similitud a los otros en cuanto a la introducción de características propias del fluido como laminar o turbulento pero que difiere en que en este caso en concreto se necesita que se incluya un diccionario llamado dynamicMeshDict en el que se determinará la forma de tratar el movimiento de la malla por parte del solver. 3.4.3. Post-proceso En este punto este punto solo quedará decir que una vez simulado el problema para las características específicas que se tengan, será el momento de analizar los resultados. Para ello se usan programas como paraView donde podemos hacer análisis a nuestro sistema. 3.5. CONCEPTO DE MASA AÑADIDA 43 0.02 0.04 0.06 0.08 U Magnitude 0 0.0953 Figura 3.8: Resultados de la simulación. 3.5. Concepto de masa añadida Entrando en más profundidad en el problema es conveniente introducir el concepto de masa añadida, aunque se estudiará con mayor detalle en apartados sucesivos. Cuando un cuerpo se mueve a través de un fluido, éste debe empujar el fluido fuera de su camino para poder desplazarse. Si el cuerpo sufre una aceleración, el fluido circundante también debe acelerarse. El cuerpo se comporta entonces como si su masa fuera mayor y el incremento aparente de la masa se conoce como masa añadida del fluido. Σ𝐹= (𝑚+𝑚𝑎)𝜕2𝑢 𝜕𝑡2(3.3) Donde 𝜕2𝑢∕𝜕𝑡2corresponde a la aceleración del cuerpo y 𝑚𝑎está referido a la masa añadida. Dicha masa es proporcional a la densidad del fluido 𝜌y al volumen del cuerpo V, quedando definida como: 𝑚𝑎=𝜌𝑉 𝐶𝑚(3.4) Donde 𝐶𝑚es el coeficiente de masa añadida. 44 3. METODOLOGÍA 3.6. Ecuación de Morison Experimentos mostrados en Massie and Journée [8] nos dicen que existe una fuerza de arrastre proporcional a la velocidad y al diámetro de la estructura, en este caso al diámetro del cilindro, D por la ecuación 3.5. 𝐹𝐷(𝑡) = 1 2𝜌 𝐶𝐷𝐷 𝑢2 0 ⋅𝑢(𝑡)|𝑢(𝑡)|(3.5) donde: 𝐶𝐷= coeficiente adimensional de arrastre. D= diámetro del cilindro. 𝑢𝑜= Amplitud u(t) velocidad dependiente del tiempo. De la ecuación de arrastre 3.5 y de la ecuación de inercia 4.42 o 4.34 (según en el caso que nos encontremos) enunciamos una ecuación semi–empírica conocida como ecuación de Morison 3.6. 𝐹(𝑡) = 𝜋 4𝜌 𝐶𝑚𝐷2𝑑2u 𝑑t2+1 2𝜌 𝐶𝑑𝐷 𝑢2 0 ⋅𝑢(𝑡)|𝑢(𝑡)|(3.6) De donde: 𝐹𝑜𝑒𝑖(𝜔𝑡+𝛿)=𝑖1 2𝜌 𝜔2𝑢2 0𝐷 𝐶𝑑+𝑖𝜋 4𝜌 𝐷2𝜔2𝑢0𝐶𝑚(3.7) Obteniendo la expresiones para calcular los coeficientes de arrastre 3.9 y de inercia 3.8. 𝐶𝑚=4𝐹0sin 𝛿 𝜋 𝜌 𝐷2𝜔2𝑢0 (3.8) 𝐶𝐷=2𝐹0cos 𝛿 𝜌 𝜔2𝑢2 0𝐷(3.9) 3.6. ECUACIÓN DE MORISON 45 3.6.1. Método de Morison Una sugerencia hecha por Morison donde se proponía calcular los coeficientes 𝐶𝑑y𝐶𝑚 sin necesidad de usar ordenadores. Su acercamiento proponía lo que se ve en las ecuaciones 3.10 y 3.11. 𝐶𝑑=2𝐹 𝜌 𝐷 𝑢0|𝑢0|(3.10) 𝐶𝑚=4𝐹 𝜋𝜌 𝐷2𝜔 𝑢0 (3.11) Para que las ecuaciones 3.11 y 3.10 puedan ser aplicadas, tal y como se indica en [8] se tiene que obtener el valor de la fuerza en función de las siguientes condiciones. Valor de 𝐹para 𝐶𝑑cuando la aceleración es de valor cero. Valor de 𝐹para 𝐶𝑚cuando la velocidad es de valor cero. Figura 3.9: Ejemplo de señales (velocidad, aceleración y fuerza) Massie and Journée [8]. 46 3. METODOLOGÍA Aunque tal y como se indica en [8] este método no es el más preciso ya que pequeñas variaciones en la velocidad causan errores significativos principalmente en el cálculo de 𝐶𝑑. No obstante se calculan estos datos con la intensión de un posterior análisis de dichos resultados. 3.6.2. Cálculo del ángulo 𝛿mediante la transformada rápida de Fourier Para aplicar las ecuaciones 3.8 y 3.9 con la finalidad de obtener dichos coeficientes se precisa despejar el valor de la incógnita 𝛿para cada uno de los casos. Este ángulo de desfase entre fuerzas y desplazamientos se calculará aplicando una transformada rápida de Fourier implementada en el código de Matlab mostrado en 3.6.2. 1 clear all; 2 clc; 3 %% Carga de datos 4 load desplazamientos; 5 load fuerza; 6 %% Busca de posiciones en los vectores de datos 7 % Los valores son tomados para cada caso 8 % Fuerza 9 fi=max(find(fx <=0.8873)); % Valor inicial para el vector de fuerzas 10 ff=max(find(fx <=23.91) -1); % Valor final para el vector de fuerzas 11 dtf=fx(2)-fx(1); % Diferencial temporal para las fuerzas 12 % Desplazamiento 13 di=max(find(dx <=0.8873)); % Valor inicial para el vector de desplazamientos 14 df=max(find(dx <=23.91) -1); % Valor final para el vector de desplazamientos 15 dtd=dx(2)-dx(1); % Diferencial temporal para los desplazamientos 16 %% Vectores con los valores de estudio 17 fx=fx(fi:ff); % Vector de fuerzas eje X 18 fy=fy(fi:ff); % Vector de fuerzas eje Y 19 npf=length(fx); % Número de puntos de los vectores fuerza 20 dx=dx(di:df); % Vector de desplazamientos eje X 21 dy=dy(di:df); % Vector de desplazamientos eje Y 22 npd=length(dx); % Número de puntos de los vectores desplazamientos 23 %% Fast Fourier Transform 24 nyquistf=pi/dtf; % Nyquist para las fuerzas 3.6. ECUACIÓN DE MORISON 47 25 nyquistd=pi/dtd; % Nyquist para los desplazamientos 26 Wf=linspace(-nyquistf ,nyquistf ,npf+1); % Vector de frecuencias de las 27 % fuerzas 28 Wd=linspace(-nyquistd,nyquistd,npd); % Vector de frecuencias de los 29 % desplazamientos 30 Yf=fft(fy); % FFT de las fuezas 31 Yd=fft(dy); % FFT de los desplazamientos 32 %% Ordenar vectores 33 Yordf=zeros(1,npf+1); 34 for i=1:npf/2-1 35 Yordf(i)=conj(Yf((npf/2)+i)); % ordena la parte negativa de la 36 % transformada 37 end 38 for i=1:npf/2+1 39 Yordf((npf/2)+i)=conj(Yf(i)); % ordena la parte positiva de la 40 % transformada 41 end 42 Yordd=zeros(1,npd); 43 for i=1:(npd -1)/2 44 Yordd(i)=conj(Yd((npd+1)/2+i)); % ordena la parte negativa de la 45 % transformada 46 end 47 for i=1:(npd+1)/2 48 Yordd(((npd -1)/2+i))=conj(Yd(i)); % ordena la parte positiva de la 49 % transformada 50 end 51 %% Plot FFTs 52 % Fuerza 53 figure; 54 subplot(3,1,1); 55 plot(fx,fy); 56 title('Señal (fuerza)'); 57 subplot (3,1,2); 58 plot(Wf/2/pi,real(Yordf)); 59 title('Parte real de la frecuencia'); 60 xlabel('Hz'); 61 xlim ([-0.8 0.8]); % Para mostrar el pico a la frecuencia de excitación 62 subplot(3,1,3); 48 3. METODOLOGÍA 63 plot(Wf/2/pi,imag(Yordf)); 64 title('Parte imaginaria de la frecuencia'); 65 xlabel('Hz'); 66 xlim ([-0.8 0.8]); % Para mostrar el pico a la frecuencia de excitación 67 % Desplazamiento 68 figure; 69 subplot(3,1,1); 70 plot(fx,fy,'r'); 71 title('Señal (desplazamiento)'); 72 subplot (3,1,2); 73 plot(Wd/2/pi,real(Yordd),'r'); 74 title('Parte real de la frecuencia'); 75 xlabel('Hz'); 76 xlim ([-0.8 0.8]); % Para mostrar el pico a la frecuencia de excitación 77 subplot(3,1,3); 78 plot(Wd/2/pi,imag(Yordd),'r'); 79 title('Parte imaginaria de la frecuencia'); 80 xlabel('Hz'); 81 xlim ([-0.8 0.8]); % Para mostrar el pico a la frecuencia de excitación Capítulo 4 Resultados 4.1. Introducción Estructuras sumergidas o semi–sumergidas (parques eólicos, plataformas Offshore, estructuras flotantes, conductos de oil and gas, boyas) en corrientes marinas están habitualmente sometidas a fuertes vibraciones que presentan un serio daño para las partes mecánicas de las estructuras. Es por ello, que con la intención de analizar y eliminar esos problemas se propone estimar el valor de los coeficientes de arrastre y de inercia sobre una estructura cilíndrica. Para ello se considera un cuerpo bidimensional oscilando con una frecuencia y amplitud conocida en un fluido que inicialmente está en reposo y cuya principal característica es que es un fluido Newtoniano incompresible. 4.2. Solución numérica del problema del cilindro oscilante Para la resolución del problema del cilindro oscilante en el seno de un fluido no ideal inicialmente en reposo se ha tomado las dimensiones y las características del problema ofrecido en Troesch and Kim [16]. Esto es, se toman las dimensiones del cilindro y características propias del mismo como son la amplitud, su frecuencia angular ya que así se podrá simular este problema y tener unos resultados de referencia con los que poder comparar. Además, interesa destacar que el fluido con el que se trabajará tiene una viscosidad cinemática “𝜈”a 20 oC de 1,003 ⋅10−6 𝑚2 𝑠. 49 50 4. RESULTADOS 4.2.1. Estructura de la malla utilizada para los estúdios numéricos En un primer lugar, para poder llevar a cabo la resolución del problema comentado con Openfoam, se hace necesario definir y realizar una geometría. Está, presentada en la figura 4.1, tiene las siguientes características: Tiene la estructura cilíndrica centrada, d = 0.203 metros. Sobre la estructura cilíndrica centrada se establece una geometría cilíndrica de diámetro D =0.406 con la que se pretende hacer un análisis más exhaustivo en esa zona. Se añaden dos franjas centrales sobre el cilindro, una transversal y otra perpendicular de las que también interesará un análisis diferente. Longitud de las mismas l = 0.287086 metros. Quedarán cuatro bloques uno por cada esquina de nuestra geometría dando lugar a una geometría de largo total de L = 1.624 metros y ancho total de L = 1.624 metros. Figura 4.1: Geometría del sistema. 4.4. SOLUCIÓN ANALÍTICA PARA UN CILINDRO EN EL SENO DE UN FLUIDO IDEAL 57 masa del cuerpo se acelera sino que la masa del fluido lo hace también. Podremos decir que la fuerza adicional para acelerar el cuerpo viene dada por: 𝑔= −𝑚𝑎 𝜕2𝑢 𝜕𝑡2(4.6) donde 𝑚𝑎está referido a la masa añadida (véase 3.5). Ahora, conocida la velocidad de las partículas de agua en un sistemas de coordenadas polar bidimensional en que cada 𝑃punto del plano se determina por una distancia 𝑟y un ángulo 𝜃, donde 𝑟es la distancia de un punto 𝑃al origen 𝑂del plano y 𝜃es el ángulo formado entre el eje polar y la recta dirigida que va del origen 𝑂al punto 𝑃, unido a la definición de potencial de velocidad 𝜙se puede enunciar la ecuación de velocidad de flujo 4.7. Esta, será una ecuación a tener en cuenta como punto de partida en este análisis en el que se marca como objetivo encontrar una ecuación para el flujo potencial genérica que se pueda aplicar a los casos: 1. Cilindro rígido en reposo en el seno de un fluido oscilante. 2. Cilindro oscilante en el seno de un fluido inicialmente en reposo.  𝑈= ∇𝜙(4.7) La ecuación de Laplace para flujo potencial viene dada por 4.8 (definición de Laplaciano en coordenadas polares). ∇2𝜙(𝑟, 𝜃, 𝑡) = 𝜕2𝜙 𝜕𝑟2+1 𝑟 𝜕𝜙 𝜕𝑟 +1 𝑟2 𝜕2𝜙 𝜕𝜃2= 0 (4.8) donde 𝜙es la velocidad potencial del flujo y ∇2es el operador Laplaciano. Seguidamente, la ecuación 4.9 en la cual se enuncia la presión del fluido. 𝑝= −𝜌𝜕𝜙 𝜕𝑡 (4.9) Seguidamente se puede enunciar la ecuación 4.10 donde 𝐹𝜃(𝜃)y𝐹𝑟(𝑟)son funciones de 𝜙aún no conocidas. 𝜙(𝑟, 𝜃, 𝑡) = 𝐹𝑟(𝑟)𝐹𝜃(𝜃)𝑒(𝑖𝜔𝑡)(4.10) Aplicando la ecuación 4.10 a la ecuación 4.8: 𝐹𝜃(𝜃)𝜕2𝐹𝑟(𝑟) 𝜕𝑟2+1 𝑟𝐹𝜃(𝜃)𝜕𝐹𝑟(𝑟) 𝜕(𝑟)+1 𝑟2𝐹𝑟(𝑟)𝜕2𝐹𝜃(𝜃) 𝜕𝜃2= 0 (4.11) 58 4. RESULTADOS Dividiendo entre 𝐹𝑟(r) y 𝐹𝜃(𝜃) 1 𝐹𝑟(𝑟) 𝜕2𝐹𝑅(𝑟) 𝜕𝑟2+1 𝑟 𝐹𝑟(𝑟) 𝜕𝐹𝑟(𝑟) 𝜕𝑟 +1 𝑟2 1 𝐹𝜃 𝜕2𝐹𝜃(𝜃) 𝜕𝜃2= 0 (4.12) Multiplicando por 𝑟2y reordenando: 𝑟2 𝐹𝑟(𝑟) 𝜕2𝐹𝑟(𝑟) 𝜕𝑟2+𝑟 𝐹𝑟(𝑟) 𝜕𝐹𝑟(𝑟) 𝜕𝑟 = − 1 𝐹𝜃 𝜕2𝐹𝜃(𝜃) 𝜕𝜃2=𝜆(4.13) Donde 𝜆es una constante de valor propio. Seguidamente, por separación de variables, se obtienen dos ecuaciones. −1 𝐹𝜃 𝜕2𝐹𝜃(𝜃) 𝜕𝜃2=𝜆⇒𝜕2𝐹𝜃(𝜃) 𝜕𝜃2+𝐹𝜃𝜆= 0 (4.14a) 𝑟2𝜕2𝐹𝑟(𝑟) 𝜕𝑟2+𝑟𝜕𝐹𝑟(𝑟) 𝜕𝑟 −𝐹𝑟(𝑟)𝜆= 0 (4.14b) Donde: La solución de 4.14a es del tipo: 𝐹𝜃=𝛼1sin √𝜆 𝜃 +𝛼2cos √𝜆 𝜃 (4.15) La solución de 4.14b es del tipo: 𝐹𝑟=𝐴1𝑒−√𝜆𝑟−1 +𝐴2𝑒√𝜆𝑟(4.16) Ahora que se ha planteado las ecuaciones 4.15 y 4.16 se debe comprobar si la ecuación 4.16 es solución. Para ello, tomando en cuenta 4.17 aplicándolo sobre 4.16 y derivando parcialmente obtenemos 4.18. cosh(𝑥) = 𝑒𝑥+𝑒−𝑥 2senh(𝑥) = 𝑒𝑥−𝑒−𝑥 2(4.17) 𝜕𝐹𝑟 𝜕𝑟 =𝐴1sinh(√𝜆 𝑙𝑛(𝑟)) √𝜆 𝑟+𝑖 𝐴2cosh (√𝜆 𝑙𝑛(𝑟))√𝜆 𝑟(4.18) 𝜕2𝐹𝑟 𝜕𝑟 =𝐴1√𝜆cosh (√𝜆 𝑙𝑛(𝑟))√𝜆− sinh(√𝜆 𝑙𝑛(𝑟)) √𝜆 𝑟2 +𝑖 𝐴2√𝜆− cosh (√𝜆 𝑙𝑛(𝑟))√𝜆+ sinh(√𝜆 𝑙𝑛(𝑟)) √𝜆 𝑟2 (4.19) 4.4. SOLUCIÓN ANALÍTICA PARA UN CILINDRO EN EL SENO DE UN FLUIDO IDEAL 59 Sustituyendo en 4.14b y resolviendo queda: 𝐴1√𝜆{√𝜆cosh[√𝜆𝑙𝑛(𝑟)] − sinh[√𝜆 𝑙𝑛(𝑟)]} + 𝑖 𝐴2√𝜆{√𝜆sinh[√𝜆𝑙𝑛(𝑟)] − cosh[√𝜆 𝑙𝑛(𝑟)]} + 𝐴1sinh [√𝜆 𝑙𝑛(𝑟)] √𝜆+𝑖 𝐴2cosh[√𝜆 𝑙𝑛(𝑟)] √𝜆− {𝜆 𝐴1cosh[√𝜆 𝑙𝑛(𝑟)] + 𝑖 𝐴2𝜆sinh[√𝜆 𝑙𝑛(𝑟)]} = 0 (4.20) Con lo que se puede expresar 𝐹𝑟(r) como: 𝐹𝑟=𝐴1cosh(√𝜆 𝑙𝑛(𝑟)) + 𝑖 𝐴2sinh(√𝜆 𝑙𝑛(𝑟)) (4.21) Sabiendo que en este caso la parte imaginaria es nula y volviendo a tener en cuenta 4.17 se propone: 𝐹𝑟=𝐴1𝑒−√𝜆𝑟−1 +𝐴2𝑒√𝜆𝑟(4.22) De aquí: 𝜕2𝐹𝑟 𝜕𝑟 = −𝐴1𝑒−√𝜆𝑟−2 +𝐴2𝑒√𝜆(4.23) 𝜕2𝐹𝑟 𝜕𝑟2= 2𝐴1𝑒−√𝜆𝑟−3 (4.24) Sustituyendo en 4.14b obtenemos 4.25. 2𝐴1𝑒−√𝜆𝑟−1 −𝐴1𝑒−√𝜆𝑟−1 +𝐴2𝑒√𝜆𝑟−𝐴1𝑒−√𝜆𝑟−1 𝜆−𝐴2𝑒√𝜆𝑟 𝜆 = 0 (4.25) 𝐴1𝑒−√𝜆𝑟−1 (1 − 𝜆) + 𝐴2𝑒√𝜆𝑟(1 − 𝜆)=0 (4.26) (𝐴1𝑒−√𝜆𝑟−1 +𝐴2𝑒√𝜆𝑟) (1 − 𝜆)=0 (4.27) Ecuación que se verifica y demuestra la ecuación 4.16 ya que de 4.27 despejando 𝜆=1 podemos obtener 4.28. 𝐹𝑟(𝑟) = 𝐴′ 1𝑟−1 +𝐴′ 2𝑟(4.28) Seguidamente introduciendo los valores de 4.14a y 4.30 en 4.10 se obtiene 4.29 𝜙= (𝐴′ 1𝑟−1 +𝐴′ 2𝑟) (𝛼1sin 𝜃+𝛼2cos 𝜃)𝑒𝑖𝜔𝑡 (4.29) Nótese que la ecuación 4.29 aparece en el documento de referencia [2] donde sin embargo existen varias erratas significativas tal y como se ha demostrado en los desarrollos anteriores. 60 4. RESULTADOS 4.4.2. Cilindro oscilante en el seno de un fluido inicialmente en reposo Para el caso que se está tratando en este apartado, cilindro oscilante en el seno de un fluido inicialmente en reposo, se enunciarán sus condiciones de contorno. De estas se conoce que en el infinito, el fluido no está perturbado y en donde existen interferencias entre el fluido y la estructura, la velocidad normal es igual que la del cilindro. 𝑢𝑟= 0 𝑝𝑎𝑟𝑎 𝑟 = ∞ 𝑢𝑟=𝑢 𝑐𝑜𝑠𝜃 𝑝𝑎𝑟𝑎 𝑟 =𝑅(4.30) Tomando en cuenta la ecuación 4.29 donde se sabe que 𝛼1debe ser cero para que se cumpla la condición de contorno 4.30 para el caso que 𝑟=𝑅y también que 𝐴′ 2= 0 para que cumpla la condición de contorno 4.30 para el caso que 𝑟= ∞. (𝐴′ 1𝑟−1 +𝐴′ 2𝑟)(𝛼2cos 𝜃)𝑒𝑖𝜔𝑡 = 0 (4.31) Sustituyendo y sacando factor común 𝐴′ 1y𝛼2queda: 𝜙=𝑎(𝑅2 𝑟+𝑟) cos 𝜃𝑒𝑖𝜔𝑡 (4.32) Aplicando 4.9 teniendo en cuenta 4.30 y 4.32 obtenemos 4.33. 𝑃= −𝜌𝜕𝜙 𝜕𝑡 = −𝑎(𝑅2 𝑟)𝑐𝑜𝑠𝜃 𝑒𝑖𝜔𝑡 𝑖𝜔 (4.33) Por tanto, teniendo en cuenta la ecuación 4.33 la fuerza del fluido actuando sobre el cilindro será: 𝑔= − ∫2𝜋 0 𝑝(𝑟, 𝜃, 𝑡)|𝑟=𝑅𝑅cos 𝜃d𝜃= − ∫2𝜋 0 𝜌 𝑅 𝑎 𝑐𝑜𝑠𝜃 𝑖𝜔𝑒(𝑖𝜔𝑡)𝑅cos 𝜃d𝜃 = − 𝜌𝑅2𝑎𝑖𝜔𝑒(𝑖𝜔𝑡) ⏟⏞⏟⏞⏟ 𝑑2u 𝑑t2 ∫2𝜋 0 cos2𝜃d𝜃 ⏟⏞⏞⏞⏞⏞⏞⏟⏞⏞⏞⏞⏞⏞⏟ 𝜋 (4.34) Tomando en cuenta 4.6 y 4.34 se llega a que: 𝑔= − 𝑚𝑎 𝑑2𝑢 𝑑𝑡2;𝑚𝑎=𝐶𝑚𝜋𝜌𝑅2(4.35) este caso que es uno de los más simples vemos que 𝐶𝑚es igual a 1. 𝐶𝑚=𝑚𝑎 𝜋𝜌𝑅2= 1 (4.36) 4.4. SOLUCIÓN ANALÍTICA PARA UN CILINDRO EN EL SENO DE UN FLUIDO IDEAL 61 4.4.3. Cilindro rígido en reposo en el seno de un fluido oscilante En este caso, el cilindro no es oscilante. De hecho está quieto siendo el fluido el que se mueve bajo un movimiento armónico. Aquí se puede decir que las condiciones de contorno nos dicen que el fluido permanece quieto cerca del cilindro pero en el infinito siguen dicho movimiento armónico. 𝑢𝑟=𝑢 𝑐𝑜𝑠𝜃 𝑝𝑎𝑟𝑎 𝑟 = ∞ 𝑢𝑟= 0 𝑝𝑎𝑟𝑎 𝑟 =𝑅(4.37) Tomando en cuenta la ecuación 4.29 se sabe que 𝛼1debe ser cero para que se cumpla la condición de contorno 4.37 para el caso que 𝑟= ∞. (𝐴′ 1𝑟−1 +𝐴′ 2𝑟)(𝛼2cos 𝜃)𝑒𝑖𝜔𝑡 = 0 (4.38) Despejando 𝐴1: −𝐴′ 1 𝑅2+𝐴′ 2= 0 𝐴′ 1=𝐴′ 2𝑅2 (4.39) Sustituyendo 4.39 en 4.38 y sacando factor común queda 4.40. 𝜙=𝑎(𝑅2 𝑟+𝑟) cos 𝜃𝑒𝑖𝜔𝑡 (4.40) Aplicando 4.9 teniendo en cuenta 4.37 y 4.40 obtenemos 4.41. 𝑃= −𝜌𝜕𝜙 𝜕𝑡 = −𝑎(𝑅2 𝑟+𝑟)𝑐𝑜𝑠𝜃 𝑒𝑖𝜔𝑡 𝑖𝜔 (4.41) Por tanto, teniendo en cuenta la ecuación 4.41 la fuerza del fluido actuando sobre el cilindro será: 𝑔= − ∫2𝜋 0 𝑝(𝑟, 𝜃, 𝑡)|𝑟=𝑅𝑅cos 𝜃d𝜃= − ∫2𝜋 0 2𝜌 𝑅 𝑎 𝑐𝑜𝑠𝜃 𝑖𝜔𝑒(𝑖𝜔𝑡)𝑅cos 𝜃d𝜃 = − 2𝜌𝑅2𝑎𝑖𝜔𝑒(𝑖𝜔𝑡) ⏟⏞⏟⏞⏟ 𝑑2u 𝑑t2 ∫2𝜋 0 cos2𝜃d𝜃 ⏟⏞⏞⏞⏞⏞⏞⏟⏞⏞⏞⏞⏞⏞⏟ 𝜋 (4.42) 62 4. RESULTADOS Tomando en cuenta 4.6 y 4.42 se llega: 𝑔= − 𝑚𝑎 𝑑2𝑢 𝑑𝑡2;𝑚𝑎= 2 𝐶𝑚𝜋𝜌𝑅2(4.43) este caso se ve que 𝐶𝑚es igual a 2. 𝐶𝑚=𝑔 𝜋𝜌𝑅2= 2 (4.44) 4.4.4. Comentarios En lo que se ha venido desarrollando en el capítulo 4.4 no se ha encontrado bibliografía alguna que demuestre la relación entre las situaciones previamente comentadas y demostradas. De hecho, en algunos documentos de referencia estos conceptos no quedan claros. Por ejemplo Troesch and Kim [16] indica en su introducción que los casos 4.4.3 y 4.4.2 son iguales, algo que se ha demostrado que no es verdad ya que los coeficientes de 𝐶𝑚según el caso son completamente diferentes. 4.5. Solución numérica del cilindro oscilante 4.5.1. Definición del problema Es el momento de pasar a realizar las simulaciones. Dichas simulaciones se realizaran para un fluido no ideal, inicialmente en reposo y que se contrastaran los resultados con Troesch and Kim [16]. Téngase en cuenta que en Troesch and Kim [16] no se especifican las propiedades especificas del fluido. Las simulaciones que se realizarán con OpenFoam tal y como se ha comentado (Véase el apéndiceA para ver la estructura de un caso de OpenFoam), se supondrá que el fluido es agua a 20 oC y con 𝜈= 1.003 10−6 𝑚2 𝑠. Además, se calcularán los valores de 4.3 ya que son necesarios para realizar las simulaciones para los valores de ?? y??. Dichos valores se introducirán en OpenFoam y se resolverá las ecuaciones de Navier–Stokes para los régimenes láminar y turbulento. Para el régimen turbulento se usará el modelo k–épsilon (apartado 2.4.2) con los resultados de las ecuaciones obtenidas de la misma. Las simulaciones que se van a llevar a cabo para el caso laminar y turbulento están clasificadas según el valor de 𝛽. Dicha clasificación es, por un lado para 𝛽igual a 23200 y con 𝐾𝑐0.2 y 0.8, y por otro lado para 𝛽igual a 48600 y con 𝐾𝑐0.2 y 0.5. Estos valores de 𝛽y de 𝐾𝑐han sido seleccionados con el mismo valor con el que se presentan los resultados 4.5. SOLUCIÓN NUMÉRICA DEL CILINDRO OSCILANTE 63 en Troesch and Kim [16] ya que así se podrá comparar los resultados de dicha publicación con los que se calculan en este trabajo. Por otro lado, a partir de las dimensiones ya dadas del cilindro (diámetro igual a 0.203 metros) se procede al cálculo de los valores de los números adimensionales llegando a los resultados que se muestran: 𝐾𝑐/ Números Re 𝐀𝟎metros 𝜔𝐫𝐚𝐝 𝐬𝐞𝐠 T (seg) 0.2 4640.0 0.00646169068953 3.54794362455 1.77093718843 0.8 18560.0 0.0258467627581 3.54794362455 1.77093718843 Tabla: 4.1: Parámetros para 𝛽= 23200 𝐾𝑐/ Números Re 𝐀𝟎metros 𝜔𝐫𝐚𝐝 𝐬𝐞𝐠 T (seg) 0.2 9720.0 0.00646169068953 7.43233017901 0.845385653738 0.5 24300.0 0.0161542267238 7.43233017901 0.845385653738 Tabla: 4.2: Parámetros para 𝛽= 48600 4.5.2. Resultados En este apartado se presentan los resultados correspondientes a los valores de 𝐶𝑚y𝐶𝑑 obtenidos. 4.5.2.1. Resultados obtenidos desde Openfoam Estos resultados han sido obtenidos a partir de las simulaciones que se han realizado en Openfoam. Con los resultados ahí obtenidos se aplican las ecuaciones 3.8 y 3.9 dando los siguientes resultados. 64 4. RESULTADOS Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 1 0.2 0.0243809268677 1.00091920741 Simulación 2 0.8 0.099934169034 0.998987211059 Tabla: 4.3: Coeficientes para régimen laminar y 𝛽= 23200 Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 3 0.2 0.234427185678 1.00097060632 Simulación 4 0.5 0.107185302706 0.99363533139 Tabla: 4.4: Coeficientes para régimen laminar y 𝛽= 48600 Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 5 0.2 1.93903948485 1.02630147321 Simulación 6 0.8 0.450366362956 0.969726732117 Tabla: 4.5: Coeficientes para régimen turbulento y 𝛽= 23200 Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 7 0.2 1.93903948485 1.01647393317 Simulación 8 0.5 0.796731149697 1.01180918256 Tabla: 4.6: Coeficientes para régimen turbulento y 𝛽= 48600 4.5.2.2. Resultados Método Morison Seguidamente se presentan los resultados obtenidos mediante la aplicación del Método de Morison. Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 1 0.2 0.037490238398782565 0.990587258852439 Simulación 2 0.8 0.10522420358501083 0.9848628359714976 Tabla: 4.7: Coeficientes obtenidos por el Método Morison para régimen laminar y 𝛽= 23200 4.5. SOLUCIÓN NUMÉRICA DEL CILINDRO OSCILANTE 65 Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 3 0.2 0.24792303026210608 0.9882409763528609 Simulación 4 0.5 0.11627596958177051 0.9815624547584448 Tabla: 4.8: Coeficientes obtenidos por el Método Morison para régimen laminar y 𝛽= 48600 Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 5 0.2 1.8985206686097085 1.0161158298005322 Simulación 6 0.8 0.6345645785755704 1.0150167178166984 Tabla: 4.9: Coeficientes obtenidos por el Método Morison para régimen turbulento y 𝛽 = 23200 Simulación 𝐊𝐜𝐶𝐷𝐂𝐌 Simulación 7 0.2 1.5610783975026201 1.0063850118108275 Simulación 8 0.5 0.7741525004277429 1.0044990223817338 Tabla: 4.10: Coeficientes obtenidos por el Método Morison para régimen turbulento y 𝛽= 48600 4.5.3. Análisis de resultados Una vez vistos los resultados se precisa de comparación con resultados obtenidos por otros autores. Es por ello que se realizará la comparativa de resultados con el artículo Troesch and Kim [16]. Análisis de 𝐶𝑚 Se observa claramente valores prácticamente iguales entre las simulaciones que se han realizado para los regímenes laminar y turbulento. Como se ha demostrado, dichos valores tienen que valer 1. En el caso de nuestras simulaciones, dicho valor es prácticamente 1 con pequeñas variaciones de decimales relacionadas con los errores propios de las simulaciones y las FFT’s. Sin embargo vemos que en Troesch and Kim [16] le salen valores próximos a 3 para el caso de 𝛽= 23200, lo cual es evidentemente 66 4. RESULTADOS una errata a la luz incluso de las expresiones que ellos presentan en el mismo trabajo para su cálculo. Además, si se tienen en cuenta los valores de 𝐶𝑚en dados por Troesch and Kim [16] para el caso de 𝛽= 48600 se ve que los valores son cercanos a dos, algo que se ha demostrado tanto analíticamente como mediante simulaciones que no corresponde con el caso del cilindro oscilando en el seno de un fluido inicialmente en reposo ya que el valor de 𝐶𝑚debe ser 1. EL valor ofrecido, 2, corresponde al caso de un cilindro quieto en el seno de un fluido vibrante armónicamente. Finalmente comparando el Método Morison con los resultados obtenidos por las simulaciones vemos que se encuentran muy cercanos. Análisis de 𝐶𝑑 Para los valores de 𝐶𝑑lo primero destacable son los valores mostrados en las figuras 4.8 y 4.9 para el régimen laminar. Éstos son del orden de los valores ofrecidos por Troesch and Kim [16]. Por otro lado, si se compara la simulación realizada para régimen turbulento con los valores empíricos de 𝐶𝑑presentados en Troesch and Kim [16], se aprecia que los resultados numéricos reproducen la tendencia y los valores empíricos, teniendo en cuenta su dispersión. Además viendo los resultados obtenidos a partir del Método Morison y comparándolos con las simulaciones realizadas se ve que prácticamente se obtienen los mismos resultados, llegando a la conclusión que dicho método ha proporcionado resultados aptos. Apéndice A Utilidades en OpenFoam A.1. Posición /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.1.1 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile {version 2.0; format ascii; class pointVectorField; location "0"; object pointMotionU; } //*************************************// dimensions [0 1 0 0 0 0 0]; internalField uniform (0 0 0); boundaryField 73 74 A. UTILIDADES EN OPENFOAM {out {type uniformFixedValue; uniformValue (0 0 0); } Top {type uniformFixedValue; uniformValue (0 0 0); } Down {type uniformFixedValue; uniformValue (0 0 0); } in {type uniformFixedValue; uniformValue (0 0 0); } cylinder {type oscillatingDisplacement; amplitude (0.0258 0 0); omega 3.5479; value uniform (0 0 0); } back {type empty; } front {type empty; } } A.2. VELOCIDAD 75 // ************************************************************************* // A.2. Velocidad /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.1.1 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile {version 2.0; format ascii; class volVectorField; object U; } //*************************************// dimensions [0 1 -1 0 0 0 0]; internalField uniform (0 0 0); boundaryField {out {type zeroGradient; } Top {type zeroGradient; } Down 76 A. UTILIDADES EN OPENFOAM {type zeroGradient; } in {type zeroGradient; } cylinder {type movingWallVelocity; value uniform (0 0 0); } back {type empty; } front {type empty; } } // ************************************************************************* // A.3. Presión /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.1.1 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile {version 2.0; A.3. PRESIÓN 77 format ascii; class volScalarField; object p; } //*************************************// dimensions [0 2 -2 0 0 0 0]; internalField uniform 0; boundaryField {out {type fixedValue; value uniform 0; } Top {type fixedValue; value uniform 0; } Down {type fixedValue; value uniform 0; } in {type fixedValue; value uniform 0; } cylinder {type zeroGradient; } 78 A. UTILIDADES EN OPENFOAM back {type empty; } front {type empty; } } // ************************************************************************* // A.4. controlDict /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.1.1 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile {version 2.0; format ascii; class dictionary; location "system"; object controlDict; } //*************************************// application pimpleDyMFoam; startFrom startTime; startTime 0; A.4. CONTROLDICT 79 stopAt endTime; endTime 50; deltaT 0.000400211793997; writeControl timeStep; writeInterval 150; purgeWrite 0; writeFormat binary; writePrecision 6; writeCompression off; timeFormat general; timePrecision 6; runTimeModifiable true; adjustTimeStep no; maxCo 0.0163108540042; functions {probes {type probes; functionObjectLibs ("libsampling.so"); outputControl timeStep; outputInterval 1; probeLocations 80 A. UTILIDADES EN OPENFOAM (( 0.15 0 0 ) ( 0.3 0 0 ) ( 0.6 0 0 ) ); fixedLocations false; fields (p U ); }forces {type forces; functionObjectLibs ( "libforces.so" ); outputControl timeStep; outputInterval 1; patches (cylinder); pName p; UName U; rhoName rhoInf; log true; rhoInf 1000;//DENSIDAD AGUA CofR (0 0 0);// NO ES UN VALOR SOLICITADO ASÍ QUE 0 } } // ************************************************************************* // A.5. fvSchemes /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.1.1 | A.5. FVSCHEMES 81 | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile {version 2.0; format ascii; class dictionary; location "system"; object fvSchemes; } //*************************************// ddtSchemes {default Euler; } gradSchemes {default Gauss linear; grad(p) Gauss linear; } divSchemes {default none; div(phi,U) Gauss linear; div((nuEff*dev(T(grad(U))))) Gauss linear; } laplacianSchemes {default none; laplacian(nu,U) Gauss linear corrected; laplacian(rAUf,pcorr) Gauss linear corrected; laplacian(rAU,pcorr) Gauss linear corrected; laplacian(rAU,p) Gauss linear corrected; laplacian(rAUf,p) Gauss linear corrected; 82 A. UTILIDADES EN OPENFOAM laplacian(diffusivity,cellMotionU) Gauss linear uncorrected; laplacian(nuEff,U) Gauss linear uncorrected; laplacian(diffusivity,cellDisplacement) Gauss linear uncorrected; } interpolationSchemes {default linear; interpolate(HbyA) linear; } snGradSchemes {default corrected; } fluxRequired {default no; pcorr ; p ; } // ************************************************************************* // A.6. fvSolution /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.1.1 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile {