Escuela Politécnica Nacional: Facultad de Ingeniería Mecánica

Descargar como pdf o txt
Descargar como pdf o txt
Está en la página 1de 93

ESCUELA POLITÉCNICA

NACIONAL

FACULTAD DE INGENIERÍA MECÁNICA

“ESTUDIO DE MALLADO ESTRUCTURADO 3D ROTATIVO EN


OPENFOAM PARA APLICACIONES EN TURBINAS TIPO
FRANCIS”

TRABAJO DE TITULACIÓN PREVIO A LA OBTENCIÓN DEL TÍTULO DE


INGENIERO MECÁNICO

VELASCO BETANCOURT MARTIN RICARDO


[email protected]

DIRECTOR:
Ing. HIDALGO DÍAZ VÍCTOR HUGO, D.Sc.
[email protected]

CODIRECTOR:
Ing. VALENCIA TORRES ESTEBAN ALEJANDRO, Ph.D.
[email protected]

Quito, mayo 2020


CERTIFICACIÓN

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:

Enfrentaste al mundo entero por mí, me diste amor y


esperanza más allá de lo que algún día pudiese describir con
palabras. Siempre fuerte y con tu semblante firme me
enseñaste todo lo necesario para seguir adelante. Este podría
ser uno de los momentos más importantes en mi vida, sin
embargo, no es exclusivo para mí. Este momento te
pertenece, porque sin ti nada de esto hubiera sido posible.
Gracias, mamá.

A mi hermano, por ser mi fuente de fortaleza y cariño


incondicional, porque gracias a ti tengo una razón más para
buscar ser mejor cada día.

A mis abuelos, por su absoluto respaldo en todo el transcurso


de mi vida, porque más que abuelos han hecho el papel de
padres.

A mi familia.

“El ruido siempre es relativo al silencio que le precede. Entre


más absoluto sea el silencio, más impactante será el
estruendo”
V.

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.

A mi director, Víctor Hugo Hidalgo Díaz, por brindarme su


apoyo incondicional, por enseñarme el valor del trabajo en
equipo, por crear oportunidades para mi crecimiento
profesional y personal; en fin, por ser el mentor que todo
estudiante espera encontrar en su camino.

A mi codirector, Esteban Alejandro Valencia Torres, por su


colaboración y guía para culminar con éxito el trabajo de
investigación.

A mi profesor, Willan Leopoldo Monar Monar, por ser mi guía


y apoyo durante toda mi vida universitaria, por sus
enseñanzas y consejos.

Al Centro de investigación y recuperación de turbinas


hidráulicas y partes industriales (CIRT) y al proyecto de
investigación PIJ-17 13, por haber facilitado información de
gran relevancia para la realización del trabajo de
investigación.

A mis amigos, quienes nunca escatimaron en dar la mano


para salir adelante todos como grupo, por las risas, los
momentos buenos y los complicados, por hacer de la
universidad un lugar más ameno ¡Son grandes, Pros! A mis
amigos del laboratorio quienes fueron fundamentales para
culminar este trabajo; además, llegaron a formar parte de mi
día a día, me enseñaron muchísimo y me apoyaron en esos
momentos de duda, los admiro y aprecio mucho.

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.

Palabras clave: algoritmo, cálculo en paralelo, malla estructurada, modelos de turbulencia,


OpenFOAM.

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.

Keywords: algorithm, OpenFOAM, parallel calculation, structured mesh, turbulence


models.

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.

Tomando en consideración estudios previos y el software utilizado en los mismos, el


presente trabajo pretende desarrollar y simular un mallado estructurado rotativo 3D en una
turbina tipo Francis con características similares a la turbina N°2 de la central hidroeléctrica
San Francisco usando el software libre y de código abierto OpenFOAM. La finalidad del
estudio es aportar con un desarrollo tecnológico de vanguardia y sin restricción comercial
que permita la optimización del diseño de turbinas hidráulicas para las condiciones hídricas
del Ecuador.

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.

Figura 1.1. Turbina Francis de la central hidroeléctrica San Francisco.


(Fuente: propia)

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].

Las turbinas hidráulicas tradicionalmente son diseñadas para un punto de operación


estable donde la eficiencia es máxima, y es conocido como BEP (Best Efficiency Point).
Este punto se caracteriza por trabajar con un flujo estable con mínima vorticidad a la salida
de la turbina, resultando en condiciones dinámicas favorables de operación con pocas
inestabilidades de flujo y pocas vibraciones en el rodete [14]. Sin embargo, durante la
condición de operación estable, cuando el rodete gira a una velocidad angular constante,
un intercambio inestable de momento toma lugar entre el borde de salida de los álabes
directrices y el borde de ingreso de los álabes del rodete. Este complejo fenómeno propaga
una onda de sonido en todas las direcciones de la turbina, y toma el nombre de interacción
rotor-estator (RSI). Se sabe que las interacciones resultantes por una distribución no
uniforme de flujo entre los álabes directrices y los álabes del rodete generan fluctuaciones
de presión con diferentes frecuencias en el dominio estático y también en el rotatorio. La
frecuencia de RSI depende principalmente de tres parámetros: el número de álabes
directrices, número de álabes del rodete y de la velocidad angular del rodete [10].

1.3. Técnicas de mallado estructurado 3D en turbinas hidráulicas


La generación de la malla o también conocido como proceso de discretización de dominio,
ha evolucionado considerablemente, hasta el punto de poder cubrir cualquier tipo de
geometría compleja con elementos hexaédricos, tetraédricos, entre otros. La malla puede
ser estructurada o no estructurada. El mallado estructurado presenta igual número de
elementos adyacentes y nodos internos, mientras que el mallado no estructurado no sigue
un patrón establecido. La malla estructurada generada presenta elementos hexaédricos
(3D) y cuadrados (2D), y debido a que ambos tipos de elementos tienen buena calidad y
un adecuado comportamiento numérico, son ampliamente usados en trabajos de
investigación que buscan aumentar su precisión en la predicción de resultados [15]. Por
ejemplo, se tiene evidencia que las predicciones de arrastre (Drag) basadas en las
ecuaciones RANS (Reynolds Average Navier-Stokes) usando un mallado estructurado
(hexaédrico) tienen una estrecha aproximación a los valores experimentales, mucho más
que las predicciones usando mallado no estructurado o mallado híbrido [16].

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.

Figura 1.2. Técnicas de mallado estructurado.


(Fuente: [16])

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].

Figura 1.3. Discretización del dominio en Multibloques.


(Fuente: [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

1.4.1. Dinámica de Fluidos Computacional (CFD) y sus etapas

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.

Figura 1.6. Diagrama de flujo del proceso de simulación CFD


(Fuente: [19])

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.

Tabla 1.1. Comparación entre modelos DNS, LES, RANS y DES.

Modelos de turbulencia Ventajas Desventajas

DNS (Direct Numerical Más preciso. No necesita Requiere un alto costo


Simulation) correlaciones empíricas. computacional (capacidad).
Capaz de caracterizar todos Difícil inclusión de
los detalles del flujo. condiciones de borde
precisas para aplicaciones
de ingeniería.

LES (Large Eddy Es capaz de capturar la Aún costoso en términos


SImulation) dinámica de los vórtices computacionales. Algunas
dominantes en el sistema. dificultades en la
Relativamente más representación del flujo en
económico que DNS. Más geometría complejas.
preciso que RANS.

RANS (Reynolds Average Adecuado para problemas Es incapaz de capturar los


Navier-Stokes) de ingeniería. El costo detalles del flujo. Depende
computacional es bajo. en gran medida de
correlaciones empíricas.

DES (Detached Eddy Adecuado para problemas Es incapaz de capturar los


Simulation) de ingeniería. Captura la detalles de flujo en regiones
inestabilidad de flujos cercanas a las paredes.
separados. Más aplicable
que RANS.

(Fuente: [8])

En la segunda etapa que corresponde al Cálculo Computacional, como su nombre lo indica,


se realizan los cálculos computacionales, y para ello el usuario debe seleccionar el solver

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

Bajo o ningún costo. Muchas veces el software es distribuido


con complejas licencias.

Código fuente está disponible. Código fuente está sujeto a cambios.

Permite ser personalizado por el usuario. Curva de aprendizaje empinada.

Soporte para el usuario mediante la Pueden aplicarse pagos para recibir


comunidad. soporte técnico.

Las simulaciones son típicamente Pre-proceso y la configuración inicial del


configuradas usando archivos de texto que caso pueden presentar un desafío al
permiten rápidamente vincular los principio (Puede o no puede disponer de
parámetros de la simulación con una Interfaz gráfica de usuario).
aplicaciones de terceros, p. ej.
herramientas de optimización.

(Fuente: [19])

Finalmente, la tercera etapa o Post-proceso involucra el tratamiento, análisis y validación


de resultados. El estudio de validación debe estar siempre inmerso en el proceso con el
objetivo de evaluar las capacidades predictivas del código, la validez de las
consideraciones de simplificación y la precisión y aplicabilidad del modelo de turbulencia.
Esta etapa es de gran importancia ya que permite generar la confiabilidad necesaria para
interpretar los resultados obtenidos en la simulación por CFD. Además, como parte de esta
etapa, el usuario puede usar los resultados de la simulación para calcular variables de
interés (coeficiente de arrastre, coeficiente de sustentación, fuerzas, torque) que pueden
ser usados como parámetros de entrada en un proceso de diseño. Ya que este proceso no

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.

1.4.2. OpenFOAM como Software libre y de código abierto de Simulación CFD

Las áreas de matemática y mecánica computacionales proveen métodos fundamentales y


herramientas para simular procesos físicos que son difíciles o muy costosos de realizarlos
en un laboratorio. Por más de treinta años se ha reconocido que la ciencia computacional
constituye una tercera rama de la ciencia, teniendo la misma jerarquía que la rama teórica
y experimental. Atravesando diferentes disciplinas al centro de la ciencia computacional, la
dinámica de fluidos constituye el núcleo de OpenFOAM [8].

OpenFOAM nació en el Imperial College de Londres, que es un centro de investigación en


CFD desde 1960. El desarrollo original empezó con el Prof. David Gosman y Dr. Radd Issa
junto con Henry Weller y el Dr. Hrvoje Jasak [20] como los principales programadores. Se
basaron en el método de volúmenes finitos (FVM), y en el enfoque de programación de
objetos orientados (object-oriented programming) en C++ para el desarrollo de un modelo
de mimetismo de ecuaciones y operaciones escalares-vectoriales-tensoriales. En 2004
Henry Weller lanzó la licencia pública general GNU del software OpenFOAM tras fundar su
compañía OpenCFD Ltd. OpenFOAM constituye una librería de herramientas CFD
desarrollada en C++ con solvers personalizados (más de sesenta) que pueden realizar
simulaciones en problemas de combustión, modelos de turbulencia, transferencia de calor,
flujo multifásico y turbomáquinas [8]. OpenFOAM es originalmente distribuido en sistemas
operativos LINUX y puede ser descargado desde la página web de OpenFOAM
Foundation.

Actualmente, existen varios paquetes de software de código abierto enfocados en la


resolución de problemas de dinámica de fluidos, sin embargo, el más usado y con mayor
soporte es OpenFOAM (OF), y su acrónimo significa Open Source Field Operation and
Manipulation. En los últimos años OpenFOAM ha tenido un vasto crecimiento a nivel
industrial y académico debido a su considerable credibilidad obtenida por estudios
validados en varios países y laboratorios de investigación. Un alto índice de universidades
y corporaciones han empezado a usar OpenFOAM, tanto de forma independiente o en
conjunto con paquetes comerciales. Empresas de gran renombre como Mercedes Benz,
BASF, BMW, Volkswagen e Intel tienen actualmente presencia en las Conferencias de
Usuarios OF [21].

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].

1.5.1. Herramientas y servicios útiles

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.

Tabla 1.3. Herramientas (utilities) útiles en el proceso de simulación CFD.

Pre-proceso

checkMesh Realiza pruebas sobre la malla y reporta los indicadores


estadísticos de calidad de malla

applyBoundaryLayer Genera un perfil de capa límite simple para los campos de


velocidad y turbulencia

Generación de malla

blockMesh Herramienta de generación de bloques y mallas base, se


recomienda para geometrías simples

snappyHexMesh Herramienta de generación automática de malla que


encaja elementos hexaédricos a la superficie que describe
la geometría. Es adecuada para geometrías complejas

Conversión de malla

ideasUnvToFoam Convierte mallas guardadas con la extensión .unv a


formato OpenFOAM, p. ej. Archivos de malla generados
en SALOME

netgenNeutralToFoam Convierte mallas guardadas en formatos de archivo


neutral generados en Netgen

14
fluentMeshToFoam Convierte mallas guardadas en extensión .msh a formato
OpenFOAM, p. ej. Archivos de malla generados en ICEM
ANSYS FLUENT

Post-proceso

yPlusRAS Calcula el valor de y+ para todos los dominios


geométricos que tienen como condición de borde wall

wallShearStress Calcula el valor del esfuerzo cortante para todos los


dominios geométricos que tienen como condición de
borde wall

(Fuente:[19])

1.5.2. Solvers disponibles

Un solver es un ejecutable que provee de un algoritmo para encontrar una solución


específica a un conjunto de ecuaciones. Hay varios solvers disponibles en OpenFOAM, sin
embargo, debido al enfoque del presente estudio sólo se tratarán los que trabajan con
fluidos incompresibles [19]. En la Tabla 1.4. se presenta algunos ejemplos con sus
características más importantes.

Tabla 1.4. Solvers disponibles en OpenFOAM para flujos incompresibles.

Solver Descripción

icoFoam Solver para comportamiento transitorio de flujos laminares


incompresibles. Puede ser aplicado a fenómenos con bajo número de
Reynolds, p. ej. Sistemas de flujo interno, UAVs y MAVs.

simpleFoam Solver para comportamiento en estado estable de flujos incompresibles


usando el algoritmo SIMPLE. Puede ser usado tanto para flujos
laminares o turbulentos (con el uso de RANS como modelo de
turbulencia).

pisoFoam Solver para flujos transitorios incompresibles usando el algoritmo


PISO. Puede ser usado en conjunto con los modelos de turbulencia
LES y RANS.

pimpleFoam Solver para flujos incompresibles usando el algoritmo PIMPLE (una


combinación de PISO y SIMPLE que permite mantener la estabilidad
numérica mientras se usan pasos de tiempo grandes). Se puede usar
en conjunto con los modelos de turbulencias LES y RANS.

(Fuente: [19])

15
1.5.3. Estructura y ajuste de un caso de simulación

La estructura de un caso de simulación, su ajuste y configuración están directamente


relacionados con el sistema de ecuaciones y las regiones del dominio que representan el
sistema físico que se busca modelar [19]. En la Figura 1.7. se muestra una estructura base
para cualquier caso de estudio. La carpeta general del caso posee tres subcarpetas, las
mismas que contienen la información de los parámetros de control de la simulación, las
condiciones de borde e iniciales, al igual que las propiedades de transporte y la información
del dominio geométrico y la malla generada.

Figura 1.7. Estructura del directorio de un caso de estudio en OpenFOAM.


(Fuente: [8])

1.5.4. OpenFOAM Turbo Tools

El grupo desarrollador de OpenFOAM notó que el área de turbomáquinas requería


herramientas más específicas para solventar ciertos problemas en la simulación CFD, por
lo que se implementaron varias utilites que cumplían con este objetivo. Como se muestra
en el trabajo de Jasak H. y Beaudoin M. [23], la simulación de turbomaquinaria representa

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.

La resolución del mallado es un requerimiento muy importante en el estudio de flujo en


turbomáquinas debido a la presencia de condiciones de presión adversa en capa límite y
efectos transitorios en la estela, que resultan en un ambiente desafiante para una
simulación CFD. Por otra parte, el diseño de máquinas rotativas involucra álabes que no
tienen un mismo punto de cabeceo (pitch) en su línea consecutiva, descartando así la
posibilidad de simplificar la geometría con una imposición de restricciones simétricas. En
respuesta a esta problemática, se han aplicado algunas simplificaciones con respecto al
flujo y a la geometría, tomando el nombre de Turbo Tools. El comportamiento rotacional
puede ser manejado de varias formas:

• Cálculo de mallas móviles (Moving mesh calculations), un componente de la malla


rotativa es transformado en movimiento de un cuerpo sólido por el movimiento de
vértices, y la presencia de rotación aparece en las ecuaciones por medio del término
de velocidad relativa. Sin embargo, este enfoque no es apropiado para
simulaciones de flujo transitorio.

• Geometrías bajo un régimen de rotación de cuerpo sólido y condiciones


estacionarias pueden ser simuladas al convertir el grupo de ecuaciones a un único
marco de referencia rotacional (single rotating frame of reference), representado
como una malla estática. Aquí el campo de velocidades que aparece en las
ecuaciones es la velocidad relativa, y la rotación se obtiene tomando en cuenta la
fuerza centrífuga y de Coriolis en términos de la ecuación de momentum.

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.

En la interfaz entre la parte estacionaria y rotacional surgen ciertos problemas, y existen


algunas técnicas para tratarlos. Para simulaciones de flujo transitorio con cambios
topológicos en su malla, la técnica a usar en su interfaz debe ser malla deslizante (sliding
mesh). Las técnicas de manejo de superficies de interfaz se listan a continuación:

• La manera más simple de tratar un problema de interfaz rotor-estator es la técnica


rotor congelado (frozen rotor), donde la sección del estator y rotor son simuladas
en una posición relativa preestablecida en estado estacionario, usando el método
MRF (Multiple Reference Frame). El comportamiento transitorio se infiere (con
aproximación) a partir de los resultados en estado estacionario. Para incrementar
la confiabilidad de esta simulación, una simulación de rotor congelado se repite
varias veces, variando la posición relativa de los elementos.

• Interfaz general de malla GGI (General Grid Interface), representa el


emparejamiento de dos dominios de manera implícita, sin ningún tipo de cambio en
la topología de malla.

• 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.

• Interfaz de superposición parcial (Partial overlap interface), es una generalización


de GGI que permite simular secciones no emparejadas geométricamente de un
estator y rotor, asumiendo su periodicidad independiente.

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.

Figura 2.1. Estructura Metodológica del caso de estudio


(Fuente: propia)

Como se muestra en la Figura 2.1. la estrategia metodológica general para el presente


estudio comienza con el desarrollo del modelo 3D de la turbina tipo Francis en el entorno
de ICEM, con el fin de obtener un mallado estructurado en todo el dominio geométrico. En
este entorno de trabajo se puede modificar tanto los puntos, líneas y superficies del modelo
con el fin de generar un modelo con transiciones suaves entre dominios. En seguida, se

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.

2.1. Modelo 3D del caso de estudio

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.

2.2.1. Modelo Geométrico

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])

En el presente estudio se tomó como punto de partida la geometría optimizada de Guascal


& Quispe [7], la misma que fue tratada en el entorno de ICEM. A pesar de ser una geometría
simplificada y con mejores características que la obtenida en el estudio de Mora [5], aún
presentaba problemas que en el presente estudio se abordaron con el afán de mejorar aún
más la geometría original, esto se puede ver en la Figura 2.4.

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)

2.2.2. Mallado Estructurado

Una vez mejorado el modelo geométrico, el siguiente paso corresponde a la generación


del mallado estructurado en todo el dominio, y para ello se utiliza la metodología de
multibloques como se puede ver en la Figura 2.6.

Figura 2.6. Metodología de multibloques.


(Fuente: [7])

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)

2.3. Importación y tratamiento del modelo 3D en OpenFOAM

Al tener cada subdominio de la turbina con su respectivo mallado estructurado, se procede


al ensamblaje final. Sin embargo, este proceso será llevado a cabo fuera del entorno de
ICEM, y para ello se necesita tener un archivo con extensión .msh de cada uno de los
subdominios. En el ANEXO I se muestra el procedimiento llevado a cabo para obtener
estos archivos.

En el presente estudio se utilizó OpenFOAM versión 7 junto con la distribución GNU/LINUX


Ubuntu 18.04.3 LTS. En primer lugar, se debe crear un directorio del caso de estudio como
se indica en la Figura 1.7., siendo las carpetas 0.orig, constant y system las que contienen
todos los archivos necesarios para generar la simulación. En la Figura 2.8. se muestra la
estructura del directorio del caso de estudio.

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)

Agregando a lo anterior, dentro de system en la carpeta de información del rodete, se debe


agregar el fichero topoSetDict. Este archivo contiene la información de los sets de celdas
que serán creados para distinguir entre dominio rotatorio y estacionario, y por ser el rodete
el elemento más importante, este debe contener en su carpeta esta información. A
continuación, se abre una terminal en la carpeta del rodete y se ejecuta el comando
topoSet, con el que se establece que el rodete será la malla esclava mientras que el resto
de la malla será la malla maestra.

En seguida, se utiliza un archivo ejecutable llamado makemergeMesh-sh, que al ser


ejecutado desde la terminal crea una carpeta llamada Final; y luego, dentro de esta carpeta
une todos los subdominios y crea la subcarpeta polyMesh con la información de mallado
del dominio total. Por último, se ejecutan los comandos setsToZones y setSet con el fin de
convertir los conjuntos de celdas en sets de zonas, ya que sin esto no se puede asignar

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.11. Dominio computacional total.


(Fuente: propia)

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

Como se menciona en el primer capítulo, la convergencia de resultados depende en gran


medida de la calidad de malla. Por lo que, verificar la calidad de malla después de haber
sido importada y ensamblada en OpenFOAM es de vital importancia para el presente
estudio. Para eso se ejecuta el comando checkMesh en la terminal, el mismo que muestra
las características más importantes con respecto a la calidad de malla y el tipo de
elementos que tiene [22]. La Tabla 2.1. muestra las características de la malla del dominio
computacional total.

Tabla 2.1. Estadísticos de la malla.

Tipo de Elementos Hexaédricos


Número de Elementos 3,7749e+06
Número de Nodos 3,9814e+06
Número de Parches 14
Max. Relación de aspecto 76,6084
Min. Área de superficie 8,6131e-06 [m2]
Min. Volumen 1,4656e-07 [m3]
Max. Oblicuidad 7,5778
(Fuente: propia)

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.

2.5. Condiciones de simulación del modelo con metodología de


interfaz de malla móvil

2.5.1. Condiciones de operación

Las condiciones de operación de la turbina tipo Francis de la Central Hidroeléctrica San


Francisco se las puede obtener de la página web de la Corporación Eléctrica del Ecuador
CELEC [24] y del estudio de Mora [5]. En la Tabla 2.2. se resume las condiciones más
importantes.

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

Temperatura del agua 16 [°C]


Velocidad de rotación 34,2716 [rad/s]
Densidad 998,9 [kg/m3] [27]
Viscosidad cinemática 1,106e -06 [m2/s] [27]
Número de Álabes
Álabes Fijos 20
Álabes Directrices 20
Álabes Rodete 13
(Fuente: propia)

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.)

Lc: Longitud característica [m].


U: Velocidad de free stream [m/s].
Ʋ: Viscosidad cinemática [m2/s].

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.

En soluciones numéricas para problemas de flujo incompresible, el enfoque de corrección


de presión es el método más utilizado en CFD [28]. El algoritmo SIMPLE (Semi-implicit
method for pressure-linked equations) propuesto por Patankar y Spalding en 1972 sigue
este enfoque, y puede ser entendido mediante el diagrama de flujo mostrado en la Figura
2.13. Además, cómo se indica en la Tabla 1.4., este algoritmo es usado principalmente
para la resolución de problemas en estado estable.

Figura 2.13. Algoritmo SIMPLE


(Fuente: [29])

En estado estable o estacionario, la simulación no presenta un movimiento relativo de malla


entre dos regiones, debido a que en este estado no existe un cambio de las propiedades
ni de la posición de sus subdominios en el tiempo. Es por eso, que se utiliza el enfoque de
rotor congelado (frozen rotor), el mismo que sigue la metodología MRF (Multiple Reference
Frame) [30]. La metodología MRF establece un marco de referencia separado para cada
región, tanto para la rotacional como para la estática. En este sentido, al rodete se lo
considera como un marco de referencia rotacional y a las zonas de álabes directrices como
un marco de referencia estacionario, así como también al tubo de descarga, los álabes fijos
y la caja espiral. En el marco de referencia rotacional la velocidad relativa de Coriolis es
agregada a las ecuaciones que gobiernan el fenómeno y el flujo es calculado a partir de la
velocidad relativa [14] [11].

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].

Tabla 2.3. Condiciones de borde globales.

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.

Tabla 2.4. Condiciones de borde.

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)

Tabla 2.5. Condiciones de borde

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.

Tabla 2.6. Superficies de transferencia de datos

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.

2.5.3. Metodología de interfaz de malla móvil o dinámica

La aproximación que se adapta mejor para el caso de estudio es el mallado corredizo o


deslizante (sliding mesh), en el cual existe un movimiento deslizante relativo entre dos
dominios. La idea básica de esta metodología es tratar a la interfase como una condición
de borde interna para cada región individual. Es decir, si existen dos regiones de celdas
que llegan a un borde en común (interfase), y a pesar de que ambas zonas tienen distinto
tamaño de elementos y número de nodos, la interfase se comportará como una condición
de borde de superficie para cada una de ellas; y, esta superficie representa una zona de
acoplamiento entre ambas regiones. El acoplamiento no es fijo, ya que se adapta y se
actualiza a medida transcurre cada paso de tiempo [31]. Esta metodología se utiliza en
estudios no estacionarios.

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.

Figura 2.14. Interfaces y superficies de transferencia.


(Fuente: propia)

34
2.5.4. Simulación con Método PIMPLE: Condiciones iniciales y de borde

El algoritmo PIMPLE (Pressure-implicit split operator) como se indica en la Tabla 1.4., es


una combinación de SIMPLE y PISO, especialmente usado en problemas de casos no
estacionarios. La mejor manera de entender el algoritmo PIMPLE es imaginarlo como si
fuese el algoritmo SIMPLE para cada paso de tiempo, donde a los correctores externos de
las ecuaciones (outer correctors) se los toma como iteraciones, y una vez que cumpla la
convergencia pasa al siguiente tiempo [33]. Tanto el algoritmo PISO como PIMPLE siguen
el mismo procedimiento de discretización de las ecuaciones, sin embargo, PIMPLE se
caracteriza por tener un bucle iterativo interno adicional que resuelve la no linealidad de las
ecuaciones de Navier-Stokes y la ecuación de la energía usando iteraciones de Picard.
Este bucle iterativo interno lo convierte en un método implícito, mientras que PISO es un
algoritmo semi implícito, ya que la corrección de la velocidad se la resuelve de forma
explícita. El proceso iterativo de PIMPLE permite realizar simulaciones con factores de
relajación que mejoran la estabilidad y la adaptabilidad de los pasos de tiempo a medida
que progresa la simulación [34]. En la Figura 2.15. se muestra el diagrama de flujo del
algoritmo PIMPLE.

Figura 2. 15. Algoritmo PIMPLE.


(Fuente: [35])

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.

2.5.5. Modelo de turbulencia

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.

Tabla 2.7. Coeficientes del modelo de turbulencia que usa OpenFOAM

Modelo de turbulencia k-ε Modelo de turbulencia k-ω SST


Cµ 0,09 α k1 0,85
C1 1,44 α k2 1
C2 1,92 α ω1 0,5
C3 0 α ω2 0,856
γ1 0,5555
γ2 0,44
β1 0,075
β2 0,0828
β* 0,09
a1 0,31
c1 10
(Fuente: propia)

2.6. Simulación en OpenFOAM

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].

En el presente estudio se usó el algoritmo SIMPLE y PIMPLE para su respectivo campo de


aplicación. En la Tabla 2.8., se muestran los ajustes del solver tanto para el caso de rotor
congelado y malla deslizante, los mismos que se encuentran en el archivo fvSolution dentro
de system.

37
Tabla 2.8. Ajustes del Solver

Metodología de Rotor congelado (frozen rotor)


Variable Solver Tolerance Rel Tol Smoother
p GAMG 1e-06 0,05 GaussSeidel
U smoothSolver 1e-07 0,1 GaussSeidel
k smoothSolver 1e-07 0,1 GaussSeidel
epsilon smoothSolver 1e-07 0,1 GaussSeidel
omega smoothSolver 1e-07 0,1 GaussSeidel
Metodología de Malla deslizante (sliding mesh)
Variable Solver Tolerance Rel Tol Smoother
Phi GAMG 1e-06 0,01 DIC
P GAMG 1e-06 0,05 GaussSeidel
U |k| omega smoothSolver 1e-06 0,1 symGaussSeidel
U |k| omega Final smoothSolver 1e-06 0 symGaussSeidel
U |k| epsilon smoothSolver 1e-06 0,1 symGaussSeidel
U |k| epsilon Final smoothSolver 1e-06 0 symGaussSeidel
(Fuente: propia)

En la metodología de malla deslizante se agrega la variable Phi para implementar el utility


de potentialFoam. El mismo que sirve para generar un campo inicial conservativo de la
variable U para casos complejos. En otras palabras, se convierte en una herramienta útil
para evitar inestabilidades debido a que crea un campo conservativo a partir de un campo
no conservativo dado por el usuario [22]. Por otra parte, los parámetros Tolerance y relTol
de la Tabla 2.8. se usan como medidas de convergencia, ya que Tolerance especifica el
valor al cual el residual inicial debe llegar, mientras que relTol especifica una relación del
valor residual actual con el valor residual final de referencia (Tabla 2.9.) que usa el solver
como medida de convergencia [11].

Tabla 2.9. Control de residuales y factores de relajación.

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.

2.6.2. Cálculo en Paralelo y características de Software/Hardware

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.

2.7. Post-proceso y Validación del modelo con resultados de


estudios previos

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.

Tabla 2. 10. Valores experimentales de rendimiento a diferentes condiciones de flujo volumétrico.

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.

3.1.1. Independencia de Malla

Debido a que la simulación recrea fenómenos complejos de flujo, es necesario comprobar


que no existe inducción de resultados debido a la malla. Por lo cual, se utilizó tres distintas
mallas para estudiar la independencia de resultados que se obtenía con cada una de ellas.
Este estudio permite encontrar el límite en el cual las predicciones de variables son
independientes con respecto al aumento o disminución del número total de elementos del
dominio. En la Figura 3.1., se muestran los valores residuales de presión en función del
número de iteración para cada una de las mallas.

Figura 3.1. Valores residuales de presión en distintas mallas.


(Fuente: propia)

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.

3.1.2. Calidad de Malla

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.

Tabla 3.1. Valores de Yplus para cada subdominio.

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)

Tabla 3.2. Número Ω para el dominio total.

N° elementos [NE] N° Nodos [ND] Ω [NE/ND]


3774986 3981484 0,948
2926033 2727438 1,072 [7]
68456549 13120572 5,211 [5]
(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.

3.3. Simulación con Método SIMPLE: Estado estable


Gracias a la información experimental (Tabla 2.10.) que se tenía del estudio de Mora [5],
se realizó una simulación en estado estable para cada una de las condiciones de caudal.
Obteniendo resultados de potencia [MW], torque [Nm], eficiencia [%] y el error de las
predicciones por CFD [%] con respecto a los valores medidos in situ. Para el cálculo de la
potencia se utilizó la Ecuación 3.1.
P = T w [MW]
(Ecuación 3.1.)
Donde:
T: Torque [Nm]

44
w: Velocidad angular [rad/s]

En la Figura 3.4., se muestra una gráfica de la potencia en función de un caudal


adimensional. La variable de caudal adimensional no es más que la división del caudal de
trabajo (Tabla 2.10.) por el caudal nominal de 58 [m3/s]. Como se puede observar en la
Figura 3.4., los valores de potencia obtenidos para cada modelo de turbulencia (Tabla 3.3.
y Tabla 3.4.) muestran una tendencia similar, a medida que las curvas se acercan al valor
de caudal adimensional máximo Q/Qn=1.07 las barras presentan menor error con respecto
al valor experimental. Siendo las predicciones del modelo de turbulencia k–ε las de mayor
precisión.

Figura 3.4. Potencia vs Q/Qn.


(Fuente: propia)

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]

En la Figura 3.6., se muestra la gráfica de Eficiencia vs Q/Qn, donde es evidente una


tendencia similar en ambos modelos de turbulencia. Y, al igual que en el caso de las curvas
de potencia, las curvas de Eficiencia a medida que se acercan al valor de caudal
adimensional máximo Q/Qn = 1,07, sus barras presentan menor error. El modelo de
turbulencia k–ε es el que mayor aproximación tiene al valor de eficiencia experimental.

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.

Tabla 3.3. Resultados de simulación con el modelo de turbulencia k–ε.

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)

Por lo tanto, debido a lo expuesto anteriormente se escogió los resultados de la simulación


numérica con caudal Qmax=62,4 [m3/s] (Tabla 2.10.) para realizar un análisis de campos de
presión, energía, velocidad y vorticidad en todo el dominio de la turbina tipo Francis.

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.17. se muestra la formación del vórtice o remolino inmediatamente después


de que el fluido sale del dominio del rodete e ingresa al tubo de descarga. Es importante
destacar que para la identificación del remolino a la salida del rodete se utilizó el método
de Q-criterion, el mismo que define a los vórtices como zonas donde la magnitud de
55
vorticidad es mayor que la magnitud del esfuerzo cortante [42]. Como se mencionó antes,
la energía cinética turbulenta tiene una relación directa con el comportamiento turbulento
del fluido, y es por ello que la forma del vórtice y del perfil de energía cinética es similar.
Por otro lado, se puede ver que el vórtice generado por la simulación que usó el modelo de
turbulencia k-ε tiene un menor tamaño con respecto al vórtice formado en la simulación
que usó el modelo de turbulencia k-ω SST. En la Figura 3.17.a) se evidencia que la cuerda
del vórtice se rompe en la zona central del codo del tubo de descarga, mientras que en la
Figura 3.17.b) el vórtice es mucho más pronunciado y no presenta un corte en su cuerda.
A pesar de ello las estructuras de vórtice formadas en ambos casos presentan una
tendencia fuerte en su zona central, y siguen el patrón helicoidal que se discute en los
estudios de Favrel et. al [40].

3.4. Simulación con Método PIMPLE: Estado inestable


La simulación en estado inestable fue configurada con el caudal de 62,4 [m3/s], debido a
que este valor de caudal presentaba la predicción más precisa en potencia con respecto a
valores experimentales, tal como se muestra en la Sección 3.3. Ya que la simulación con
método PIMPLE siguió la metodología de malla móvil (sliding mesh), los resultados son
presentados para cuatro estados: i) ¼ de revolución (0.045 segundos), ii) ½ de revolución
(0.091 segundos), iii) ¾ de revolución (0.137 segundos) y iv) 1 revolución (0.183
segundos). Los resultados de los campos de variables se presentan a continuación.

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)

En la Figura 3.20 y Figura 3.21., se muestra el perfil de presión en el dominio de álabes


directrices y rodete. La distribución de presión sobre los álabes directrices es mucho más
pronunciada que en el caso de los resultados de la simulación con método SIMPLE,
también se observa que el borde de ataque de los álabes del rodete presenta valores de

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)

En la Figura 3.24. se muestra la formación del vórtice o remolino inmediatamente después


de que el fluido abandona el dominio del rodete e ingresa al tubo de descarga. El vórtice
generado sigue la forma del perfil de energía cinética turbulenta en el cual la mayor
cantidad de flujo turbulento se encuentra en las paredes del tubo de descarga, y cuando
cumple una revolución (iv), se puede observar que el comportamiento es mucho más
caótico que en el caso del vórtice obtenido en la simulación con método SIMPLE. Por otra
parte, la Figura 3.24., no muestra los cuatro estados de la simulación en estado inestable,
ya que el perfil del vórtice formado en ½ revolución y ¾ de revolución se asemeja mucho
al perfil del vórtice cuando cumple 1 revolución. También, se puede observar que no hay
una cuerda de vórtice definida debido al perfil caótico presente, sin embargo, es notorio
que la Figura 3.24.b) presenta estructuras de vórtice mucho más abundantes que en el
caso de la Figura 3.24.a). La tendencia de mayor producción de estructuras de vórtices en

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.

Se confirmó la ventaja de utilizar mallado estructurado 3D para simulaciones de dominios


computacionales complejos como son las turbinas hidráulicas. Ya que, este tipo de mallado
permite el uso de la metodología de malla deslizante (sliding mesh) para simulaciones en
estado inestable, debido a la homogeneidad en la distribución de celdas en las interfaces
de dominios. Además, el proceso de interpolación de datos entre celdas en un mallado
estructurado presenta un comportamiento estable y de una adecuada precisión. También,
esta malla permite el uso de la metodología MRF para simulaciones en estado estable.

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 concluye que el modelo de turbulencia si influye en la predicción de estructuras de


vórtice. Se evidenció que el modelo de turbulencia k-ω SST predice estructuras de vórtice
que se asemejan a estudios previos, mientras que el modelo de turbulencia k-ε genera
estructuras de vórtice más débiles en la zona central.

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.

Los resultados obtenidos de las simulaciones en estado estable reproducen


adecuadamente el fenómeno en función de los datos experimentales disponibles. Sin
embargo, las simulaciones en estado inestable no llegan a cumplir con los perfiles de
velocidad, energía cinética turbulenta y vorticidad, debido a que estas simulaciones se
desarrollan en una revolución. Es por ello que, con base en información de estudios previos
llevados a cabo por diferentes autores, se concluye que la simulación en estado inestable
debería cumplir por lo menos tres revoluciones para reproducir el fenómeno en su totalidad.

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.

4.2. Trabajos Futuros


Tras haber definido la metodología para estudiar turbinas tipo Francis con mallado
estructurado rotativo en el entorno de OpenFOAM, tanto para simulaciones en estado
estable como inestable. Es necesario realizar estudios en OpenFOAM para el caso de una
turbina tipo Francis con las mismas características, pero con el uso de la herramienta de
mallado snappyHexMesh, la misma que genera un mallado híbrido hexadominante que
puede reducir el tiempo de simulación drásticamente sin reducir la calidad de resultados.

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

[1] «Balance Nacional de Energía Eléctrica – ARCONEL».


https://www.regulacionelectrica.gob.ec/balance-nacional/ (accedido abr. 22, 2020).
[2] Ecuatran, «Conoce cuáles son las 9 Centrales Hidroeléctricas que existen en
Ecuador», Ecuatran. https://www.ecuatran.com/blog/conoce-cuales-son-las-9-
centrales-hidroelectricas-que-existen-en-ecuador/ (accedido abr. 22, 2020).
[3] «CELEC EP - Hidráulicos». https://www.celec.gob.ec/generacion/hidraulicos.html
(accedido abr. 22, 2020).
[4] V. Hidalgo, X. Luo, A. Peña, E. Valencia, R. Soto, y A. Yu, «Benefits of hydropower
research in Ecuador using OpenFOAM based on CFD technology», Congr. Cienc.
Tecnol. ESPE, vol. 9, n.o 1, pp. 123-127, may 2014, doi: 10.24133/cctespe.v9i1.93.
[5] C. Mora, «Análisis de la eficiencia de una turbina tipo francis con características
similares a la de la central hidroeléctrica San Francisco-Ecuador», jul. 2018, Accedido:
abr. 22, 2020. [En línea]. Disponible en:
http://bibdigital.epn.edu.ec/handle/15000/19561.
[6] H. Zhang y L. Zhang, «Numerical simulation of cavitating turbulent flow in a high head
Francis turbine at part load operation with OpenFOAM», Procedia Eng., vol. 31, pp.
156-165, ene. 2012, doi: 10.1016/j.proeng.2012.01.1006.
[7] E. J. Guascal Sanguña y P. A. Quispe Quispe, «Desarrollo y estudio de un mallado
estructurado para un Turbina tipo Francis con validación experimental», sep. 2019,
Accedido: abr. 22, 2020. [En línea]. Disponible en:
http://bibdigital.epn.edu.ec/handle/15000/20484.
[8] G. Chen, Q. Xiong, P. J. Morris, E. G. Paterson, A. Sergeev, y Y. C. Wang, «OpenFOAM
for computational fluid dynamics», Not. Am. Math. Soc., vol. 61, n.o 4, pp. 354-363, abr.
2014, doi: 10.1090/noti1095.
[9] P. Breeze, Power Generation Technologies, Edición: 2. Amsterdam: Newnes, 2014.
[10] R. Goyal, «FLOW FIELD IN A HIGH HEAD FRANCIS TURBINE DRAFT TUBE
DURING TRANSIENT OPERATIONS», 2017, Accedido: abr. 22, 2020. [En línea].
Disponible en: http://urn.kb.se/resolve?urn=urn:nbn:se:ltu:diva-66297.
[11] H. Kapoor, «Open Source Mesh Generation and CFD Simulations for Francis Turbine»,
2014, Accedido: abr. 22, 2020. [En línea]. Disponible en:
https://odr.chalmers.se/handle/20.500.12380/207128.
[12] B. Aakti, O. Amstutz, E. Casartelli, G. Romanelli, y L. Mangani, «On the performance
of a high head Francis turbine at design and off-design conditions», vol. 579, p. 012010,
ene. 2015, doi: 10.1088/1742-6596/579/1/012010.

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.

Figura A.1. Submenú ICEM para escoger solver


(Fuente: propia)

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.

Figura A.2. Ventana emergente para guardar el archivo .msh.


(Fuente: propia)

69
ANEXO II.
ARCHIVO EJECUTABLE PARA LA UNIÓN DEL DOMINIO TOTAL,
MANIPULACIÓN DEL DOMINIO DENTRO DE OPENFOAM Y
FICHERO TOPOSETDICT

El archivo ejecutable usado para la unión de los subdominios de la turbina se le dio el


nombre de makemergeMesh-sh. A continuación, se muestra lo que contiene el archivo.
rm -r final
if [[ "$1" -eq "" ]] ; then
numProcs=8
else
numProcs=$1
fi
echo numProcs
echo $numProcs
echo "Merging meshes..."
mkdir final
cp -rv rodete/constant final/ > log.finalMesh-cp 2>&1
cp -rv rodete/system final/ >> log.finalMesh-cp 2>&1
mergeMeshes final tuboDescarga -overwrite > log.finalMesh-mergeMeshes 2>&1
mergeMeshes final cajaEspiral -overwrite > log.finalMesh-mergeMeshes 2>&1
mergeMeshes final alabesDirectrices -overwrite > log.finalMesh-mergeMeshes 2>&1
cd final
cd ..
echo "Copying mesh into case directory..."
cp -rv mesh/final/constant/polyMesh constant > log.finalMesh-cpMeshToCase 2>&1

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.

A continuación, se muestra el contenido del archivo topoSetDict, el mismo que permite


identificar las zonas que serán tomadas como estáticas y el dominio rotacional.

/*--------------------------------*- C++ -*----------------------------------*\


| ========= | |

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

Los ficheros de información de condiciones de borde son los que se encuentran en la


carpeta 0.orig, y son los que contienen información sobre la velocidad de entrada, presión
en la entrada y salida del dominio total, el valor inicial de energía cinética turbulenta, valor
inicial de ε y valor inicial de ω. A continuación, se muestra el contenido de estos ficheros.
/*--------------------------------*- C++ -*----------------------------------*\
71
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.5-dev |
| \\ / A nd | Revision: exported |
| \\/ M anipulation | Web: http://www.OpenFOAM.org |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
location "0";
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
//- Set patchGroups for constraint patches
#includeEtc "caseDicts/setConstraintTypes"
BLADE //wall
{
type movingWallVelocity;
value uniform (0 0 0);
}
BLADE1 //wall
{
type movingWallVelocity;
value uniform (0 0 0);
}
WALL_RUNNER_ //wall
{
type movingWallVelocity;
value uniform (0 0 0);
}

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
}
// ************************************************************************* //

/*--------------------------------*- C++ -*----------------------------------*\


| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.5-dev |
| \\ / A nd | Revision: exported |
| \\/ M anipulation | Web: http://www.OpenFOAM.org |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
location "0";
object nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 2 -1 0 0 0 0];
internalField uniform 1.106e-6;
boundaryField
{
//- Set patchGroups for constraint patches
#includeEtc "caseDicts/setConstraintTypes"
BLADE //wall
{
type nutUSpaldingWallFunction;
value uniform 0;
}
BLADE1 //wall
{
type nutUSpaldingWallFunction;
value uniform 0;
}
WALL_RUNNER_ //wall
{
type nutUSpaldingWallFunction;
value uniform 0;
}
WALL_DRAFT_TUBE //wall
{
type nutUSpaldingWallFunction;

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
}
// ************************************************************************* //

/*--------------------------------*- 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;
object k;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 2 -2 0 0 0 0];
internalField uniform 1.168949;

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
}
// ************************************************************************* //

/*--------------------------------*- C++ -*----------------------------------*\


| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 4.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/

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
}
// ************************************************************************* //

/*--------------------------------*- C++ -*----------------------------------*\


========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox

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
}
// ************************************************************************* //

El archivo que contiene la información de malla deslizante tiene el nombre de


dynamicMeshDict, y a continuación, se muestra su contenido.

/*--------------------------------*- C++ -*----------------------------------*\


========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant";
object dynamicMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dynamicFvMesh dynamicMotionSolverFvMesh;
motionSolverLibs ("libfvMotionSolvers.so");
motionSolver solidBody;
cellZone rotor;
solidBodyMotionFunction rotatingMotion;
origin (0.0001245737075805664 -4.5239925384521484e-5 0.29139767587184906);
axis (0 0 1);
omega 34.271634; // rad/s
// ************************************************************************* //

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.

/*--------------------------------*- C++ -*----------------------------------*\


========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org

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

Para iniciar la simulación se ejecuta el comando stk2, y a continuación se muestra su


contenido.
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
# Source tutorial run functions
. $WM_PROJECT_DIR/bin/tools/RunFunctions
# SEGUNDO PASO
cp -r 0.orig 0
runApplication decomposePar
runParallel patchSummary
runParallel simpleFoam
#runApplication reconstructPar
#------------------------------------------------------------------------------

Mientras que para el cálculo en paralelo se necesita del fichero decomposeParDict, el


mismo que contiene la información sobre cuantos núcleos se va a usar y la forma de
división del dominio para el cálculo en paralelo. A continuación, se muestra el contenido de
este fichero.

/*--------------------------------*- C++ -*----------------------------------*\


| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 5 |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object decomposeParDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
numberOfSubdomains 6;
//method hierarchical;
method simple;
//globalFaceZones (AMI_IMPELLER_INLET AMI_IMPELLER_VOLUTE);
simpleCoeffs
{
n (2 3 1);
delta 0.001;
}
// ************************************************************************* //

81

También podría gustarte