Escuela Politécnica Nacional: Facultad de Ingeniería Mecánica
Escuela Politécnica Nacional: Facultad de Ingeniería Mecánica
Escuela Politécnica Nacional: Facultad de Ingeniería Mecánica
NACIONAL
DIRECTOR:
Ing. HIDALGO DÍAZ VÍCTOR HUGO, D.Sc.
[email protected]
CODIRECTOR:
Ing. VALENCIA TORRES ESTEBAN ALEJANDRO, Ph.D.
[email protected]
Certifico que el presente trabajo fue desarrollado por MARTIN RICARDO VELASCO
BETANCOURT, bajo mi supervisión.
_____________________
Ing. Víctor Hugo Hidalgo Díaz D.Sc.
DIRECTOR DE PROYECTO
_____________________
Ing. Esteban Alejandro Valencia Torres Ph.D.
CODIRECTOR DE PROYECTO
i
DECLARACIÓN
Yo, MARTIN RICARDO VELASCO BETANCOURT, declaro bajo juramento que el trabajo
aquí descrito es de mi autoría; que no ha sido previamente presentado para ningún grado
o calificación profesional; y, que he consultado las referencias bibliográficas que se
incluyen en este documento.
A través de la presente declaración cedo mis derechos de propiedad intelectual
correspondiente a este trabajo, a la Escuela Politécnica Nacional, según lo establecido por
la Ley de Propiedad Intelectual, por su Reglamento y por la normativa institucional vigente.
_____________________
Martin Ricardo Velasco Betancourt
ii
DEDICATORIA
A mi madre:
A mi familia.
iii
AGRADECIMIENTO
A la Escuela Politécnica Nacional, por haber proporcionado la
formación necesaria para cumplir mi labor profesional con
confianza y determinación.
iv
ÍNDICE DE CONTENIDO
CERTIFICACIÓN ................................................................................................................ i
DECLARACIÓN ................................................................................................................. ii
DEDICATORIA ................................................................................................................. iii
AGRADECIMIENTO ......................................................................................................... iii
ÍNDICE DE CONTENIDO .................................................................................................. v
ÍNDICE DE FIGURAS ...................................................................................................... vii
ÍNDICE DE TABLAS ......................................................................................................... ix
RESUMEN......................................................................................................................... x
ABSTRACT ...................................................................................................................... xi
INTRODUCCIÓN .............................................................................................................. 1
Pregunta de Investigación ................................................................................................ 3
Objetivo general................................................................................................................ 3
Objetivos específicos ........................................................................................................ 3
1. MARCO TEÓRICO .................................................................................................. 4
1.1. Antecedentes........................................................................................................... 4
1.2. Turbinas tipo Francis y sus aplicaciones .................................................................. 5
1.3. Técnicas de mallado estructurado 3D en turbinas hidráulicas.................................. 6
1.4. Software libre y de código abierto de simulación de Dinámica de Fluidos
Computacional .................................................................................................................10
1.4.1. Dinámica de Fluidos Computacional (CFD) y sus etapas .......................................10
1.4.2. OpenFOAM como Software libre y de código abierto de Simulación CFD ..............13
1.5. Herramientas de OpenFOAM como software de CFD ............................................14
1.5.1. Herramientas y servicios útiles ...............................................................................14
1.5.2. Solvers disponibles .................................................................................................15
1.5.3. Estructura y ajuste de un caso de simulación .........................................................16
1.5.4. OpenFOAM Turbo Tools.........................................................................................16
2. METODOLOGÍA .....................................................................................................19
2.1. Modelo 3D del caso de estudio ...............................................................................20
2.2. Desarrollo y generación del mallado estructurado ..................................................21
2.2.1. Modelo Geométrico ................................................................................................21
2.2.2. Mallado Estructurado ..............................................................................................23
2.3. Importación y tratamiento del modelo 3D en OpenFOAM .......................................24
2.4. Verificación de la calidad de malla ..........................................................................29
2.5. Condiciones de simulación del modelo con metodología de interfaz de malla
móvil…… .........................................................................................................................29
2.5.1. Condiciones de operación ......................................................................................29
2.5.2. Simulación con Método SIMPLE: Condiciones iniciales y de borde. .......................31
2.5.3. Metodología de interfaz de malla móvil o dinámica .................................................33
v
2.5.4. Simulación con Método PIMPLE: Condiciones iniciales y de borde ........................35
2.5.5. Modelo de turbulencia ............................................................................................36
2.6. Simulación en OpenFOAM .....................................................................................37
2.6.1. Solver .....................................................................................................................37
2.6.2. Cálculo en Paralelo y características de Software/Hardware ..................................39
2.7. Post-proceso y Validación del modelo con resultados de estudios previos .............40
3. RESULTADOS Y DISCUSIÓN ...............................................................................41
3.1. Mallado ...................................................................................................................41
3.1.1. Independencia de Malla ..........................................................................................41
3.1.2. Calidad de Malla .....................................................................................................42
3.2. Residuales ..............................................................................................................43
3.3. Simulación con Método SIMPLE: Estado estable ...................................................44
3.4. Simulación con Método PIMPLE: Estado inestable.................................................56
4. CONCLUSIONES ...................................................................................................62
4.1. Conclusiones ..........................................................................................................62
4.2. Trabajos Futuros.....................................................................................................64
Referencias Bibliográficas ...............................................................................................65
ANEXOS..........................................................................................................................69
ANEXO I. PROCEDIMIENTO DE OBTENCIÓN DEL ARCHIVO DE MALLA .MSH DESDE
ICEM ANSYS...................................................................................................................69
ANEXO II. ARCHIVO EJECUTABLE PARA LA UNIÓN DEL DOMINIO TOTAL,
MANIPULACIÓN DEL DOMINIO DENTRO DE OPENFOAM Y FICHERO
TOPOSETDICT…. ...........................................................................................................70
ANEXO III. FICHEROS CON INFORMACIÓN DE CONDICIONES DE BORDE Y
CONDICIÓN DE MALLA DESLIZANTE ...........................................................................71
ANEXO IV. FICHERO DE CONTROL DE SIMULACIÓN CONTROLDICT ......................79
ANEXO V. ARCHIVO EJECUTABLE PARA INICIAR LA SIMULACIÓN Y FICHERO
USADO PARA EL CÁLCULO EN PARALELO .................................................................81
vi
ÍNDICE DE FIGURAS
Figura 1.1. Turbina Francis de la central hidroeléctrica San Francisco. ............................ 5
Figura 1.2. Técnicas de mallado estructurado................................................................... 7
Figura 1.3. Discretización del dominio en Multibloques. .................................................... 7
Figura 1.4. Estrategias de multibloques ............................................................................ 8
Figura 1.5. Comparación de diferentes tipos de mallas para una simulación RANS ......... 9
Figura 1.6. Diagrama de flujo del proceso de simulación CFD .........................................10
Figura 1.7. Estructura del directorio de un caso de estudio en OpenFOAM. ....................16
Figura 2.1. Estructura Metodológica del caso de estudio .................................................19
Figura 2.2. Rodete Turbina Francis y casa de máquinas de la Central San Francisco .....20
Figura 2.3. Modelo Geométrico del Rodete a) Mora y b) Guascal & Quispe. ...................21
Figura 2.4. Geometría mejorada del codo del tubo de descarga ......................................22
Figura 2.5. Diferencia en la discretización de la entrada del tubo de descarga a) Optimizada:
polígono de 12 lados, b) Guascal & Quispe: polígono de 6 lados. ...................................23
Figura 2.6. Metodología de multibloques. ........................................................................23
Figura 2.7. Transición de malla entre subdominios a) Guascal & Quispe: 1,2e+05
elementos, b) Optimizada: 9,7e+05 elementos. ...............................................................24
Figura 2.8. Estructura del caso de estudio. ......................................................................25
Figura 2.9. Diagrama de procesos para obtener la carpeta polyMesh con la información de
la malla final. ....................................................................................................................26
Figura 2.10. Esquema de carpeta de subdominio en el directorio creado para el tratamiento
del mallado. .....................................................................................................................27
Figura 2.11. Dominio computacional total. .......................................................................28
Figura 2.12. Subdominios a) Rodete, b) Tubo de descarga, c) Álabes directrices, d) Caja
espiral y Álabes Fijos. ......................................................................................................28
Figura 2.13. Algoritmo SIMPLE ........................................................................................31
Figura 2.14. Interfaces y superficies de transferencia. .....................................................34
Figura 2. 15. Algoritmo PIMPLE. ......................................................................................35
Figura 3.1. Valores residuales de presión en distintas mallas. .........................................41
Figura 3.2. Valores residuales para la simulación con método SIMPLE (estado estable): a)
modelo de turbulencia k-ε, b) modelo de turbulencia k-ω SST. ........................................43
Figura 3.3. Valores residuales para la simulación con método PIMPLE (estado inestable):
a) modelo de turbulencia k-ε, b) modelo de turbulencia k-ω SST. ....................................44
Figura 3.4. Potencia vs Q/Qn............................................................................................45
Figura 3.5. Potencia vs Q/Qn vs Error. .............................................................................46
Figura 3.6. Eficiencia vs Q/Qn. ........................................................................................47
vii
Figura 3.7. Perfil de velocidad a la salida del rodete: a) modelo de turbulencia k-ε, b) modelo
de turbulencia k-ω SST....................................................................................................48
Figura 3.8. Perfil de Velocidad en el plano de corte YZ, a) modelo de turbulencia k-ε, b)
modelo de turbulencia k-ω SST. ......................................................................................49
Figura 3.9. Perfil de Presión en el plano de corte YZ, a) modelo de turbulencia k-ε, b)
modelo de turbulencia k-ω SST. ......................................................................................50
Figura 3.10. Vista superior del perfil de Presión en álabes directrices y rodete: a) modelo
de turbulencia k-ε, b) modelo de turbulencia k-ω SST. ....................................................51
Figura 3.11. Vista inferior del perfil de presión en la salida del rodete: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST. .........................................................51
Figura 3.12. Perfil de Presión en las paredes del rodete y en sus álabes: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST. .........................................................52
Figura 3.13. Comparación de daños en el material con los campos de alta presión de álabes
directrices y álabes del rodete obtenidos en la simulación numérica. ..............................52
Figura 3.14. Perfil de Energía cinética turbulenta k en el plano de corte YZ: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST. .........................................................53
Figura 3.15. Vista en detalle del perfil de energía cinética turbulenta: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST. .........................................................54
Figura 3.16. Perfil de energía cinética turbulenta a la entrada del tubo de descarga: a)
modelo de turbulencia k-ε, b) modelo de turbulencia k-ω SST. ........................................54
Figura 3.17. Cuerda de vórtice obtenida mediante una iso superficie de Q-criterion = 64: a)
modelo de turbulencia k-ε, b) modelo de turbulencia k-ω SST. ........................................55
Figura 3.18. Perfil de presión en el plano de corte YZ: a) modelo de turbulencia k-ε, b)
modelo de turbulencia k-ω SST. ......................................................................................57
Figura 3.19. Perfil de Presión en el rodete a) modelo de turbulencia k-ε, b) modelo de
turbulencia k-ω SST.........................................................................................................57
Figura 3.20. Vista superior del Perfil de presión en álabes directrices y rodete: a) modelo
de turbulencia k-ε, b) modelo de turbulencia k-ω SST. ....................................................58
Figura 3.21. Vista superior de perfil de Presión en álabes directrices y paredes del rodete
sin contar el cubo de salida del rodete: a) modelo de turbulencia k-ε, b) modelo de
turbulencia k-ω SST.........................................................................................................58
Figura 3.22. Perfil de Energía cinética turbulenta en el plano de corte YZ: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST. .........................................................59
Figura 3.23. Perfil de velocidad en el plano de corte YZ: a) modelo de turbulencia k-ε, b)
modelo de turbulencia k-ω SST. ......................................................................................60
viii
Figura 3.24. Cuerda de vórtice obtenida mediante una iso superficie de Q-criterion = 64: a)
modelo de turbulencia k-ε, b) modelo de turbulencia k-ω SST. ........................................61
Figura 3.25. Perfil de presión en la entrada del tubo de descarga: a) modelo de turbulencia
k-ε, b) modelo de turbulencia k-ω SST. ...........................................................................62
Figura A.1. Submenú ICEM para escoger solver .............................................................69
Figura A.2. Ventana emergente para guardar el archivo .msh. ........................................69
ÍNDICE DE TABLAS
Tabla 1.1. Comparación entre modelos DNS, LES, RANS y DES....................................11
Tabla 1.2. Ventajas y desventajas del uso de paquetes de software libre y de código abierto
CFD. ................................................................................................................................12
Tabla 1.3. Herramientas (utilities) útiles en el proceso de simulación CFD. .....................14
Tabla 1.4. Solvers disponibles en OpenFOAM para flujos incompresibles. ......................15
Tabla 2.1. Estadísticos de la malla...................................................................................29
Tabla 2.2. Condiciones de operación. ..............................................................................30
Tabla 2.3. Condiciones de borde globales. ......................................................................32
Tabla 2.4. Condiciones de borde. ....................................................................................32
Tabla 2.5. Condiciones de borde .....................................................................................32
Tabla 2.6. Superficies de transferencia de datos .............................................................33
Tabla 2.7. Coeficientes del modelo de turbulencia que usa OpenFOAM..........................37
Tabla 2.8. Ajustes del Solver ...........................................................................................38
Tabla 2.9. Control de residuales y factores de relajación. ................................................38
Tabla 2. 10. Valores experimentales de rendimiento a diferentes condiciones de flujo
volumétrico. .....................................................................................................................40
Tabla 3.1. Valores de Yplus para cada subdominio. ........................................................42
Tabla 3.2. Número Ω para el dominio total. ......................................................................42
Tabla 3.3. Resultados de simulación con el modelo de turbulencia k–ε. ..........................47
Tabla 3.4. Resultados de simulación con modelo de turbulencia k–ω SST. .....................48
ix
RESUMEN
En el presente trabajo se desarrolló una metodología para estudiar las turbinas tipo Francis
utilizando el software CFD libre y de código abierto OpenFOAM, junto con un marco de
referencia múltiple y una malla móvil estructurada para casos de estado estable e inestable,
respectivamente. OpenFOAM ofrece claras ventajas sobre otros paquetes CFD como, por
ejemplo, un mejor proceso de cálculo en paralelo, bajo costo, un mejor control de los
parámetros de simulación y una optimización de los recursos computacionales. La malla
estructurada usada fue obtenida de estudios previos, sin embargo, se efectuaron algunas
modificaciones en su geometría debido a la presencia de puntos de desequilibrio en la
simulación numérica. Se realizaron dos tipos de simulaciones numéricas junto con los
modelos de turbulencia k-ε y k-ω SST, primero usando el algoritmo SIMPLE para casos de
estado estable y el algoritmo PISO + SIMPLE (PIMPLE) para casos de estado inestable.
Los resultados obtenidos para las simulaciones de estado estable reproducen el fenómeno
con una gran aproximación. Se estudiaron diez puntos de operación diferentes en las
simulaciones en estado estable, obteniendo un error de 5.02% para el caso de potencia al
usar el caudal máximo. Debido a que las simulaciones en estado inestable se desarrollaron
en una revolución, se capturó parcialmente el fenómeno.
x
ABSTRACT
The development of a methodology to study Francis turbines using the open source CFD
software OpenFOAM along with a multiple reference frame and a structured moving mesh
for steady and unsteady state cases respectively, was carried out in this research work.
OpenFOAM offers clear advantages over others CFD packages, such as a better parallel
calculation process, low cost, better control of simulation parameters and optimization of
computational resources. The structured mesh used was obtained from previous studies,
nevertheless, some modifications were done due to the presence of imbalance points in the
numerical simulation. Two types of numerical simulations were performed together with
turbulence models k-ε and k-ω SST, first using the SIMPLE algorithm for steady state cases,
and the PISO + SIMPLE (PIMPLE) algorithm for unsteady state cases. The results obtained
for the steady state simulations reproduce the phenomenon with a great approximation. Ten
different operating points were studied in the steady state simulations, obtaining an error of
5.02% for the power when using the maximum value of flow. Due to the one revolution
approach of the unsteady state simulations, the phenomenon was partially captured.
xi
ESTUDIO DE MALLADO ESTRUCTURADO 3D ROTATIVO EN
OPENFOAM PARA APLICACIONES EN TURBINAS TIPO FRANCIS
INTRODUCCIÓN
En la última década, el gobierno ecuatoriano ha ejecutado varios proyectos estratégicos
con el objetivo de obtener un mayor aprovechamiento energético. Teniendo en cuenta que
las principales fuentes de energía renovable son hídricas, actualmente el país cuenta con
una capacidad efectiva (MW) en energía hidráulica de 62,59% según el informe emitido por
la ARCONEL a octubre de 2019 [1]. La mayor cantidad de proyectos hidroeléctricos que
se encuentran en funcionamiento están en la cuenca oriental, entre los que se puede
mencionar Coca Codo Sinclair, Quijos, La Sopladora, Mazar-Dudas y Delsitanisagua. Por
otra parte, existen proyectos hidroeléctricos que se encuentran en la cuenca occidental
tales como San Francisco y Toachi Pilatón [2]. De las tecnologías disponibles en
generación de energía mediante turbinas hidráulicas, las centrales La Sopladora, San
Francisco, Toachi Pilatón y Quijos utilizan turbinas tipo Francis [3]. Es decir, las turbinas
tipo Francis representan un pilar fundamental en la generación de energía eléctrica en el
Ecuador.
A pesar de ello, este tipo de turbinas sufren un gran desgaste y pérdida de rendimiento por
las condiciones de operación a las que están sometidas al trabajar en las fuentes hídricas
del país [4]. Por lo que, un adecuado estudio de este tipo de turbinas puede conllevar a la
mitigación de distintos tipos de problemas. En el caso específico de la central San
Francisco, existen problemas de cavitación-erosión en el rodete debido a la interacción de
micropartículas de sedimentos de origen volcánico con la turbina [5]; tomando en cuenta
que, la turbina no fue diseñada para soportar este tipo de impactos. En este contexto el
presente estudio forma parte del proyecto interno PIJ 17-13 dentro del cual se investiga y
evalúa modelos de cavitación-erosión aplicados a prevención de daños en turbinas
hidráulicas de las centrales hidroeléctricas del Ecuador.
En torno a este aspecto, desmontar una turbina para su estudio representaría una pérdida
económica y de generación eléctrica para el país, por lo que no es una metodología viable.
Sin embargo, gracias al creciente desarrollo que ha tenido la tecnología CFD (dinámica de
fluidos computacional) y el poder computacional de los ordenadores, es posible realizar
simulaciones numéricas complejas de turbinas hidráulicas con el fin de predecir su
comportamiento bajo condiciones reales de trabajo [6]. Existen trabajos previos realizados
por Mora [5] y Guascal & Quispe [7], en los que se llevó a cabo un estudio numérico con
1
malla no estructurada y un estudio numérico con malla estructurada respectivamente.
Ambos estudios fueron realizados en el entorno de software comercial, y su dominio
computacional fue una turbina tipo Francis con características similares a la turbina tipo
Francis N°2 de la central hidroeléctrica San Francisco.
No obstante, existen varios paquetes de software en los que se puede realizar el estudio
numérico de una turbina hidráulica. En el caso de un software comercial es necesario
contar con la licencia de la empresa que lo comercializa, la misma que tiene un alto costo.
Por otro lado, el uso de software libre y de código abierto tiene ventajas económicas y
analíticas. OpenFOAM ha tenido un vasto crecimiento a nivel industrial y académico en los
últimos años, siendo el software más utilizado en varios estudios involucrados con
mecánica de fluidos, transferencia de calor, etc. [8], por lo que resulta conveniente su uso
frente a otros paquetes de software comerciales.
2
Pregunta de Investigación
¿Cómo influye el mallado estructurado rotativo 3D en la simulación de una turbina tipo
Francis en el entorno de OpenFOAM?
Objetivo general
Estudiar el mallado estructurado 3D rotativo en OpenFOAM para aplicaciones en turbinas
tipo Francis.
Objetivos específicos
• Revisar el estado del arte de estudios de turbinas tipo Francis y mallado
estructurado rotativo.
• Desarrollar un modelo 3D de la geometría de estudio con base en escaneo 3D.
• Definir una metodología para realizar el mallado estructurado 3D con interfaces
móviles.
• Simular el modelo obteniendo resultados para validación y estudio numérico.
• Comparar los resultados obtenidos en estudios previos usando Software comercial
con los obtenidos usando OpenFOAM.
• Validar el estudio con datos que se encuentran en el estado del arte.
3
1. MARCO TEÓRICO
1.1. Antecedentes
En el año 2018, bajo el aval de la Escuela Politécnica Nacional, se da inicio al proyecto
interno junior PIJ 17-13 titulado “Investigación y evaluación de modelos de cavitación-
erosión aplicados a la prevención de daños en turbinas hidráulicas de las centrales
hidroeléctricas del Ecuador”. Esta iniciativa permite el estudio de problemas que se
evidencian en la operación diaria de las centrales hidroeléctricas del país para que la
academia aporte con soluciones prácticas a la industria.
Los estudios previos se han centrado en analizar los problemas que presentan las turbinas
tipo Francis que operan en algunas centrales del país. El primer estudio lo llevó a cabo el
Ing. Christian Andrés Mora Sánchez M. Sc en el año 2018 [5], cuyo objetivo principal era
analizar la eficiencia de una turbina tipo Francis de la central hidroeléctrica San Francisco,
ubicada entre la cuenca media y baja del Río Pastaza, en la ciudad de Baños de Agua
Santa, provincia de Tungurahua. Para ello, se realizó un proceso de escaneo 3D de todos
los elementos de la turbina y con métodos de ingeniería inversa se obtuvo las geometrías
finales en archivos STL y CAD, para posteriormente ser usadas en el software de
simulación CFD ANSYS CFX. El conjunto del dominio computacional fue discretizado en
una malla no estructurada con 68,5 millones de elementos. Como resultados del estudio
se logró estimar la eficiencia general de la turbina con un error menor al 10%, y se
determinó que la holgura entre álabes directores y placas de desgaste tiene una relación
directa con la fuga de caudal, generando mayores pérdidas hidráulicas y con ello una
descompensación en la eficiencia.
Posteriormente, los ingenieros Edison Javier Guascal Sanguña y Alexander Pedro Quispe
Quispe en el año 2019 [7], realizaron un estudio en el cual se optimizó la geometría
realizada por Mora, al eliminar líneas y superficies excedentes y discretizar el dominio
computacional en una malla estructurada con 2,9 millones de elementos. Esto permitió
tener un mayor control en los resultados de los perfiles de velocidad y presión, obteniendo
cambios menos bruscos de presión en las secciones de transición como son el rodete y los
álabes directores. Además, gracias a la optimización de mallado el tiempo de cálculo
computacional decreció en un 33% y se obtuvo un estimado de la eficiencia general con
un error menor al 5%. Estos resultados mostraron la ventaja que se tiene al trabajar con
una malla estructurada, sin embargo, el procedimiento para obtener dicha malla es más
complicado.
4
1.2. Turbinas tipo Francis y sus aplicaciones
Las turbinas tipo Francis son las turbinas de reacción más utilizadas en la industria de
generación eléctrica porque operan en un rango moderado de cabeza hidráulica, entre 50
y 700 [m], sin embargo, su máxima eficiencia se obtiene entre 100 y 300 [m] [9]. Además,
debido a su elevada eficiencia (~90-95%) bajo un amplio rango de condiciones de trabajo,
contribuyen con cerca del 60% de capacidad efectiva de producción de energía hidráulica
en el mundo [10] [11]. De ahí que, los dos proyectos hidroeléctricos más grandes a nivel
mundial, Itaipú (Brasil-Paraguay) y Las Tres Gargantas (China), poseen turbinas tipo
Francis con una capacidad de generación de 700 [MW] cada una [9]. La turbina está
ubicada entre una fuente de agua de alta presión y una salida de agua de baja presión,
generalmente en la base de la represa. Los principales componentes de una turbina tipo
Francis son el caracol o caja espiral, álabes fijos, álabes guías o directrices, rodete y tubo
de descarga, y pueden ser visualizados en la Figura 1.1.
El agua ingresa al dominio de la turbina a través de la caja espiral, siguiendo hacia los
álabes fijos, los mismos que aseguran la integridad estructural de la caja espiral y dirigen
el flujo hacia los álabes directrices. La forma en espiral del caracol permite tener una
distribución uniforme de velocidad en los álabes fijos [13]. Los álabes directrices llevan el
agua con una dirección tangencial hacia los álabes internos del rodete, y el impacto del
flujo con el borde de ataque de los álabes es lo que produce el movimiento de rotación.
Nótese que los álabes directrices son ajustables en ángulo para permitir la operación de la
turbina con suaves cambios de flujo. El tubo de descarga es un difusor que conecta la
salida del rodete con la salida del dominio de la turbina, y su principal objetivo es recuperar
la energía perdida en el rodete en forma de un incremento de presión para así elevar la
5
cabeza hidráulica efectiva de la turbina [10]. Es esencial que la salida del tubo de descarga
se encuentre sumergida para asegurar que la turbina se mantenga llena de agua [13].
Los algoritmos que se emplean para generar este mallado generalmente poseen
complejas técnicas iterativas de suavizado que favorecen el emparejamiento de los
elementos con los bordes o límites de la geometría. Los diferentes softwares de generación
de mallado estructurado son comúnmente usados en el área de CFD, donde el
6
alineamiento de los elementos es requerido para poder capturar ciertos fenómenos físicos
[15]. Para la generación de malla existen varias técnicas ampliamente usadas, y se pueden
clasificar en cuatro tipos principales que se muestran en la Figura 1.2.
La descomposición del dominio en subdominios o bloques ha sido una de las técnicas más
usadas para la generación de mallado estructurado en geometrías complejas. La estrategia
de multibloques fue formulada por Lee et al. [17], y consiste en segmentar el dominio
computacional o bloque general en bloques más pequeños para posteriormente generar la
malla para cada bloque individual. La Figura 1.3. muestra la idea esquemática de la
descomposición de un bloque en tres subdominios alrededor de un perfil aerodinámico. En
este ejemplo la combinación de bloques 1 y 3 usa una combinación de malla H-H y una
combinación de malla C-H es usada en el bloque 2 [16].
7
Figura 1.4. Estrategias de multibloques
(Fuente: [16])
Dentro de la técnica de mallado por multibloques existen algunas variantes o métodos que
dependen de la continuidad o discontinuidad de las líneas de la malla a través de los límites
de bloque. La Figura 1.4.a) muestra el método de Superposición (Overset method), el cual
no intenta unir las mallas de un bloque con su bloque vecino ni los límites de los bloques,
llegando a ser una forma fácil de generar la malla. La gran ventaja del uso de este método
es que, gracias a su independencia de coincidencia de las mallas, permite tratar problemas
de mallas móviles. Sin embargo, para una transferencia adecuada de datos entre mallas
es necesario la formulación de una región de interfaz. La Figura 1.4.b) muestra el método
Parchado (Patched multiblock), el cual tampoco intenta unir las líneas de malla de los
bloques, pero sí los límites de los bloques vecinos, asegurando una continuidad de límites
de bloque más no la continuidad de líneas de malla. La transferencia de datos entre mallas
presenta menos dificultades que en el caso de Superposición. Este método ofrece la
ventaja de poder refinar la malla en regiones específicas sin necesidad de imponer un
refinamiento en secciones innecesarias. La Figura 1.4.c) muestra el método Compuesto
(Composite multiblock), que se puede decir es un caso especial del método Parchado, con
la diferencia que esta variante busca una continuidad de líneas de malla a través de los
límites de bloques, asegurando un mismo número de nodos en los límites para los bloques
vecinos. No obstante, su mayor desventaja es la dependencia que posee el refinamiento
8
en los límites de malla, ya que, si en un bloque se refina en mayor grado un sector limítrofe,
esto induce un refinamiento en el bloque adyacente debido a su correspondencia de líneas
de malla en la interfaz. En consecuencia, la mayoría de los softwares comerciales han
optado por usar la metodología multibloques Compuesta en sus generadores de malla [16].
En este contexto, a pesar del gran desarrollo de herramientas amigables para el usuario,
la descomposición en multibloques para un mallado estructurado se mantiene como una
actividad que consume un elevado recurso de tiempo. Además, se han llevado a cabo
varios estudios computacionales para determinar qué método de mallado conlleva a
soluciones numéricas más precisas con el uso de RANS [16]. En este contexto, la Figura
1.5. muestra que la metodología de multibloques lleva a obtener resultados más precisos,
así como también ofrece mayor facilidad de uso.
Figura 1.5. Comparación de diferentes tipos de mallas para una simulación RANS
(Fuente: [16])
Finalmente, debido a lo expuesto, quedan en evidencia las ventajas que posee la técnica
de mallado estructurado, especialmente con la metodología de multibloques. Es por ello
que, en la última década muchos de los estudios de vanguardia que se han llevado a cabo
en el tema de turbinas hidráulicas y turbomaquinaria en general, utilizan esta técnica de
mallado [18].
9
1.4. Software libre y de código abierto de simulación de Dinámica
de Fluidos Computacional
Dinámica de Fluidos Computacional o CFD, por sus siglas en inglés, es una herramienta
usada para predecir el desempeño de sistemas y subsistemas ingenieriles, donde el flujo
de un fluido en específico juega un papel importante. En consecuencia, no es una sorpresa
que CFD sea ampliamente usado en un vasto rango de áreas de ingeniería tales como
automotriz, hidráulica o aeronáutica. También, la simulación CFD es usada como una
herramienta complementaria de diseño (y en algunos casos es la herramienta principal) en
las metodologías tradicionales, las mismas que incurren en altos costos experimentales y
de prototipado. Una de las mayores ventajas para incorporar simulación CFD en el proceso
ingenieril de diseño es la posibilidad de comparar varios diseños de forma paralela a un
costo reducido. CFD es esencialmente un método para la resolución de un conjunto de
ecuaciones diferenciales parciales que describen a un fluido como sistema, las mismas
que representan el principio de conservación de masa, momentum y energía. Además,
pueden agregarse ecuaciones de transporte adicionales para modelar otras propiedades
como la energía cinemática turbulenta k y la tasa de disipación ε, cuando se use un modelo
de turbulencia k-ε standard [19]. Para realizar un estudio adecuado usando simulación CFD
se debe seguir una metodología establecida, tal como se muestra en la Figura 1.6.
10
La primera etapa llamada Pre-proceso es considerada la más desafiante. Aquí es donde el
usuario evalúa la naturaleza y las características dominantes del problema, y si es posible
se definen consideraciones de simplificación para casos muy complejos. Una vez que el
fenómeno físico dominante es identificado, el usuario selecciona el enfoque de simulación
(Tabla 1.1.) que será usado, ya sea DNS , DES, LES o RANS, y esta decisión generalmente
se basa en los recursos computacionales disponibles. El siguiente paso es la discretización
espacial del dominio, es decir, construir la malla en concordancia con el enfoque de
simulación escogido anteriormente. Ciertamente, es necesario enfatizar que el desarrollo
del mallado es el paso más crítico y de mayor demanda de tiempo en una simulación típica
de CFD; y, como se menciona en la Sección 1.3, la precisión de los resultados finales
depende en gran medida de la calidad del mallado. Y para finalizar esta etapa, se definen
las condiciones de borde del problema basándose en el fenómeno físico a estudiar, el
enfoque de simulación y las restricciones físicas-geométricas del dominio.
(Fuente: [8])
11
más apropiado y el software CFD. Esta selección estará basada en varios factores como
son: costo de software, limitaciones de licencia, complejidad de uso, solvers disponibles,
herramientas de pre-proceso y post-proceso, paralelización, etc. En especial el usuario
debe escoger entre usar un paquete de software de simulación comercial o un software
libre y de código abierto. Ahora bien, la ventaja de un paquete de software comercial es el
acceso a soporte técnico que tendría el usuario junto con una detallada documentación.
Sin embargo, para usar este tipo de paquetes es necesario una licencia de un elevado
costo [8]. Por otra parte, un paquete de software libre y de código abierto presenta algunas
ventajas, como no presentar la restricción económica referente a la adquisición de la
licencia. En la Tabla 1.2., se muestran las ventajas y desventajas de usar paquetes de
software CFD libre y de código abierto.
Tabla 1.2. Ventajas y desventajas del uso de paquetes de software libre y de código abierto CFD.
Ventajas Desventajas
(Fuente: [19])
12
es lineal, sino más bien cíclico, si el usuario detecta que los resultados no concuerdan con
el fenómeno físico estudiado, se redefine el caso de estudio y se regresa a la etapa de Pre-
proceso nuevamente.
13
1.5. Herramientas de OpenFOAM como software de CFD
OpenFOAM es un paquete de software libre y de código abierto que permite crear archivos
ejecutables o aplicaciones que pueden ser de dos categorías: solvers, que son diseñados
para la resolución de problemas específicos, y utilities que son usados para el manejo y
manipulación de información. Una gran ventaja es la variedad de solvers/utilities que
posee, ya que el usuario podrá modificar y ajustar cualquier aplicación que tenga
naturaleza similar para su caso de estudio [22].
OpenFOAM posee varios servicios y herramientas (utilities) que ayudan al usuario a través
del proceso de simulación CFD, desde el pre-proceso, generación de malla, conversión de
formato de malla, simulación y post-proceso de resultados [19]. Un resumen de algunas
herramientas que posee OpenFOAM en cada una de las etapas de la simulación CFD se
muestra en la Tabla 1.3.
Pre-proceso
Generación de malla
Conversión de malla
14
fluentMeshToFoam Convierte mallas guardadas en extensión .msh a formato
OpenFOAM, p. ej. Archivos de malla generados en ICEM
ANSYS FLUENT
Post-proceso
(Fuente:[19])
Solver Descripción
(Fuente: [19])
15
1.5.3. Estructura y ajuste de un caso de simulación
16
una importante parte en la simulación CFD contemporánea. La simulación de máquinas
multi-álabe conlleva varios desafíos más allá de los estándares en aerodinámica y flujos
turbulentos. En primer lugar, la complejidad de las geometrías multi-álabe y multi-etapa en
equipos rotativos se empareja con la complejidad física del flujo inestable rotacional y la
interacción álabe-estela. Siguiendo en esta línea, la predicción de rendimiento de una
máquina rotativa requiere una simulación precisa que considere el desprendimiento de la
capa límite turbulenta en presencia de un gradiente adverso de presiones. Esto implica que
la malla debe tener un refinamiento adecuado. Además, la rotación implica necesariamente
tener un movimiento relativo entre un dominio llamado rotor y otro llamado estator. A pesar
de tener efectos claramente variantes en el tiempo, para su uso práctico en ingeniería se
toman condiciones cuasi estáticas que buscan tener un costo computacional bajo (en
tiempo de respuesta), lo que permite realizar estudios de condiciones de funcionamiento
fuera de diseño y optimización de geometría.
17
• Para simulaciones multi-etapa, el marco de referencia rotacional múltiple MRF
(multiple rotating frame) es usado, una vez más representado en una malla estática.
La rotación se toma en cuenta con la restricción que la condición de borde externa
de cada región rotativa forme una superficie de rotación. La ecuación de momentum
se formula en términos de velocidad absoluta o relativa teniendo en cuenta la
interfaz (superficie de rotación) rotor-estator, lo que complica al modelo.
• GGI cíclico (GGI cyclic), es usado para establecer un manejo adecuado de límites
que no tienen emparejamiento periódico. Se usa en simulaciones de estado estable,
pero presenta un factor limitante en la calidad de malla.
18
2. METODOLOGÍA
En este capítulo se describe la metodología utilizada durante el desarrollo del presente
estudio. Este estudio busca evidenciar el beneficio que se puede obtener al utilizar
OpenFOAM como software de simulación CFD para turbinas hidráulicas tipo Francis
usando mallado estructurado en su dominio.
19
procede a importar el modelo al entorno de OpenFOAM, donde se manipulan los distintos
dominios del modelo para posicionarlo tomando las coordenadas del rodete como la
referencia principal de ensamble. Luego, se realiza el ensamble de los dominios y se define
los sets de celdas, y al finalizar esta etapa se verifica la calidad de malla.
El siguiente paso es definir las condiciones de borde y las iniciales, tomando en cuenta las
condiciones de operación de la turbina tipo Francis de la Central hidroeléctrica San
Francisco. En función de las condiciones de borde e iniciales se selecciona el tipo de solver,
el modelo de turbulencia y la metodología de malla móvil. Así es como se llega a la etapa
de simulación, donde son llevados a cabo los cálculos computacionales de los campos de
las variables de interés. Los resultados obtenidos son visualizados en ParaView 5.6.3,
donde se muestran los campos de velocidad y presión en todo el dominio de la turbina, así
como también las zonas de turbulencia. Finalmente, se realiza la validación tanto de
resultados como de la malla estructurada con estudios previos, siendo este un paso
necesario para identificar la confiabilidad de los resultados obtenidos en OpenFOAM.
El modelo 3D usado en el presente estudio es el de una turbina tipo Francis de eje vertical
de la central hidroeléctrica San Francisco, ubicada entre la cuenca media y baja del Río
Pastaza, en la ciudad de Baños de Agua Santa, provincia de Tungurahua.
Específicamente, la unidad Francis a estudiar es la primera que entró en funcionamiento
(Unidad 02) en mayo del 2007. Esta central tiene dos unidades de eje vertical con
capacidad de generación de 115 [MW] cada una [24].
Figura 2.2. Rodete Turbina Francis y casa de máquinas de la Central San Francisco
(Fuente: propia)
20
2.2. Desarrollo y generación del mallado estructurado
El modelo geométrico y la construcción del mallado estructurado para cada subdominio fue
llevado a cabo en su totalidad en el software comercial ICEM, debido principalmente a la
versatilidad y robustez que posee para las actividades de pre-proceso, brindando un fácil
manejo mediante su amigable interfaz gráfica.
La turbina tipo Francis posee una geometría bastante compleja que requiere un alto nivel
de detalle en su dimensionamiento, y por esto recrear una turbina tipo Francis a partir de
los planos de taller y montaje sería una tarea muy exhaustiva y fuera del alcance del
presente trabajo. Gracias a los estudios previos realizados por Mora [5] y Guascal & Quispe
[7], se obtuvo en el primer caso el dominio geométrico mediante un proceso de ingeniería
inversa con escaneado 3D, mientras que en el segundo caso, se optimizó la geometría
obtenida al eliminar líneas, puntos y superficies excedentes como se puede ver en la Figura
2.3. Cabe mencionar que la malla fue construida teniendo en cuenta una apertura de álabes
directrices de 24° con respecto a su posición cerrada, siendo esta la condición de apertura
máxima de los álabes directrices.
Figura 2.3. Modelo Geométrico del Rodete a) Mora y b) Guascal & Quispe.
(Fuente: [7])
21
Figura 2.4. Geometría mejorada del codo del tubo de descarga
(Fuente: propia)
El área que más problemas presentó era el codo del tubo de descarga, ya que no tenía
una continuidad en las líneas de superficie, debido a que las curvas fueron aproximadas a
líneas pendientes tangentes a dichas curvas. Para solucionar este particular se crearon
más puntos en el modelo geométrico con el fin de aproximar una curva y eliminar las líneas
pendientes como se muestra en la Figura 2.4. También, se discretizó el dominio del tubo
de descarga con más cantidad de bloques con el fin de asegurar que la malla generada se
acople adecuadamente a la modificación geométrica del codo del tubo. Esto permitió tener
más control sobre la cantidad de elementos en la malla final del tubo de descarga.
Por otra parte, se tenía un inconveniente en la superficie de entrada del tubo de descarga,
ya que esta superficie circular estaba discretizada como un hexágono, mientras que la
superficie de salida del rodete (ambas superficies estaban en contacto) tenía una
discretización de su superficie circular de un polígono de 13 lados, teniendo así mayor
cantidad de puntos para transferencia de datos. En la Figura 2.5. se puede ver la
modificación de esta superficie circular en la entrada del tubo de descarga, la misma que
se discretizó usando un polígono de 12 lados.
22
Figura 2.5. Diferencia en la discretización de la entrada del tubo de descarga a) Optimizada:
polígono de 12 lados, b) Guascal & Quispe: polígono de 6 lados.
(Fuente: propia)
En este proceso es importante tener en cuenta la densidad de malla que posee cada uno
de los subdominios de la turbina con la finalidad de evitar transiciones bruscas entre mallas
23
y asegurar el éxito de la metodología de mallas móviles [25]. Ya que, si no se asegura una
densidad adecuada los cálculos computacionales pueden presentar inducción de
resultados o peor aún, una falla integral de cálculo que lleve a un corte abrupto en la etapa
de simulación. La Figura 2.7. muestra la mejora de mallado estructurado en el caso del
tubo de descarga asegurando una transición de malla suave entre el rodete y la entrada
del tubo de descarga.
Figura 2.7. Transición de malla entre subdominios a) Guascal & Quispe: 1,2e+05 elementos, b)
Optimizada: 9,7e+05 elementos.
(Fuente: propia)
24
Figura 2.8. Estructura del caso de estudio.
(Fuente: propia)
Una vez creado el directorio general del caso, es necesario tener la información del dominio
computacional que será simulado. En esta sección se detalla cómo fue el tratamiento de
los subdominios de la turbina en el entorno de OpenFOAM para obtener la subcarpeta
polyMesh, la misma que contiene la información del modelo geométrico y del mallado
estructurado del dominio total. En la Figura 2.9. se muestra el procedimiento llevado a cabo
para la obtención de la información del modelo.
25
Figura 2.9. Diagrama de procesos para obtener la carpeta polyMesh con la información de la malla
final.
(Fuente: propia)
Para empezar, se creó un directorio con el fin de manipular únicamente los archivos con
información de la malla estructurada y poder realizar el ensamblaje final. Primero, dentro
del nuevo directorio se creó una carpeta para cada subdominio en el cual se incluían las
subcarpetas constant y system en cada una. Los archivos .msh de cada subdominio son
copiados a la carpeta específica con su nombre. Luego, se abre la terminal dentro de la
carpeta y se ejecuta el comando fluentMeshToFoam, este comando sobrescribe la
subcarpeta polyMesh, que está dentro de constant, y captura toda la información de la
geometría y el mallado estructurado, y los ficheros que crea son mostrados en la Figura
2.10. Se debe tomar en cuenta que el hecho de transformar información de malla y
geometría de un software comercial a un software libre y de código abierto conlleva una
pérdida de información mínima que debe ser verificada.
Con la información total del mallado de cada subdominio lo que prosigue es verificar la
adecuada orientación de cada uno. Ya que, el elemento principal es el rodete, este será el
elemento con respecto al cuál los demás subdominios serán posicionados. Para ello se
ejecutan los comandos transformPoints -translate y transformPoints -rollPitchYaw, los
mismos que sirven para trasladar o rotar cierto subdominio en función de un vector para
que se sitúe en la ubicación adecuada [22]. Por ejemplo, existía un problema de alineación
del tubo de salida con respecto al rodete (desviación hacia la izquierda del eje de rotación
del rodete y una superposición en la zona inferior) y lo mismo ocurría con el dominio de los
álabes directrices.
26
Figura 2.10. Esquema de carpeta de subdominio en el directorio creado para el tratamiento del
mallado.
(Fuente: propia)
27
las propiedades de flujo a las celdas. Estos comandos modifican la subcarpeta polyMesh
creando ficheros donde se guardan los diferentes sets de zonas de celdas. Por ejemplo,
uno de sus ficheros llamado rotor contiene las zonas de celdas que pertenecen al rodete,
y estos archivos juegan un importante papel en la metodología de malla móvil. En el
ANEXO II se describe el código del archivo ejecutable, así como también la sintaxis de los
comandos y ficheros antes mencionados. En la Figura 2.11. y Figura 2.12. se muestra el
dominio computacional total y los subdominios que componen a la turbina.
Figura 2.12. Subdominios a) Rodete, b) Tubo de descarga, c) Álabes directrices, d) Caja espiral y
Álabes Fijos.
(Fuente: propia)
28
2.4. Verificación de la calidad de malla
Por otra parte, se usan dos números adimensionales, el primer número adimensional Ω
propuesto por Hidalgo V. et al [26], el cual permite evaluar la calidad de malla, así como
también el tiempo de simulación que ocupará el caso de estudio. El número Ω es la relación
entre el número total de elementos NE y el número total de nodos ND. El segundo número,
es el parámetro conocido como Yplus, que permite verificar si la malla coincide con el rango
de valores permitidos para el uso del modelo de turbulencia seleccionado. Los valores
obtenidos de Yplus y el número Ω obtenidos en el presente caso de estudio se observan
en la Sección 3.1.2.
29
Tabla 2.2. Condiciones de operación.
Condición Valor
Fabricante VATECH
Potencia Nominal 115 [MW]
Caudal Nominal de trabajo 58 [m3/s]
Altura aguas abajo 6,34 [m]
Rendimiento 95,5 %
Diámetro entrada Caja 3 [m]
Espiral
Gracias a los datos obtenidos de las condiciones de operación se puede estimar el rango
de turbulencia que tiene el flujo dentro del dominio, y para ello se calcula el número de
Reynolds. El número de Reynolds es un valor adimensional que se obtiene como una
relación de fuerzas inerciales y fuerzas viscosas [27]. La Ecuación. 2.1., muestra el cálculo
de este parámetro adimensional.
Lc ∗ U 3 ∗ 8.204
𝑅e = = = 2.225 e + 07
υ 1,106e − 06
(Ecuación 2.1.)
De acuerdo con Fox et. al [27], un valor de Reynolds de 2300 representa el límite en el cual
el flujo pasa de ser laminar a turbulento; por lo tanto, el flujo del caso de estudio es
altamente turbulento.
30
2.5.2. Simulación con Método SIMPLE: Condiciones iniciales y de borde.
Las condiciones de borde globales se muestran en la Tabla 2.3., las mismas que fueron
usadas en estudios previos por Mora [5] y Guascal & Quispe [7]. El cálculo de la velocidad
de entrada se lo consigue mediante la relación de caudal Q y área transversal A, V = Q/A.
31
Por otro lado, el valor de la presión debe ser normalizado, es decir, el valor de presión se
lo divide para el valor de la densidad del fluido a la temperatura de trabajo. Esto debido a
que OpenFOAM trabaja con presiones que se miden en unidades de [m2/s2].
Condición Valor
Inlet (In_Spiral_Case) 8,827 [m/s]
Outlet (Out_Draft_Tube) 80,088 [m2/s2]
(Fuente: propia)
A continuación, en la Tabla 2.4. y Tabla 2.5. se muestran las condiciones de borde usadas
para los campos de variables. Las partes o subdominios mostrados en ambas tablas son
aquellos considerados como paredes (wall), y también la entrada y salida del dominio total
considerados como patch.
Campos
Parte U p nut
Blade (wall) movingWallVelocity zeroGradient nutUSpaldingWallFunction
Blade1 (wall) movingWallVelocity zeroGradient nutUSpaldingWallFunction
Wall_Runner movingWallVelocity zeroGradient nutUSpaldingWallFunction
Wall_Draft_Tube fixedValue zeroGradient nutUSpaldingWallFunction
Wall_Guide_Vane fixedValue zeroGradient nutUSpaldingWallFunction
Wall_Spiral_Case fixedValue zeroGradient nutUSpaldingWallFunction
In_Spiral_Case fixedValue zeroGradient calculated (Tabla2.1.)
(patch) (Tabla2.2.)
Out_Draft_Tube inletOutlet fixedValue calculated (Tabla2.1.)
(patch) (Tabla2.2.)
(Fuente: propia)
Campos
Parte omega k epsilon
Blade (wall) omegaWallFunction kqRWallFunction epsilonWallFunction
Blade1 (wall) omegaWallFunction kqRWallFunction epsilonWallFunction
Wall_Runner omegaWallFunction kqRWallFunction epsilonWallFunction
Wall_Draft_Tube omegaWallFunction kqRWallFunction epsilonWallFunction
Wall_Guide_Vane omegaWallFunction kqRWallFunction epsilonWallFunction
Wall_Spiral_Case omegaWallFunction kqRWallFunction epsilonWallFunction
In_Spiral_Case fixedValue fixedValue fixedValue
(patch) (Ecuación 2.3 ) (Ecuación 2.2) (Ecuación 2.4)
Out_Draft_Tube zeroGradient zeroGradient zeroGradient
(patch)
(Fuente: propia)
32
Las superficies de transferencia de datos siguen la metodología AMI, explicada en la
Sección 2.5.3., y las mismas se muestran en la Tabla 2.6.
Parte Tipo
In_Runner cyclicAMI
Out_Runner cyclicAMI
In_Draft_Tube cyclicAMI
Out_Guide_Vane cyclicAMI
In_Guide_Vane cyclicAMI
Out_Spiral_Case cyclicAMI
(Fuente: propia)
La condición de cyclicAMI permite que dos superficies sean tratadas como si estas
estuviesen físicamente conectadas; se usan para geometrías repetitivas. Cada par de
superficies de conexión debe tener un área similar que cumpla con una tolerancia dada en
el fichero de boundary que se encuentra dentro de polymesh [22].
Por otra parte, se debe ajustar dos ficheros que se encuentran en la carpeta constant, con
el fin de dar las características del fluido y el enfoque de solución para el mismo. El primero,
transportProperties, contiene la información de las propiedades del fluido como, por
ejemplo, el valor de la viscosidad cinemática del fluido. El segundo, turbulenceProperties,
contiene información sobre el carácter turbulento o laminar del fluido y el modelo de
turbulencia que será usado en la simulación numérica.
Para conseguir que esta metodología funcione se utiliza el enfoque AMI (Arbitrary Mesh
Interface) de transferencia de datos. El enfoque AMI está basado en la creación de una
“supermalla”, la misma que está definida por elementos de intersección de dos mallas
33
consecutivas en un paso de tiempo. Este enfoque, aplicado para la interpolación de mallas
adaptivas, tiende a reducir errores de continuidad, mejorando así la eficiencia numérica del
cálculo. El concepto de “supermalla” le da una eficiencia mayor que el enfoque GGI
(General Grid Interface). La mejora que presenta AMI radica en que las integrales que se
resuelven no lo hacen en áreas superpuestas, sino más bien se calculan en la “supermalla”
que aumenta la capacidad de precisión en la interpolación y con ello la mejor estimación
de propiedades de conservación. En resumen, AMI permite realizar simulaciones de
dominios adyacentes desconectados al crear una malla intermedia llamada “supermalla”,
lo que mejora notablemente el proceso de interpolación y lo hace más simple [32]. En la
Figura 2.14. se muestran las superficies de transferencia de datos del dominio total de la
turbina.
34
2.5.4. Simulación con Método PIMPLE: Condiciones iniciales y de borde
El estado inestable o no estacionario difiere con respecto al estacionario en que los valores
que toman los campos cambian con el tiempo. Además, como se vio en la Sección 2.5.3.
el estado inestable o no estacionario trabaja con interfaces de malla móviles. Puesto que,
el dominio y las condiciones de borde no cambian, las condiciones establecidas en las
Tablas 2.3., 2.4., 2.5., y 2.6. se mantienen. Sin embargo, existe un fichero que debe ser
incluido en este caso para indicarle a OpenFOAM cuál dominio debe rotar, cuál debe ser
estático y sobre qué punto debe girar. El fichero tiene el nombre de dynamicMeshDict y se
encuentra dentro de la subcarpeta constant, ya que este fichero contiene información de
las superficies que tendrán movimiento deslizante relativo. Además, contiene la
35
información de la velocidad de giro del rodete (Tabla 2.1.), eje de rotación y las
coordenadas del origen. En el ANEXO III se muestra la sintaxis interna del fichero
dynamicMeshDict al igual que los ficheros que contienen las condiciones de borde para el
estado estable e inestable.
Debido al flujo altamente turbulento que sale del rodete en condiciones de operación, es
de vital importancia usar un modelo de turbulencia apropiado que pueda predecir la
eficiencia y las propiedades del flujo con una adecuada precisión. Este tipo de flujo no
puede ser modelado completamente solo con estudios en estado estable, sin embargo, el
uso de modelos de turbulencia adecuados permite una buena predicción de las
propiedades medias (promedio) del flujo y las pérdidas hidráulicas en esta zona de estudio.
Los modelos de turbulencia conocidos como modelos de dos ecuaciones han sido
ampliamente usados y validados para estudios de turbomaquinaria [14]. En el presente
estudio, tanto para la simulación en estado estable e inestable se usó el modelo k–ε y k–ω
SST.
Los valores iniciales de los parámetros k (energía cinética turbulenta), ω (tasa de disipación
turbulenta específica) y ε (tasa de disipación turbulenta) son calculados usando las
expresiones de Wilcox. D [36], así como también el valor de IT (intensidad turbulenta %).
No obstante, debido a que el flujo es altamente turbulento, se puede usar un valor igual o
mayor a 10% para la intensidad turbulenta IT [14]. Los valores iniciales mostrados fueron
calculados usando el caudal experimental de 62,4 [m3/s] mostrado en la Tabla 2.10.
3 3 m2
k= (U 𝐼𝑇 )2 = (8,827 ∗ 0,1)2 = 1,1689 [ 2 ]
2 2 s
(Ecuación 2.2.)
1 1 1 √1,1689
− √k − √k
ω = Cμ 4 = Cμ 4 = 0,09−4 = 9,39978 [ s−1 ]
lturb 0,07 Dh 0,07 ∗ 3
(Ecuación 2.3.)
3 3 3
3 3
𝑘2 𝑘2 3 1,0095 2 m2
𝜀= Cμ4 = Cμ4 = 0,09 4 = 0,9889 [ 3 ]
lturb 0,07 Dh 0,07 ∗ 3 s
(Ecuación 2.4.)
2 k 2 1,1689
𝐼𝑇 % = 100 ∗ √ 2
= 100 ∗ √ ∗ = 10
3U 3 8,8272
(Ecuación 2.5.)
36
La constante del modelo de turbulencia Cµ tiene un valor aproximado de 0,09 [36], U es la
velocidad de entrada en [m/s], Dh es el diámetro hidráulico de entrada en [m] y lturb es la
longitud turbulenta en [m]. En la Tabla 2.7. se muestran los coeficientes para cada uno de
los modelos de turbulencia.
2.6.1. Solver
La mayoría de los solvers de OpenFOAM para dinámica de fluidos usan los algoritmos
PISO y SIMPLE o una combinación llamada algoritmo PIMPLE (PISO + SIMPLE). Estos
algoritmos son procedimientos iterativos para ecuaciones de acoplamiento de momentum
y conservación de masa, PISO y PIMPLE son usados para problemas no estacionarios,
mientras que SIMPLE para problemas en estado estacionario [22]. Varios estudios se han
llevado a cabo y se ha validado el uso de estos algoritmos [6] [11] [14].
37
Tabla 2.8. Ajustes del Solver
PIMPLE SIMPLE
p 1e-03 1e-03
U 1e-04 1e-04
k 1e-05 1e-05
omega 1e-05 1e-05
epsilon 1e-05 1e-05
Factores de relajación
p 0,3
U 0,7
k 0,7
omega 0,7
epsilon 0,7
(Fuente: propia)
38
La Tabla 2.9. muestra los valores de control de residuales tanto para la solución con el
algoritmo SIMPLE y PIMPLE. Los valores de residuales no son la solución, sin embargo,
según Aakti et. al [12], residuales que tengan un valor menor o igual a 10e-03 para
simulaciones en estado inestable muestran una señal de convergencia. Mientras que, para
simulaciones en estado estable el límite superior se estima en 10e-02. Los residuales se
puede decir que son variables que miden el desequilibrio o error. Por lo general, el valor
residual monitoreado es el valor global, es decir, la suma de los residuales locales de cada
nodo de la malla. La convergencia entonces se consigue cuando el valor de residual
cumple una tolerancia R ≤ e, donde e es un valor asignado conocido como tolerancia de
convergencia para un sistema de ecuaciones [37].
Por otra parte, se debe ajustar el fichero controlDict, ya que este archivo es el que controla
el paso temporal entre iteraciones de la solución numérica, y especifica el algoritmo de
solución que se desea usar. Además, posee características adicionales tales como el
control del número de Courant, el procesamiento de valores de Yplus, residuales,
Qcriterion, etc. En el ANEXO IV se muestra el contenido del fichero controlDict y su sintaxis.
A continuación, se presenta el cálculo del tiempo total requerido de la simulación en estado
inestable para que el rodete haga un recorrido de una revolución.
−1
rad 180o
(34,2716 [ ]∗ ) ∗ 360o = 0,1833 [s]
s π
(Ecuación 2.6.)
En conclusión, el fichero controlDict es el archivo que permite al usuario tener control sobre
el estudio numérico que se está realizando. En el presente estudio el paso de tiempo usado
en la simulación en estado inestable es de 1e-05 [s], y se lo denota como deltaT.
Debido a la gran cantidad de elementos que tiene la malla de la turbina tipo Francis, llevar
a cabo la simulación numérica del caso de estudio en núcleos en serie conllevaría una
saturación de procesos y una pérdida de tiempo. Es por ello que, la simulación numérica
se la configura para que realice sus cálculos de forma paralela. Para lo cual, se requieren
ciertos ficheros de texto que indiquen al ordenador cómo dividir el dominio y en cuántas
partes debe dividirlo. En primer lugar, se debe tener un fichero llamado decomposeParDict
dentro de la carpeta de system. Este archivo indica el método de descomposición del
dominio y la cantidad de subdominios que se desea; el usuario puede manipular estos
valores según el beneficio de su estudio. Para el caso de estudio se utiliza la metodología
39
de descomposición simple y se divide en 6 subdominios. El objetivo de este procedimiento
es poder dividir el dominio con un mínimo esfuerzo pero de una manera que se garantice
una solución económica en función de recursos computacionales [22].
El ordenador utilizado para el estudio es un CPU Dell OPTIPLEX 7040 Intel core i7-6700,
3.40 GHZ con memoria RAM DDR4 de 16 GB y un disco duro SATA de capacidad de 2 TB
a 7200 RPM que pertenece al Laboratorio de Mecánica-Informática ICB 207 de la Facultad
de Ingeniería Mecánica de la Escuela Politécnica Nacional. El ordenador posee ocho
núcleos físicos, sin embargo, los cálculos se los llevó a cabo usando seis núcleos. El
sistema operativo es LINUX/UBUNTU 18.04.3 LTS y el software CFD es OpenFOAM
versión 7. Adicionalmente, el post-proceso se llevó a cabo usando GNUPLOT y ParaView
5.6.3.
Finalmente, se creó un fichero ejecutable con el nombre de stk2 con los comandos
necesarios para iniciar a la simulación en paralelo. En el ANEXO V se muestra la sintaxis
de este archivo junto con la del fichero decomposeParDict.
Una vez obtenidos los resultados numéricos del presente caso de estudio, se realizó el
post-procesamiento de los mismos para ser comparados y validados con los valores de
potencia experimentales que se registran en el informe de rendimiento que realizó CELEC
EP en la central San Francisco en el año 2016, como se muestra en la Tabla 2.10.
Rendimiento
Pot E [MW] Pot H [MW]
Q [m3/s] P1 [Pa] P2 [Pa] Turbina [%]
Medido (P1-P2)*Q
Pot E/Pot H
17.0 2100000 71000 29.3 34.5 84.94
25.9 2090000 72000 40.3 52.3 77.11
33.4 2070000 75000 55.5 66.6 83.29
33.2 2060000 70000 56.2 66.1 85.06
40.3 2050000 75000 70.6 79.6 88.70
45.1 2050000 80000 79.2 88.8 89.14
50.8 2000000 80000 90.9 97.5 93.20
56 2000000 80000 100.1 107.5 93.10
59.7 1950000 80000 105.8 111.6 94.77
62.4 1950000 80000 109.5 116.7 93.84
(Fuente: [5])
40
Además, como se menciona en la Sección 2.6.2., para el post-proceso se utilizó
GNUPLOT, ya que es una poderosa herramienta para graficar resultados en ejes
coordenados, como por ejemplo graficar residuales en función del número de iteración.
Esta gráfica permite observar el comportamiento de las variables más importantes y a su
vez identificar posibles errores. Por otro lado, se obtuvo imágenes de los campos de
presión, velocidad y vorticidad usando ParaView 5.6.3.
3. RESULTADOS Y DISCUSIÓN
3.1. Mallado
La malla estructurada usada en el presente estudio es una malla mejorada a partir de la
obtenida en el estudio de Guascal & Quispe [7]. Como se ha mencionado antes, el mallado
es un elemento sensible en la simulación numérica por CFD.
41
Como se muestra en la Figura 3.1., la malla intermedia de 3,77 millones de elementos
presenta una tendencia convergente mayor frente a los otros dos casos, y por ello esta fue
la escogida para seguir con el estudio numérico propuesto. Además, debido a la
complejidad de su geometría, existen áreas en las que al refinarlas demasiado los
elementos se distorsionan a tal punto que se generan pirámides en lugar de hexaedros.
Esto genera errores en la simulación ya que el proceso de interpolación no puede realizarse
adecuadamente entre celdas contiguas. Por consiguiente, se escogió la malla intermedia
que presenta una buena calidad de elementos sin comprometer la integridad de la
simulación.
La calidad de malla se puede evaluar con base en los estadísticos obtenidos usando la
herramienta checkMesh de OpenFOAM, los mismos que se muestran en la Tabla 2.1. En
esta tabla se evidencia un valor mínimo de volumen de celda, así como también un valor
de área mínima de cara de celda, los mismos que indican que no existen volúmenes
negativos dentro del dominio. De esta manera se puede asegurar la estabilidad de malla
para el estudio numérico. Por otro lado, la calidad de malla puede ser analizada usando
parámetros adimensionales como los que se muestran en la Tabla 3.1. y Tabla 3.2.
Yplus
Dominio
k-ε k-ω SST
Rodete 72,31 73,19
Tubo de Descarga 22,7 26,23
Álabes Directrices 208,27 66,79
Caja Espiral 43,42 47,00
Álabes Fijos 43,42 47,00
(Fuente: propia)
En la Tabla 3.1., se evidencia que los valores de Yplus cumplen con los rangos establecidos
para cada uno de los modelos de turbulencia [38] [39]. Mientras que, el valor Ω para el
presente estudio (Tabla 3.2.) muestra una tendencia muy cercana a 1. El valor de Ω
implícitamente estima el recurso computacional que ocupa la simulación, es por ello que
42
este valor a medida que se acerca a 1 representa el punto de equilibrio entre tamaño de
malla y recurso computacional requerido. Por ejemplo, en el estudio de Mora [5] el número
Ω es 5 veces mayor al valor obtenido por Guascal & Quispe [7], llegando a tener un tiempo
total de simulación mucho más elevado que en el segundo caso. En el presente estudio se
obtuvo un valor Ω de 0.948, lo que implica una optimización en el tiempo y recurso
computacional.
3.2. Residuales
En la Sección de 2.6.1., se menciona que los residuales no son la solución numérica, sin
embargo, son un indicador del comportamiento convergente o divergente de la simulación.
En la Figura 3.2. y Figura 3.3., se muestran los valores residuales para cada método de
simulación.
Figura 3.2. Valores residuales para la simulación con método SIMPLE (estado estable): a) modelo
de turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
En la Figura 3.2., se puede ver que las curvas de residuales para las tres componentes de
la velocidad, en ambos casos, presentan mayor estabilidad que las demás variables. A
pesar de eso, la curva de residuales para la presión tiene la tendencia convergente más
pronunciada alcanzando valores de 10e-06. Por otro lado, en la Figura 3.3., se evidencia
un comportamiento más inestable en los valores residuales de la variable presión, sin
embargo, sus valores oscilan con una amplitud constante alrededor de 10e-05. Además,
como en el caso anterior, las componentes de la velocidad presentan estabilidad en sus
curvas durante la simulación llegando a obtener valores residuales de 10e-09.
43
Figura 3.3. Valores residuales para la simulación con método PIMPLE (estado inestable): a)
modelo de turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
Para el caso de la simulación con método SIMPLE (Figura 3.2.), con base en los resultados
obtenidos para cada variable, se estableció que la convergencia se alcanzó en la iteración
200. El rango de simulación entre la iteración 200 y 270 confirma la convergencia de los
resultados, ya que los valores obtenidos en la iteración final difieren muy poco con los
obtenidos en la iteración 200. En la Figura 3.3., la simulación con método PIMPLE recrea
una revolución del rodete que corresponde a 0,183 [s] con un paso temporal de 1e-05 [s]
para cada iteración, llegando a tener 183300 iteraciones en total. Es evidente que a partir
de 0.06 [s], en ambos casos, se muestran curvas estables de valores residuales (Figura
3.3.). Es importante señalar que en el estudio con método SIMPLE cada simulación tomó
alrededor de 55 minutos, mientras que en el caso del método PIMPLE cada simulación
tomó alrededor de 100 horas.
44
w: Velocidad angular [rad/s]
Para tener una visualización más clara del error para cada modelo de turbulencia se obtuvo
la Figura 3.5., en donde las curvas de Error [%] vs Q/Qn muestran una tendencia
decreciente a medida que se acercan al valor de caudal adimensional máximo Q/Qn = 1,07.
Estas curvas muestran un error relativo de 5,02 % y 6,49 % para el modelo de turbulencia
k–ε y k–ω SST respectivamente, cuando el caudal es el máximo. Siendo la predicción del
modelo de turbulencia k–ε la que presenta el menor error con respecto al valor
experimental. En la Figura 3.5., también se muestran las curvas de potencia [MW] vs Q/Qn,
con el fin de evidenciar la relación que existe entre el aumento de caudal y la reducción de
error en la predicción de potencia.
45
Figura 3.5. Potencia vs Q/Qn vs Error.
(Fuente: propia)
Por otra parte, se calcula la eficiencia de la turbina en función de los resultados obtenidos
de la simulación numérica, y para ello se utiliza la Ecuación 3.2.
𝑇𝑤
𝜂= [%]
(𝑝𝑖𝑛 − 𝑝𝑜𝑢𝑡 ) 𝑄
(Ecuación 3.2.)
Donde:
pin, out: Presión [Pa]
Q: Caudal [m3/s]
T: Torque [Nm]
w: Velocidad angular [rad/s]
46
Figura 3.6. Eficiencia vs Q/Qn.
(Fuente: propia)
Las Figuras 3.4., 3.5., y 3.6., fueron construidas con base en los resultados expuestos en
la Tabla 3.3. y 3.4., donde se muestran los valores obtenidos de la simulación numérica y
el error relativo con los datos experimentales.
k–ε
Experimental (Tabla 2.10.) Simulación Numérica
Pot E Torque Error Pot
Q/Qn Eficiencia Pot [MW] Eficiencia
[MW] [Nm] [%]
0,293 29,3 0,849 28863,26 0,98 96,62 0,028
0,446 40,3 0,771 225588,87 7,73 80,81 0,147
0,575 55,5 0,832 458590,53 15,71 71,68 0,235
0,572 56,2 0,850 433951,83 14,87 73,53 0,225
0,694 70,6 0,887 891418,02 30,55 56,72 0,383
0,777 79,2 0,891 1253066,2 42,94 45,77 0,483
0,875 90,9 0,932 1759705,8 60,30 35,65 0,602
0,965 100,1 0,931 2293576,9 78,60 21,47 0,731
1,029 105,8 0,947 2711054,3 92,91 12,18 0,832
1,075 109,5 0,938 3034581,8 104,00 5,022 0,891
(Fuente: propia)
47
Tabla 3.4. Resultados de simulación con modelo de turbulencia k–ω SST.
k–ω SST
Experimental (Tabla 2.10.) Simulación Numérica
Pot E Torque Error Pot
Q/Qn Eficiencia Pot [MW] Eficiencia
[MW] [Nm] [%]
0,293 29,3 0,849 126489,88 4,33 85,20 0,125
0,446 40,3 0,771 57290,062 1,96 95,12 0,037
0,575 55,5 0,832 390374,7 13,37 75,89 0,200
0,572 56,2 0,850 434840,55 14,90 73,48 0,225
0,694 70,6 0,887 905563,06 31,03 56,03 0,389
0,777 79,2 0,891 1242683,4 42,58 46,22 0,479
0,875 90,9 0,932 1751251,9 60,01 35,65 0,599
0,965 100,1 0,931 2267675,6 77,71 22,36 0,722
1,029 105,8 0,947 2671228,6 91,54 13,47 0,820
1,075 109,5 0,938 2987501,9 102,38 6,49 0,877
(Fuente: propia)
En la Figura 3.7. se presenta el perfil de velocidad a la salida del rodete, y en la Figura 3.8.
se presenta el perfil de velocidad en el plano de corte YZ. Es evidente que el perfil de
velocidad difiere dependiendo el modelo de turbulencia utilizado en la simulación. Sin
embargo, ambos gráficos muestran una tendencia similar en la formación del remolino o
vórtice a la salida del rodete. Y también, ambos casos muestran una velocidad cercana a
30 [m/s] alcanzada por el fluido dentro del rodete de la turbina.
Figura 3.7. Perfil de velocidad a la salida del rodete: a) modelo de turbulencia k-ε, b) modelo de
turbulencia k-ω SST.
(Fuente: propia)
48
Figura 3.8. Perfil de Velocidad en el plano de corte YZ, a) modelo de turbulencia k-ε, b) modelo de
turbulencia k-ω SST.
(Fuente: propia)
La Figura 3.9., muestra el perfil de presión en el plano de corte YZ. Los valores de presión
obtenidos tienen la tendencia esperada y concuerdan con los estudios previos [5]. No
obstante, los resultados obtenidos muestran una transición de presión más suave entre el
dominio de los álabes fijos, álabes directrices y rodete. Se puede observar que el perfil de
presión tiene un comportamiento similar para ambos modelos de turbulencia, con la única
diferencia en la formación del remolino o vórtice justo al momento en que el fluido sale del
rodete e ingresa al tubo de descarga. Para el modelo de turbulencia k-ε, Figura 3.9.a), la
formación del remolino presenta un menor tamaño que en el caso del modelo de
turbulencia k-ω SST, Figura 3.9.b). La presencia del vórtice en el centro del tubo de
descarga inmediatamente a la salida del rodete concuerda con estudios llevados a cabo
por Favrel et. al [40], donde se establece que el fenómeno de cavitación aparece cuando
los niveles de presión en el tubo de descarga son bajos, promoviendo la generación de un
vórtice con forma helicoidal.
49
Figura 3.9. Perfil de Presión en el plano de corte YZ, a) modelo de turbulencia k-ε, b) modelo de
turbulencia k-ω SST.
(Fuente: propia)
En la Figura 3.10. y Figura 3.12., se muestra el perfil de presión para el dominio de álabes
directrices y rodete. En la Figura 3.10., se puede ver claramente que la zona de más alta
presión está en los álabes directrices. A continuación, sigue una zona de presión media en
los álabes del rodete y una zona de presión baja en el cubo de salida del rodete. Debido a
las altas presiones presentes en los álabes directrices, generalmente es aquí donde se
encuentra la mayor cantidad de daños en el material de los álabes tal cómo se muestra en
la Figura 3.13. Sin embargo, los bordes de ataque de los álabes del rodete al ser los
elementos que se encuentran inmediatamente después de los álabes directrices presentan
también valores de alta presión, tal como se muestra en la Figura 3.12. Como se puede
observar en la Figura 3.13., los bordes de ataque de los álabes del rodete evidencian un
desgaste mucho mayor que los bordes de salida.
50
Figura 3.10. Vista superior del perfil de Presión en álabes directrices y rodete: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
Los resultados de presión para ambos modelos de turbulencia en la Figura 3.10. y Figura
3.12. concuerdan con los valores obtenidos en estudios previos [5]. Además, puede
observarse que el comportamiento de ambos modelos de turbulencia es similar en la
formación del campo de presión, con ligeras variaciones debido a la modelación
matemática de cada uno de ellos. En la Figura 3.11. se observa el perfil de presión a la
salida del rodete y entrada del tubo de descarga, donde es evidente la zona de formación
del vórtice (zona de baja presión).
Figura 3.11. Vista inferior del perfil de presión en la salida del rodete: a) modelo de turbulencia k-ε,
b) modelo de turbulencia k-ω SST.
(Fuente: propia)
51
Figura 3.12. Perfil de Presión en las paredes del rodete y en sus álabes: a) modelo de turbulencia
k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
Figura 3.13. Comparación de daños en el material con los campos de alta presión de álabes
directrices y álabes del rodete obtenidos en la simulación numérica.
(Fuente: propia)
52
Por otra parte, en la Figura 3.14. se muestra el perfil de energía cinética turbulenta en el
plano de corte YZ. La energía cinética turbulenta k se mide en [m2/s2] y tiene una relación
directa con la intensidad turbulenta que presenta un fluido, y es por ello que es una variable
importante a estudiar en problemas de dinámica de fluidos [41]. Se observa que la forma
del perfil de energía cinética turbulenta es semejante a la forma del flujo turbulento
(remolino) que se genera a la salida del rodete y entrada del tubo de descarga, tal como
se muestra en la Figura 3.17. Adicionalmente, en este caso es evidente que el perfil de
energía cinética turbulenta cambia con respecto al modelo de turbulencia utilizado en la
simulación numérica. Por ejemplo, en el caso de la simulación que usa el modelo de
turbulencia k-ε se presenta una sobreestimación de energía cinética turbulenta con
respecto al perfil obtenido en la simulación con modelo de turbulencia k-ω SST.
Figura 3.14. Perfil de Energía cinética turbulenta k en el plano de corte YZ: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
53
Figura 3.15. Vista en detalle del perfil de energía cinética turbulenta: a) modelo de turbulencia k-ε,
b) modelo de turbulencia k-ω SST.
(Fuente: propia)
Figura 3.16. Perfil de energía cinética turbulenta a la entrada del tubo de descarga: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
Una mejor visualización del perfil de energía cinética turbulenta se presenta en la Figura
3.15., la misma que es una vista en detalle de la zona comprendida por la caja espiral, los
álabes fijos, los álabes directrices, el rodete y tubo de descarga. En esta gráfica es
54
apreciable de mejor manera la cantidad de energía turbulenta producida, y se puede
observar cómo esta variable cambia dentro de todo el dominio de la turbina. A pesar de la
diferencia que existe en los perfiles de energía cinética turbulenta para ambas
simulaciones, los dos modelos de turbulencia logran captar el fenómeno turbulento
generado a la salida del rodete. En esta zona los dos modelos de turbulencia representan
de manera muy similar al remolino formado. En la Figura 3.16., se muestra el perfil de
energía cinética turbulenta en la entrada del tubo de descarga, y se puede observar que el
valor más alto se encuentra en el mismo lugar donde se genera el remolino a la salida del
rodete.
Figura 3.17. Cuerda de vórtice obtenida mediante una iso superficie de Q-criterion = 64: a) modelo
de turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
En la Figura 3.18., se muestra el perfil de presión en el plano de corte YZ. Se puede notar
que el perfil de presión al cumplir una revolución presenta valores similares en ambos
casos. También, se puede observar que la generación de vórtice en la simulación que usó
el modelo de turbulencia k-ω SST tiene un mayor tamaño que en el caso de la simulación
que usó el modelo de turbulencia k-ε. Sin embargo, el perfil del vórtice generado es mucho
menor que en la simulación con método SIMPLE. Además, en el perfil de ¼ de revolución
(i), se puede notar que el campo de presión es más pronunciado en la simulación que usó
el modelo de turbulencia k-ω SST que en el caso de la simulación que usó el modelo de
turbulencia k-ε.
Por otro lado, el perfil de presión en el caso del método PIMPLE difiere con respecto a la
simulación con el método SIMPLE en los valores máximos de presión alcanzados. No
obstante, ambas simulaciones presentan un patrón similar en la formación del campo de
presión.
56
Figura 3.18. Perfil de presión en el plano de corte YZ: a) modelo de turbulencia k-ε, b) modelo de
turbulencia k-ω SST.
(Fuente: propia)
En la Figura 3.19., se muestra el perfil de presión en los álabes del rodete, y se puede
observar que el perfil generado por las simulaciones de ambos modelos de turbulencia
presenta un patrón similar. Sin embargo, el perfil de presión en la Figura 3.19.b) presenta
un valor máximo superior al obtenido en el caso de la Figura 3.19.a).
Figura 3.19. Perfil de Presión en el rodete a) modelo de turbulencia k-ε, b) modelo de turbulencia
k-ω SST.
(Fuente: propia)
57
Figura 3.20. Vista superior del Perfil de presión en álabes directrices y rodete: a) modelo de
turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
Figura 3.21. Vista superior de perfil de Presión en álabes directrices y paredes del rodete sin
contar el cubo de salida del rodete: a) modelo de turbulencia k-ε, b) modelo de turbulencia k-ω
SST.
(Fuente: propia)
58
presión mucho más altos. Mientras que, los valores de presión en los bordes de salida de
los álabes del rodete se asemejan a los valores obtenidos en la simulación con método
SIMPLE. Es importante destacar que el perfil de presión en la Figura 3.20 y Figura 3.21.,
muestra una tendencia similar para ambos modelos de turbulencia.
Como se mencionó en la Sección 3.3, la energía cinética turbulenta tiene una relación
directa con el comportamiento turbulento del fluido. En este sentido, en la Figura 3.22. se
muestra el perfil de energía cinética turbulenta durante una revolución, donde se evidencia
el cambio del perfil de energía con respecto al uso de un modelo de turbulencia diferente.
Se puede observar en la Figura 3.22.a) que la simulación que usó el modelo de turbulencia
k-ε presenta una sobreestimación de la energía cinética turbulenta con respecto al modelo
de turbulencia k-ω SST.
El perfil de energía no presenta una zona marcada del vórtice de salida del rodete como
en el caso de la simulación con método SIMPLE. Sin embargo, las dos simulaciones
mostradas en la Figura 3.22. indican una tendencia similar de comportamiento del vórtice
formado al alcanzar una revolución, tal como se muestra en la Figura 3.24. Es importante
destacar que el flujo turbulento tiende a ir hacia afuera del eje de rotación del rodete,
alcanzando sus valores máximos en las paredes del tubo de descarga.
Figura 3.22. Perfil de Energía cinética turbulenta en el plano de corte YZ: a) modelo de turbulencia
k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
En contraste con el perfil de energía cinética, la Figura 3.23. muestra el perfil de velocidad
en el plano de corte YZ. En esta gráfica se puede notar que el cono de velocidad de salida
59
del rodete se asemeja al perfil de energía cinética turbulenta de la Figura 3.22. Pero, el
perfil de velocidad obtenido difiere con respecto a los resultados de la simulación con
método SIMPLE. Sin embargo, los valores de velocidad en el dominio de caja espiral,
álabes fijos, álabes directrices y rodete son similares a los obtenidos en la simulación con
método SIMPLE, mostrando un valor aproximado de 30 [m/s] en el rodete al cumplir una
revolución. Además, las simulaciones llevadas a cabo con ambos modelos de turbulencia
muestran una tendencia similar en el perfil obtenido y en los valores máximos de velocidad.
Figura 3.23. Perfil de velocidad en el plano de corte YZ: a) modelo de turbulencia k-ε, b) modelo de
turbulencia k-ω SST.
(Fuente: propia)
60
el caso del modelo de turbulencia k-ω SST se cumple en ambos métodos de simulación
(SIMPLE y PIMPLE). Como se pudo evidenciar en los resultados mostrados anteriormente,
los valores de presión para el caso del modelo de turbulencia k-ω SST son mayores que
los valores obtenidos con el modelo de turbulencia k-ε para las estructuras de vórtices.
Figura 3.24. Cuerda de vórtice obtenida mediante una iso superficie de Q-criterion = 64: a) modelo
de turbulencia k-ε, b) modelo de turbulencia k-ω SST.
(Fuente: propia)
En la Figura 3.25. se indica el perfil de presión en la entrada del tubo de descarga junto
con líneas de contorno de presión, con el fin de evidenciar el movimiento del vórtice hasta
cumplir una revolución. Este gráfico permite entender de mejor manera cómo se da la
formación y movimiento del vórtice, ya que se puede observar cómo el punto de presión
más baja (zona de formación de vórtice) cambia su posición en el tiempo. Esta metodología
de identificación de vórtices fue utilizada por Qian et. al [43], con el fin de entender la
influencia de la inyección de aire en la entrada del tubo de descarga sobre las pulsaciones
de presión.
En la Figura 3.25. a medida que transcurre el tiempo, el foco de formación del vórtice se
establece en la posición central de la entrada del tubo de descarga, tal y como se ve en el
tiempo (iii) y (iv). Es decir que, al cumplir una revolución el vórtice llega a ubicarse en el
centro de la entrada del tubo de descarga. Además, es importante mencionar que el perfil
61
de presión difiere con respecto al modelo de turbulencia usado. La simulación que utilizó
el modelo de turbulencia k-ω SST presenta una tendencia mucho más pronunciada del
vórtice desde ¼ de revolución (i), a comparación del perfil obtenido con la simulación que
utilizó el modelo de turbulencia k-ε.
Figura 3.25. Perfil de presión en la entrada del tubo de descarga: a) modelo de turbulencia k-ε, b)
modelo de turbulencia k-ω SST.
(Fuente: propia)
4. CONCLUSIONES
4.1. Conclusiones
Se desarrolló una metodología que permite estudiar y evaluar el comportamiento de una
turbina tipo Francis bajo condiciones de operación reales, mediante el software de
simulación CFD de código abierto OpenFOAM.
62
Con el estudio de independencia de malla se determinó que una malla intermedia de
alrededor de 3,7 millones de elementos cumple con la estabilidad y calidad de malla
necesaria para el estudio numérico propuesto presentando una buena confiabilidad de
resultados. Por otra parte, la simulación en estado inestable con 4,9 millones de elementos
tuvo un fallo integral en el cálculo computacional, debido a la distorsión de celdas en ciertos
puntos críticos de la geometría. Mientras que, la simulación con 2,9 millones de elementos
no brindaba la precisión necesaria en la predicción de resultados, especialmente en la
salida del rodete y entrada del tubo de descarga.
Se llevó a cabo dos tipos de simulaciones, la primera en estado estable usando el método
SIMPLE, y la segunda en estado inestable usando el método PIMPLE. Para el caso estable
se realizaron simulaciones para las diferentes condiciones de operación de la turbina,
mientras que para la simulación en estado inestable se optó por utilizar la condición de
operación que presentaba menos error en las predicciones del método SIMPLE.
Se validó el uso de los modelos de turbulencia k-ε y k-ω SST en simulaciones en estado
estable, obteniendo resultados de perfil de velocidad y presión similares a estudios previos.
Ambos modelos de turbulencia presentan resultados que concuerdan con la predicción de
potencia generada y la eficiencia de la turbina. Sin embargo, el modelo de turbulencia k-ε
presentó la mejor predicción de eficiencia y potencia, teniendo un error del 5,02 %.
Se determinó que la mejora de la malla en la zona del codo del tubo de descarga evitó la
generación de estructuras de vórtice en las paredes. Por lo que, discretizar una curva del
dominio con una sucesión de rectas tangentes genera zonas de inducción de resultados.
63
Se concluye que el uso de software libre y de código abierto presentó varias ventajas frente
a paquetes de software CFD comerciales. Primero, el control de la simulación mediante
archivos de texto facilita la configuración de simulaciones con distintas condiciones de
operación. Segundo, presenta una gran versatilidad para el cálculo en paralelo. Tercero,
permite entender de mejor manera cómo es la discretización de las ecuaciones que serán
resueltas en los cálculos computacionales. Por lo tanto, el software libre y de código abierto
presenta una gran ventaja para estudiar problemas de turbomaquinaria frente a paquetes
de software comercial.
Para finalizar, se concluye que el presente estudio cumplió con el objetivo principal y logró
generar una metodología válida para estudios posteriores.
No obstante, el estudio más importante sería realizar simulaciones con método PIMPLE y
PISO para tres o más revoluciones y comparar con los resultados obtenidos en el presente
estudio.
Además, se podría hacer un estudio de la influencia del ángulo en el que se encuentran los
álabes directrices, con el fin de obtener resultados numéricos para distintas condiciones de
flujo.
64
Referencias Bibliográficas
65
[13] S. L. Dixon y C. A. Hall, Fluid Mechanics and Thermodynamics of Turbomachinery,
Edición: 6. Burlington, MA: Butterworth-Heinemann, 2010.
[14] L. Stoessel, «Numerical simulations of the flow in the Francis-99 turbine», 2014,
Accedido: abr. 22, 2020. [En línea]. Disponible en:
https://odr.chalmers.se/handle/20.500.12380/204724.
[15] I. Sadrehaghighi, «Mesh Generation in CFD - A Review», CFD Open Ser., p. 262, 2017.
[16] T. J. Baker, «Mesh generation: Art or science?», Prog. Aerosp. Sci., vol. 41, n.o 1, pp.
29-63, ene. 2005, doi: 10.1016/j.paerosci.2005.02.002.
[17] K. D. H. Lee, «Grid generation for general three-dimensional configurations», ene.
1980, Accedido: abr. 22, 2020. [En línea]. Disponible en:
https://ntrs.nasa.gov/search.jsp?R=19810006200.
[18] Z. Ali y P. G. Tucker, «Multiblock Structured Mesh Generation for Turbomachinery
Flows», en Proceedings of the 22nd International Meshing Roundtable, Cham, 2014,
pp. 165-182, doi: 10.1007/978-3-319-02335-9_10.
[19] H. J. Medina, A. Beechook, J. Saul, S. Porter, S. Aleksandrova, y S. Benjamin, «Open
source Computational Fluid Dynamics using OpenFOAM», presentado en Light Aircraft
Design: Methods and Tools, 2015, Accedido: abr. 22, 2020. [En línea]. Disponible en:
https://pureportal.coventry.ac.uk/en/publications/open-source-computational-fluid-
dynamics-using-openfoam-2.
[20] H. G. Weller, G. Tabor, H. Jasak, y C. Fureby, «A tensorial approach to computational
continuum mechanics using object-oriented techniques», Comput. Phys., vol. 12, n.o 6,
pp. 620-631, nov. 1998, doi: 10.1063/1.168744.
[21] «Comparing CFD Software - Part 2: Open Source CFD Software Packages», Resolved
Analytics CFD Consulting and Simulation Strategy.
https://www.resolvedanalytics.com/theflux/comparing-cfd-software-part-2-open-
source-cfd-software-packages (accedido abr. 22, 2020).
[22] C. Greenshields, «OpenFOAM User Guide: CFD Direct, Architects of OpenFOAM»,
CFD Direct, mar. 02, 2017. https://cfd.direct/openfoam/user-guide/ (accedido abr. 22,
2020).
[23] H. Jasak y M. Beaudoin, «OpenFOAM Turbo Tools: From General Purpose CFD to
Turbomachinery Simulations», presentado en ASME-JSME-KSME 2011 Joint Fluids
Engineering Conference, may 2012, pp. 1801-1812, doi: 10.1115/AJK2011-05015.
[24] E. CELEC, «San Francisco».
https://www.celec.gob.ec/hidroagoyan/index.php/centrales/pucara/15-
centrales/sanfrancisco (accedido abr. 22, 2020).
66
[25] H. Weizhang y R. D. Russell, Adaptive Moving Mesh Methods. New York: Springer-
Verlag, 2011.
[26] V H Hidalgo, X W Luo, A Yu, y R Soto, «CAVITATING FLOW SIMULATION WITH
MESH DEVELOPMENT USING SALOME OPEN SOURCE SOFTWARE», 2014, doi:
10.13140/2.1.2423.0402.
[27] Robert. W. Fox, A. T. McDonald, y P. J. Pritchard, Introduction to Fluid Mechanics, 6ta
ed. Wiley, 2003.
[28] X. L. Liu, W. Q. Tao, y Y. L. He, «A simple method for improving the SIMPLER algorithm
for numerical simulations of incompressible fluid flow and heat transfer problems», Eng.
Comput., vol. 22, n.o 8, pp. 921-939, ene. 2005, doi: 10.1108/02644400510626488.
[29] B. Indraratna, C. Kumara, S.-P. Zhu, y S. W. Sloan, «Mathematical Modeling and
Experimental Verification of Fluid Flow through Deformable Rough Rock Joints», 2015,
doi: 10.1061/(ASCE)GM.1943-5622.0000413.
[30] N. Tonello, Y. Eude, B. de L. de Meux, y M. Ferrand, «Frozen Rotor and Sliding Mesh
Models Applied to the 3D Simulation of the Francis-99 Tokke Turbine with
Code_Saturne», ene. 2017, doi: 10.1088/1742-6596/782/1/012009.
[31] J. McNaughton, I. Afgan, D. D. Apsley, S. Rolfo, T. Stallard, y P. K. Stansby, «A simple
sliding-mesh interface procedure and its application to the CFD simulation of a tidal-
stream turbine», Int. J. Numer. Methods Fluids, vol. 74, n.o 4, pp. 250-269, 2014, doi:
10.1002/fld.3849.
[32] F. O. M. Carneiro, L. F. M. Moura, P. A. Costa Rocha, R. J. Pontes Lima, y K. A. R.
Ismail, «Application and analysis of the moving mesh algorithm AMI in a small scale
HAWT: Validation with field test’s results against the frozen rotor approach», Energy,
vol. 171, pp. 819-829, mar. 2019, doi: 10.1016/j.energy.2019.01.088.
[33] «CFD: PIMPLE Algorithm», SimScale CAE Forum, ene. 10, 2018.
https://www.simscale.com/forum/t/cfd-pimple-algorithm/81418 (accedido abr. 22,
2020).
[34] C. Introini, S. Lorenzi, A. Cammi, D. Baroli, B. Peters, y S. Bordas, «A Mass
Conservative Kalman Filter Algorithm for Computational Thermo-Fluid Dynamics»,
Materials, vol. 11, n.o 11, nov. 2018, doi: 10.3390/ma11112222.
[35] O. Winter y P. Sváček, «Numerical simulation of flow induced airfoil vibrations with large
amplitudes», en Programs and Algorithms of Numerical Mathematics 19, abr. 2019, pp.
186-194, doi: 10.21136/panm.2018.20.
[36] D. C. Wilcox, Turbulence Modeling for CFD, Edición: 3rd. Canada, Calif: D C W
Industries, 2006.
67
[37] Computational Fluid Dynamics: A Practical Approach, Edición: 3. Butterworth-
Heinemann, 2018.
[38] M. Ariff, S. M. Salim, y S. C. Cheah, «Wall y+ approach for dealing with turbulent flows
over a surface mounted cube: part 1 - low Reynolds number», Proc. Seventh Int. Conf.
CFD Miner. Process Ind. Melb. Aust., 2009, Accedido: abr. 30, 2020. [En línea].
Disponible en: https://researchportal.hw.ac.uk/en/publications/wall-y-approach-for-
dealing-with-turbulent-flows-over-a-surface-m-2.
[39] F. R. Menter, M. Kuntz, y R. Langtry, Turbulence, heat, and mass transfer 4:
proceedings of the fourth international symposium on turbulence, heat and mass
transfer, Antalya, Turkey, 12-17 October, 2003. New York: Begell House, 2003, pp.
625-632.
[40] A. T. Favrel, A. Müller, C. K. Landry, K. Yamamoto, y F. Avellan, «Study of the vortex-
induced pressure excitation source in a Francis turbine draft tube by particle image
velocimetry», 2015, doi: 10.1007/s00348-015-2085-5.
[41] R. B. Stull, «Turbulence Kinetic Energy, Stability and Scaling», en An Introduction to
Boundary Layer Meteorology, R. B. Stull, Ed. Dordrecht: Springer Netherlands, 1988,
pp. 151-195.
[42] V. Holmén, «Methods for Vortex Identification», Master’s Theses Math. Sci., 2012, doi:
http://lup.lub.lu.se/student-papers/record/3241710/file/3241720.pdf.
[43] Z. Qian, J. Yang, y W. Huai, «Numerical simulation and analysis of pressure pulsation
in francis hydraulic turbine with air admission», J. Hydrodyn. Ser B, vol. 19, n.o 4, pp.
467-472, ago. 2007, doi: 10.1016/S1001-6058(07)60141-3.
68
ANEXOS
ANEXO I.
PROCEDIMIENTO DE OBTENCIÓN DEL ARCHIVO DE MALLA
.MSH DESDE ICEM ANSYS
El procedimiento para obtener los archivos .msh desde ICEM comienza con tener el
mallado listo, es decir se debe tener abierto el archivo que se desea exportar. El siguiente
paso es abrir el submenú Output Mesh y darle clic en el ícono de una caja de herramientas
roja, tal como se muestra en la Figura A.1.
Al darle clic en este icono en la parte inferior izquierda se despliega un submenú donde se
indica el solver que se va a usar para el estudio. En este caso se debe seleccionar “ANSYS
Fluent” y darle clic en “Ok”. A continuación, se debe dar clic en el ícono de un diskette para
guardar. Entonces se acepta guardar el proyecto y en breves momentos aparece una
nueva ventana donde indica con qué formato y nombre se quiere guardar el archivo, tal
como se muestra en la Figura A.2. Para Finalizar se le da clic en el ícono “Done” y el archivo
se guarda en formato .msh para posteriormente llevarlo al entorno de OpenFOAM.
69
ANEXO II.
ARCHIVO EJECUTABLE PARA LA UNIÓN DEL DOMINIO TOTAL,
MANIPULACIÓN DEL DOMINIO DENTRO DE OPENFOAM Y
FICHERO TOPOSETDICT
Este archivo copia la carpeta rodete, y la convierte la carpeta principal en función de la cuál
las demás mallas van a ser unidas. El comando mergeMeshes cumple la función de dar la
orden para que se unan las mallas y a la vez se escriba un archivo log con la información
de todo este proceso. Una vez unidas todas las mallas, se copia la carpeta polymesh a la
carpeta llamada Final.
Los comandos que sirven para manipular la posición de los subdominios son:
transformPoints -translate ‘(x y z)’ y transformPoints -rollPitchYaw ‘(x y z)’, y dentro de los
paréntesis se escribe el vector para traslación o en el caso de querer girar la malla se
escribe la dirección sobre la cual se quiere que gire.
70
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 7 |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "mesh/component2/system";
object topoSetDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
actions
(
{
name rotorSet;
type cellSet;
action clear;
}
{
name rotorSet;
type cellSet;
action invert;
}
{
name rotor;
type cellZoneSet;
action new;
source setToCellZone;
sourceInfo
{
set rotorSet;
}
}
);
// ************************************************************************* //
ANEXO III.
FICHEROS CON INFORMACIÓN DE CONDICIONES DE BORDE Y
CONDICIÓN DE MALLA DESLIZANTE
WALL_DRAFT_TUBE //wall
{
type fixedValue;
value uniform (0 0 0);
}
WALL_GUIDE_VANE //WALL
{
type fixedValue;
value uniform (0 0 0);
}
WALL_SPIRAL_CASE //WALL
{
type fixedValue;
value uniform (0 0 0);
}
IN_SPIRAL_CASE //PATCH INLET
72
{
type fixedValue;
value uniform (0 8.82779 0);
}
OUT_DRAFT_TUBE //Patch OUTLET
{
type inletOutlet;
inletValue uniform (0 0 0);
value uniform (0 0 0);
}
#include "include/ami" //ami 00 ami01 //rotor //estator
}
// ************************************************************************* //
73
value uniform 0;
}
WALL_GUIDE_VANE //WALL
{
type nutUSpaldingWallFunction;
value uniform 0;
}
WALL_SPIRAL_CASE //WALL
{
type nutUSpaldingWallFunction;
value uniform 0;
}
IN_SPIRAL_CASE //PATCH INLET
{
type calculated;
value $internalField;
}
OUT_DRAFT_TUBE //Patch OUTLET
{
type calculated;
value $internalField;
}
#include "include/ami" // ami00 ami01 rotor estator
}
// ************************************************************************* //
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 4.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
location "0";
object epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 2 -3 0 0 0 0];
internalField uniform 0.9889084;
boundaryField
{
//- Set patchGroups for constraint patches
#includeEtc "caseDicts/setConstraintTypes"
BLADE //wall
{
type epsilonWallFunction;
value $internalField;
}
BLADE1 //wall
74
{
type epsilonWallFunction;
value $internalField;
}
WALL_RUNNER_ //wall
{
type epsilonWallFunction;
value $internalField;
}
WALL_DRAFT_TUBE //wall
{
type epsilonWallFunction;
value $internalField;
}
WALL_GUIDE_VANE //WALL
{
type epsilonWallFunction;
value $internalField;
}
WALL_SPIRAL_CASE //WALL
{
type epsilonWallFunction;
value $internalField;
}
IN_SPIRAL_CASE //PATCH INLET
{
type fixedValue;
value $internalField;
}
OUT_DRAFT_TUBE //Patch OUTLET
{
type zeroGradient;
}
#include "include/ami" // ami00 ami01 rotor estator
}
// ************************************************************************* //
75
boundaryField
{
//- Set patchGroups for constraint patches
#includeEtc "caseDicts/setConstraintTypes"
BLADE //wall
{
type kqRWallFunction;
value $internalField ;
}
BLADE1 //wall
{
type kqRWallFunction;
value $internalField ;
}
WALL_RUNNER_ //wall
{
type kqRWallFunction;
value $internalField ;
}
WALL_DRAFT_TUBE //wall
{
type kqRWallFunction;
value $internalField ;
}
WALL_GUIDE_VANE //WALL
{
type kqRWallFunction;
value $internalField ;
}
WALL_SPIRAL_CASE //WALL
{
type kqRWallFunction;
value $internalField ;
}
IN_SPIRAL_CASE //PATCH INLET
{
type fixedValue;
value $internalField ;
}
OUT_DRAFT_TUBE //Patch OUTLET
{
type zeroGradient;
}
#include "include/ami" // ami00 ami01 rotor estator
}
// ************************************************************************* //
76
FoamFile
{
version 2.0;
format binary;
class volScalarField;
location "0";
object p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;
boundaryField
{
//- Set patchGroups for constraint patches
#includeEtc "caseDicts/setConstraintTypes"
BLADE //wall
{
type zeroGradient;
}
BLADE1 //wall
{
type zeroGradient;
}
WALL_RUNNER_ //wall
{
type zeroGradient;
}
WALL_DRAFT_TUBE //wall
{
type zeroGradient;
}
WALL_GUIDE_VANE //WALL
{
type zeroGradient;
}
WALL_SPIRAL_CASE //WALL
{
type zeroGradient;
}
OUT_DRAFT_TUBE //Patch OUTLET
{
type fixedValue;
value uniform 80.088096; //PRESION NORMALIZADA
}
IN_SPIRAL_CASE //PATCH INLET
{
type zeroGradient;
#include "include/ami" //ami 00 ami01 //rotor //estator
}
// ************************************************************************* //
77
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 6
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object omega;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 0 -1 0 0 0 0];
internalField uniform 9.3997846;
boundaryField
{
//- Set patchGroups for constraint patches
#includeEtc "caseDicts/setConstraintTypes"
BLADE //wall
{
type omegaWallFunction;
value $internalField;
}
BLADE1 //wall
{
type omegaWallFunction;
value $internalField;
}
WALL_RUNNER_ //wall
{
type omegaWallFunction;
value $internalField;
}
WALL_DRAFT_TUBE //wall
{
type omegaWallFunction;
value $internalField;
}
WALL_GUIDE_VANE //WALL
{
type omegaWallFunction;
value $internalField;
}
WALL_SPIRAL_CASE //WALL
{
type omegaWallFunction;
value $internalField;
}
IN_SPIRAL_CASE //PATCH INLET
{
type fixedValue;
value $internalField;
}
OUT_DRAFT_TUBE //Patch OUTLET
78
{
type zeroGradient;
}
#include "include/ami" // ami00 ami01 rotor estator
}
// ************************************************************************* //
ANEXO IV.
FICHERO DE CONTROL DE SIMULACIÓN CONTROLDICT
El fichero controlDict es uno de los más importantes en toda la configuración del caso de
estudio en OpenFOAM, ya que este archivo posee los valores de control con respecto a
las iteraciones totales que debe cumplir el estudio. Posee además control de funciones
para el post-proceso como, por ejemplo, q-criterion, Yplus, residuales, etc. A continuación,
se muestra el contenido de este archivo.
79
\\ / A nd | Version: 6
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application simpleFoam;
startFrom latestTime;
startTime 0;
stopAt endTime;
endTime 270;
deltaT 1;
writeControl timeStep;
writeInterval 10;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
adjustTimeStep yes;
maxCo 1;
maxAlphaCo 0.5;
maxDeltaT 1e-5;
functions
{
forces
{
type forces;
libs ("libforces.so");
patches (BLADE BLADE1);
p p;
U U;
rho rhoInf;
rhoInf 998.9;
CofR (0.0016646385192871094 0.013781309127807617 0.2870127260684967);
}
#includeFunc residuals
#includeFunc Q
#includeFunc yPlus
#includeFunc wallShearStress
}
// ************************************************************************* //
80
ANEXO V.
ARCHIVO EJECUTABLE PARA INICIAR LA SIMULACIÓN Y
FICHERO USADO PARA EL CÁLCULO EN PARALELO
81