Modelado CFD
Modelado CFD
Autor
Carlos Moreno Montagud
Tutor
Marcos Carreres Talens
Cotutor
Mario Belmar Gil
En la sociedad actual existe cada vez mayor preocupación por las emisiones con-
taminantes, especialmente el 𝑁 𝑂𝑥 , y eso se traduce en la creciente concienciación de
la población sobre dicho problema así como en leyes cada vez más estrictas. En este
contexto surgen estrategias de inyección y combustión diferentes que tratan de corregir
este problema, como la Lean Direct Injection (LDI) para motores de flujo contínuo.
III
Resum
En la societat actual hi ha cada vegada una major preocupació per les emissions
contaminants, especialment el 𝑁 𝑂𝑥 , fet que es tradueix en la concienciació de la po-
blació sobre aquest problema així com en lleis cada vegada més estrictes. En aquest
context, sorgeixen estratègies d’injecció i combustió diferents que tracten de corregir el
problema, com la Lean Direct Injection (LDI) per a motors de flux continu.
V
Abstract
In this Final Master’s Project, we intend to make progress in the study of this
strategy through numerical resolution with CFD codes. The steps to follow are: analy-
zing the fuel injection, its atomization and subsequent evaporation, its interaction with
the turbulent flow in the combustion chamber, the global efficiency and the emissions
produced in this process. Based on experimental results, obtained from the LDI com-
bustion chamber designed in CORIA for the KIAI project, an attempt will be made
to characterize the turbulent velocity field that is produced internally. To do this, non-
reactive cases will be simulated and compared, premixed (injecting a poor mixture of
air and fuel) and not premixed (injecting gaseous fuel directly into the turbulent air
flow). The CFD codes used are CONVERGE (with adaptive mesh) and OpenFOAM
(with traditional mesh), using RANS and LES approaches.
VII
Great works are performed not by strength
but by perseverance.
Samuel Johnson
En primer lugar, al tutor Marcos Carreres y cotutor Mario Belmar por guiarme
durante este proyecto, así como en mi estancia en el CMT. Me han brindado las he-
rramientas necesarias para sacarlo adelante, me han ayudado con sus acertadas co-
rrecciones y gracias a ellos recordaré siempre con cariño la experiencia vivida y los
conocimientos adquiridos en ella.
En segundo lugar, a mis compañeros de carrera Jose Manuel Sellés e Iván Olmeda
que han hecho más llevaderos estos años de universidad, los trabajos en grupo, el
período de prácticas, los descansos a mediodía para comer y las horas de gimnasio.
Y finalmente, pero no menos importante, a mis padres Diego y Ana, así como a
mi hermano Diego, por su dedicación e insistencia para que pueda llegar lo más lejos
posible en la vida.
XI
Índice general
Resumen III
Resum V
Abstract VII
Agradecimientos XI
Documento I Memoria 1
1 Introducción 7
1.1. Contexto histórico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.2. Objetivos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
4 Metodología 47
4.1. Configuración experimental . . . . . . . . . . . . . . . . . . . . . . . . 47
4.2. Software empleado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
XIII
ÍNDICE GENERAL
5 Resultados 69
5.1. Caso premezclado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
5.2. Caso no premezclado . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
Bibliografía 95
1 Introducción 183
XIV
ÍNDICE GENERAL
1 Introducción 213
XV
Índice de figuras
1.1. Emisiones de 𝐶𝑂2 del sector de aviación en seis escenarios diferentes según
estimaciones del uso de combustible . . . . . . . . . . . . . . . . . . . . . . 8
1.2. Número de vuelos realizados por la industria global de aerolíneas desde 2004
hasta 2018 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.3. Distribución geográfica de flujos NOx en el año 2000 . . . . . . . . . . . . . 9
A.1. Contorno de energía cinética turbulenta del caso premezclado LES_17_OF. 105
A.2. Contornos de velocidades del caso premezclado LES_17_OF. . . . . . . . 106
A.3. Contorno de energía cinética turbulenta del caso premezclado RANS_17_CG.107
A.4. Contornos de velocidades del caso premezclado RANS_17_CG. . . . . . . 108
A.5. Contorno de energía cinética turbulenta del caso premezclado
URANS_18_OF. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109
A.6. Contornos de velocidades del caso premezclado URANS_18_OF. . . . . . 110
XIX
Documento I
MEMORIA
Índice de documento
1 Introducción
1.1. Contexto histórico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.1.1. Impacto Ambiental De Las Aeronaves . . . . . . . . . . . . . . . 7
1.1.2. Reducción De Emisiones Contaminantes . . . . . . . . . . . . . 10
1.2. Objetivos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
3
Índice de documento
3.5.1. Pre-procesado . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.5.2. Solver . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
3.5.3. Post-procesado . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
4 Metodología
4.1. Configuración experimental . . . . . . . . . . . . . . . . . . . . . . . . 47
4.1.1. Geometría Del Quemador KIAI . . . . . . . . . . . . . . . . . . 47
4.1.2. Condiciones Iniciales Y De Contorno . . . . . . . . . . . . . . . 49
4.2. Software empleado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
4.2.1. OpenFOAM . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
4.2.2. CONVERGE . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4.3. Cálculo de errores y desviación . . . . . . . . . . . . . . . . . . . . . . 55
4.4. Configuración numérica . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
4.4.1. Geometría En OpenFOAM . . . . . . . . . . . . . . . . . . . . . 57
4.4.2. Geometría En CONVERGE . . . . . . . . . . . . . . . . . . . . 58
4.4.3. Modelado De La Turbulencia . . . . . . . . . . . . . . . . . . . 60
4.4.4. Variables De Estudio . . . . . . . . . . . . . . . . . . . . . . . . 60
4.5. Rutinas de posprocesado . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.6. Estudio de independencia de malla . . . . . . . . . . . . . . . . . . . . 63
4.6.1. Malla OpenFOAM . . . . . . . . . . . . . . . . . . . . . . . . . 63
4.6.2. Malla CONVERGE . . . . . . . . . . . . . . . . . . . . . . . . . 65
5 Resultados
5.1. Caso premezclado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
5.1.1. Simulaciones En OpenFOAM . . . . . . . . . . . . . . . . . . . 69
5.1.2. Simulaciones En CONVERGE . . . . . . . . . . . . . . . . . . . 72
5.1.3. Perfiles De Velocidades . . . . . . . . . . . . . . . . . . . . . . . 76
5.1.4. Contornos De Velocidad Y Energía Cinética Turbulenta . . . . 79
5.2. Caso no premezclado . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
5.2.1. Simulaciones En CONVERGE . . . . . . . . . . . . . . . . . . . 82
5.2.2. Perfiles De Velocidades Y Especies . . . . . . . . . . . . . . . . 85
5.2.3. Contornos De Velocidad Y Energía Cinética Turbulenta . . . . 90
Bibliografía
4
Índice de documento
5
Capítulo 1
Introducción
7
Memoria
8
1. INTRODUCCIÓN
9
Memoria
1.2. OBJETIVOS
10
Capítulo 2
Fundamentos de combustión
en motores aeroderivados
11
Memoria
Prompt NOx . Es una forma de thermal NOx que se forma en lugares cercanos
al frente de llama cuando se oxidan productos intermedios de la combustión como el
HCN.
Mecanismo Prompt
𝑁2 + 𝐶𝐻 ↔ 𝐻𝐶𝑁 + 𝑁
𝑁 + 𝑂2 ↔ 𝑁 𝑂 + 𝑂
𝐻𝐶𝑁 + 𝑂𝐻 ↔ 𝐶𝑁 + 𝐻2 𝑂
𝐶𝑁 + 𝑂2 ↔ 𝑁 𝑂 + 𝐶𝑂
12
2. FUNDAMENTOS DE COMBUSTIÓN EN MOTORES AERODERIVADOS
13
Memoria
Finalmente, en la zona de mezcla pobre (0.5 < 𝜑 < 0.7) se oxidan las altas concen-
traciones de hidrocarburos y otros productos intermedios de la combustión, tratando
de no sobrepasar temperaturas de 1900K que conllevarían una elevadísima formación
de NOx . No obstante, se debe asegurar una temperatura lo suficientemente alta que
permita consumir los restos de CO, UHC y hollín que no pueden ser expulsados en el
escape [17].
14
2. FUNDAMENTOS DE COMBUSTIÓN EN MOTORES AERODERIVADOS
Cada una está diseñada para optimizar aspectos concretos del rendimiento de la cámara
de combustión.
15
Memoria
mantenerla estable. Esta capacidad para cambiar de un tipo de operación a otro in-
troduce complejidad en el diseño y aumenta los costes. Además existen varios peligros
potenciales, nada deseables en aviación general y menos todavía en aviación comercial:
un límite de estabilidad estrecho, alta generación de ruido y posibilidad de auto-ignición
y retroceso de llama en la zona premezclada. [23]
El swirler (sección 2.4.2) toma una importancia notable ya que el proceso de mezcla
y combustión depende casi en su totalidad del campo de velocidades generado. El ver-
dadero desafío de los quemadores LDI es conseguir que la mezcla sea lo más homogénea
posible, evitando picos de temperatura locales en la cámara que disparen las emisiones
de NOx .
swirler que crea un alto grado de turbulencia para atomizar los chorros de combustible
en gotas muy finas. La dispersión de éstas y las zonas de cizalladura resultantes de los
swirlers adyacentes contribuyen a la homogeneización de la mezcla, al contrario que en
los LDI de un único inyector limitados por la pared de la cámara. La configuración tam-
bién permite controlar los inyectores individualmente para redistribuir el combustible
y reducir inestabilidades. [36]
17
Memoria
Concretamente, en el caso de los motores de aviación el peso y tamaño son muy im-
portantes, mientras que para motores industriales lo es la capacidad para funcionar con
diferentes combustibles. Además, independientemente de su aplicación se consideran
primordiales las bajas emisiones contaminantes.
En este apartado se dará una visión general de los elementos básicos de la cámara
de combustión en turbinas de gas. Existen otros elementos a parte de los mencionados
pero que no son relevantes en referencia a los objetivos que se plantean en el presente
trabajo.
2.4.1. Difusor
La pérdida de presión en la cámara de combustión ocurre fundamentalmente cuando
se añade calor a un flujo de gas según la fórmula 2.1, donde 𝑇3 y 𝑇4 son respectivamente
las temperaturas a la entrada y salida de la cámara. Como ésta es proporcional al
cuadrado de la velocidad del aire, con velocidades a la salida del compresor de 170
m/s y un ratio de temperaturas de 2.5, las pérdidas originadas pueden llegar a suponer
el 25 % de la presión alcanzada. Para reducirlas a un nivel aceptable (hasta la quinta
parte) se utiliza el difusor [11].
𝑇4
(︂ )︂
Δ𝑃ℎ𝑜𝑡 = 0.5𝜌𝑈 2 −1 (2.1)
𝑇3
En su forma más simple, el difusor es un pasaje divergente por donde el flujo se
decelera, convirtiendo la velocidad en presión estática. La eficiencia de este proceso es
muy importante ya que las pérdidas ocurridas se manifiestan en una pérdida de presión
total, y por tanto, en el rendimiento del ciclo.
Por una parte, en difusores largos con ángulo de divergencia pequeño, la pérdida de
presión se debe principalmente a la fricción con las paredes. Por otro lado, los difusores
cortos con mayor ángulo de divergencia, tienen menor longitud y fricción, pero sufren
principalmente las pérdidas por separación de la capa límite. Lo ideal en el diseño de
estos componentes es un equilibrio entre ambas geometrías, ya que para una relación
de áreas determinada existe un ángulo de divegencia óptimo con pérdidas mínimas que
suele estar entre 6-12°.
18
2. FUNDAMENTOS DE COMBUSTIÓN EN MOTORES AERODERIVADOS
2.4.2. Swirler
Para posibilitar la combustión es necesario un flujo recirculatorio de aire y productos
quemados que cumple dos funciones: evitar el desprendimiento de llama (por el campo
de velocidades, figura 2.5) y mantener una fuente contínua de ignición para la mezcla
entrante [37]. Esto puede obtenerse mediante varias soluciones, pero concretamente en
este trabajo se optará por el swirler, porque es el que dispone la geometría a estudiar
y ofrece ventajas como regiones de fuertes tensiones cortantes, alta turbulencia y tasas
de mezcla rápidas.
Las características del flujo de aire producido por el swirler y sus aplicaciones prác-
ticas han sido estudiadas en [38, 39], siendo la última de especial interés porque se dirige
directamente a los tipos de swirlers más relevantes en turbinas de gas, proporcionando
información valiosa para el diseño de estos componentes.
19
Memoria
𝐺𝜑
𝑆= (2.2)
𝐺𝑥 𝑟𝑡
∫︁ 𝑟𝑡
𝐺𝜑 = (𝑊 𝑟)𝜌𝑈 (2𝜋𝑟)𝑑𝑟
𝑟ℎ
∫︁ 𝑟𝑡
= (𝑈 𝑡𝑎𝑛𝜃)2𝜋𝜌𝑈 𝑟2 𝑑𝑟 (2.3)
𝑟ℎ
2
= 𝜋𝜌𝑈 2 (𝑟𝑡3 − 𝑟ℎ3 )
3
∫︁ 𝑟𝑡 ∫︁ 𝑟𝑡
𝐺𝑥 = 𝜌𝑈 2 (2𝜋𝑟)𝑑𝑟 + 𝑝(2𝜋𝑟)𝑑𝑟 (2.4)
𝑟ℎ 𝑟ℎ
𝐺𝜑
𝑆= (2.5)
𝐺′𝑥 𝑟𝑡
∫︁ 𝑟𝑡
𝐺′𝑥 = 𝜌𝑈 2 (2𝜋𝑟)𝑑𝑟
𝑟ℎ
∫︁ 𝑟𝑡
= 2𝜋𝜌𝑈 2 𝑟𝑑𝑟 (2.6)
𝑟ℎ
2 1 − (𝑟ℎ /𝑟𝑡 )3
[︃ ]︃
𝑆 = 𝑡𝑎𝑛𝜃 (2.7)
3 1 − (𝑟ℎ /𝑟𝑡 )2
donde 𝜃 corresponde al ángulo que forman las palas del swirler respecto a la dirección
axial. En los swirlers axiales de álabes planos este ángulo es constante, mientras que
en aquellos con álabes torsionados el ángulo es nulo en la sección de entrada y 𝜃 en la
20
2. FUNDAMENTOS DE COMBUSTIÓN EN MOTORES AERODERIVADOS
Basándose en el número de swirl, se puede dividir el flujo en dos regiones: débil (S <
0.6) y fuerte (S > 0.6). En el flujo rotatorio débil los gradientes de presión axiales son
insuficientes para causar recirculación. Opuestamente, el flujo rotatorio fuerte tiene
gradientes altos de presión en dirección radial y axial cerca de la salida del swirler,
consiguiendo originar dicha zona central de recirculación.
Tipos De Swirler. Existen varios tipos de swirler, pero los dos más comunes
son axial y radial, cuyos esquemas se muestran en la figura 2.6 y se suelen encontrar
dispuestos de manera individual o concéntrica.
Existen estudios relevantes sobre swirlers axiales como [40], donde emplean 6 palas
con un ángulo 60 grados en combinación con tobera convergente-divergente.
En el caso de los swirlers radiales hay estudios como [24, 25, 26, 27] con aplicación
directa en quemadores LDI. Sus investigaciones incluyen el efecto de uno y varios
swirlers contrarrotantes, y de la inyección del combustible en las emisiones de NOx ,
aunque utilizan gas natural en lugar de combustible líquido.
21
Memoria
2.4.3. Inyector
Aunque en la sección 2.4.2 se ha hablado de ambos tipos de combustible (gas y
líquido), la gran mayoría de los motores turbina de gas en la industria aeronáutica
operan con combustibles líquidos debido a su elevada densidad energética; el más utili-
zado es el keroseno por su alto poder calorífico y bajo punto de fusión. El combustible
se inyecta en forma de spray, proceso que puede tener lugar en la corriente del flujo de
recirculación (inyección directa) o en una zona anterior a la entrada de la cámara de
combustión (premezclado).
22
2. FUNDAMENTOS DE COMBUSTIÓN EN MOTORES AERODERIVADOS
Según los procesos mecánicos mediante los cuales se realizan la atomización primaria
y secundaria, existen diversas tecnologías de inyectores:
Figura 2.7: Variante del Pressure Swirl Atomiser más sencilla, co-
nocida como SIMPLEX. [44]
23
Memoria
cónica sobre una pared cilíndrica donde se encuentra con un primer chorro de aire
torbellinado que lo rompe en ligamentos y lo arrastra hasta el orificio de salida.
En la sección de salida del inyector un segundo chorro de aire torbellinado en
sentido contrario al primero genera una fuerte zona de cortadura y completa la
atomización secundaria. Por tanto, el combustible inyectado en la cámara de com-
bustión se encuentra totalmente atomizado y parcialmente mezclado, mejorando
así las prestaciones respecto a los inyectores pressure swirl.
24
Capítulo 3
Comenzó a utilizarse en los años 60, enfocado sobre todo para el diseño y fabricación
en el sector de la industria aeroespacial, pero rápidamente se extendió a otras áreas de
investigación durante los años 80-90 llegando a establecerse como parte esencial en el
desarrollo tecnológico e industrial.
Se espera que en un futuro éstas sean capaces de reemplazar a los métodos expe-
rimentales debido a que normalmente son muy costosos, a menudo son irrepetibles y
generalmente son muy dependientes de la instrumentación utilizada. Opuestamente,
una simulación puede repetirse innumerables veces, se puede tomar medidas en cual-
quier punto del dominio y se puede generar múltiples visualizaciones muy útiles para
ingenieros y diseñadores.
25
Memoria
décadas, las simulaciones de CFD necesitan manejar grandes cantidades de datos (la
geometría del problema se divide en millones de elementos, con sus respectivas variables
físicas y las ecuaciones complejas que se deben resolver en cada uno de ellos) haciendo
que lleguen a durar semanas, meses e incluso años.
La presión tiene signo negativo por convenio, al tratarse de una fuerza de compre-
sión. En estas ecuaciones se tienen en cuenta las tensiones superficiales y, además, la
presencia de términos fuente (𝑆𝑀 𝑥 , 𝑆𝑀 𝑦 y 𝑆𝑀 𝑧 ), que representan la acción debida a las
fuerzas de volumen como podría considerarse la gravedad.
𝐷𝑢 𝜕𝑝
𝜌 =− + 𝑑𝑖𝑣(𝜇 𝑔𝑟𝑎𝑑 𝑢) + 𝑆𝑀 𝑥 (3.5)
𝐷𝑡 𝜕𝑥
𝐷𝑣 𝜕𝑝
𝜌 =− + 𝑑𝑖𝑣(𝜇 𝑔𝑟𝑎𝑑 𝑣) + 𝑆𝑀 𝑦 (3.6)
𝐷𝑡 𝜕𝑦
𝐷𝑤 𝜕𝑝
𝜌 = − + 𝑑𝑖𝑣(𝜇 𝑔𝑟𝑎𝑑 𝑤) + 𝑆𝑀 𝑧 (3.7)
𝐷𝑡 𝜕𝑧
3.2.3. Energía
Según la Primera Ley de la Termodinámica, la tasa de variación de energía en una
partícula fluida equivale a la suma de tasa de adición de calor y el trabajo realizado por
la misma. Por un lado, la tasa de trabajo aportado a la partícula fluida ejercido por una
fuerza de superficie es igual al producto de fuerza y velocidad. Además hay que tener
en cuenta el fenómeno de conducción de calor hasta la partícula fluida considerada, así
como la tasa de adición de calor que este mismo fenómeno produce. La ecuación de la
energía queda como se muestra en 3.8, donde el término 𝑆𝐸 representa una fuente de
energía por unidad de volumen y tiempo.
27
Memoria
𝑝 = 𝜌𝑅𝑇 (3.9)
𝑖 = 𝐶𝑣 𝑇 (3.10)
28
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
La velocidad y longitud característica de los remolinos grandes son del mismo or-
den de magnitud que la escala de velocidad y longitud del flujo medio, estando su
comportamiento dominado, por tanto, por los efectos inerciales. La estructura de estos
remolinos es altamente anisótropa, presentando fluctuaciones diferentes en función de
la dirección y afectadas por las condiciones de contorno del problema.
No obstante, los remolinos pequeños presentan unas escalas de longitud del orden
de 0.1-0.01 mm y unas frecuencias en torno a los 10 kHz, predominando los efectos
viscosos. En estas escalas la energía asociada a los remolinos más pequeños se disipa en
forma de energía térmica interna, así que la estructura de los remolinos más pequeños
solamente depende de la tasa de disipación de la energía cinética turbulenta 𝜀 y de la
viscosidad cinemática del fluido 𝜈.
29
Memoria
3.3.1. RANS
Los modelos de turbulencia RANS (Reynolds-Averaged Navier-Stokes) están cen-
trados únicamente en capturar el flujo medio y los efectos que la turbulencia origina en
las propiedades medias del flujo. Las ecuaciones de Navier-Stokes se promedian tempo-
ralmente antes de aplicar los métodos numéricos, descartando los detalles del flujo en
fluctuaciones instantáneas y dando lugar a la aparición de términos extra como conse-
cuencia de las interacciones entre fluctuaciones turbulentas, los esfuerzos de Reynolds:
−𝜌𝑢′2 , −𝜌𝑣 ′2 −𝜌𝑤′2 , −𝜌𝑢′ 𝑣 ′ , −𝜌𝑢′ 𝑤′ y −𝜌𝑣 ′ 𝑤′ .
De esta manera, para ser capaz de calcular flujos turbulentos con las ecuaciones
RANS es necesario desarrollar modelos de turbulencia capaces de predecir los esfuerzos
de Reynolds y los términos de transporte escalar para cerrar el sistema. Los modelos de
turbulencia más comunes se clasifican en base al número de ecuaciones de transporte
adicionales que se necesitan resolver junto a las ecuaciones RANS y se enumeran en la
tabla 3.1.
30
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
Se puede obtener una ecuación para la energía cinética media 𝐾 multiplicando cada
una de las componentes de la ecuación de Reynolds por las correspondientes del vector
de velocidad media [48], obteniendo la ecuación 3.12.
⎛ ⎞
𝜕(𝜌𝐾)
+ 𝑑𝑖𝑣(𝜌𝐾 𝑈
⃗ ) = 𝑑𝑖𝑣 ⎜
⎝−𝑃 ⃗ + 2𝜇𝑈
𝑈 ⃗ 𝑆𝑖𝑗 − 𝜌𝑈 ⎠ − 2𝜇𝑆𝑖𝑗 · 𝑆𝑖𝑗 + 𝜌𝑢′𝑖 𝑢′𝑗 · 𝑆𝑖𝑗
⃗ 𝑢′𝑖 𝑢′𝑗 ⎟
𝜕𝑡 ⏟ ⏞ ⏟ ⏞
𝑐
⏟ ⏞ ⏟ ⏞ ⏟ ⏞ ⏟ ⏞
𝑏 𝑒 𝑔
⏟ ⏞
𝑎 𝑑 𝑓
(3.12)
Resulta sencillo demostrar que en flujos con alto Reynolds los términos turbulentos
(e) y (g) son siempre varios ordenes de magnitud mayores que sus correspondientes
términos viscosos (d) y (f).
31
Memoria
𝜕(𝜌𝜀) 𝜇𝑡 𝜀 𝜀2
(︂ )︂
+ 𝑑𝑖𝑣(𝜌𝜀𝑈 ) = 𝑑𝑖𝑣
⃗ 𝑔𝑟𝑎𝑑(𝜀) + 𝐶1𝜀 2𝜇𝑡 𝑆𝑖𝑗 · 𝑆𝑖𝑗 − 𝐶2𝜀 𝜌 (3.16)
𝜕𝑡 𝜎𝜀 𝑘 𝑘
Las ecuaciones contienen cinco constantes 𝐶𝜇 , 𝜎𝑘 , 𝜎𝜀 , 𝐶1𝜀 y 𝐶2𝜀 que deben ser ajus-
tadas mediante experimentos para un amplio rango de flujos turbulentos. Los términos
correspondientes al transporte turbulento se representan mediante la difusión del gra-
diente de la magnitud transportada, donde los números de Prandtl 𝜎𝑘 y 𝜎𝜀 conectan las
difusividades de 𝑘 y 𝜀 con la viscosidad turbulenta 𝜇𝑡 . Además, los términos de produc-
ción y destrucción de energía cinética turbulenta se encuentran fuertemente acoplados,
presentando un valor grande de tasa de disipación 𝜀 allí donde la producción de 𝑘 es
elevada.
2 𝑘 3/2
𝑘 = (𝑈𝑟𝑒𝑓 𝐼)2 𝜀 = 𝐶𝜇3/4 𝑙 = 0.07𝐿 (3.17)
3 𝑙
La principal ventaja de este modelo es que no hay necesidad de ajustar las constan-
tes, y su éxito está demostrado y validado para multitud de flujos turbulentos. Ocurre
incluso en aquellos con zonas de recirculación [52] como el del presente trabajo, donde
la predicción del transporte turbulento es imprescindible para reproducir correctamente
el proceso de combustión.
32
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
componente de swirl, o frente a aquellos que presenten una deformación rápida como
sucede en capas límites con un elevado grado de curvatura o en secciones divergentes.
𝐷𝑅𝑖𝑗 𝜕𝑅𝑖𝑗
= + 𝐶𝑖𝑗 = 𝑃𝑖𝑗 + 𝐷𝑖𝑗 − 𝜀𝑖𝑗 + Π𝑖𝑗 + Ω𝑖𝑗 (3.18)
𝐷𝑡 𝜕𝑡
Los términos correspondientes a la convección, producción y rotación se definen
como:
𝜕(𝜌𝑈𝑘 𝑢′𝑖 𝑢′𝑗 )
𝐶𝑖𝑗 = 𝜕𝑥𝑘
= 𝑑𝑖𝑣(𝜌𝑢′𝑖 𝑢′𝑗 𝑈
⃗) (3.19)
𝑃𝑖𝑗 = −(𝑅𝑖𝑚 𝑥𝑈𝑚𝑗 + 𝑅𝑗𝑚 𝑥𝑈𝑚𝑖 ) (3.20)
Ω𝑖𝑗 = −2𝜔𝑘 (𝑢′𝑗 𝑢′𝑚 𝑒𝑖𝑘𝑚 + 𝑢′𝑖 𝑢′𝑚 𝑒𝑗𝑘𝑚 ) (3.21)
(︃ )︃
𝜕 𝜐𝑡 𝜕𝑅𝑖𝑗 𝜐𝑡
(︂ )︂
𝐷𝑖𝑗 = = 𝑑𝑖𝑣 𝑔𝑟𝑎𝑑(𝑅𝑖𝑗 ) (3.22)
𝜕𝑥𝑚 𝜎𝑘 𝜕𝑥𝑚 𝜎𝑘
2
𝜀𝑖𝑗 = 𝜀𝛿𝑖𝑗 (3.23)
3
Finalmente, el término de interacciones presión-deformación es el más difícil de
modelar con precisión, ya que el efecto que causan se debe a dos fenómenos físicos in-
dependientes: un proceso lento que reduce la anisotropía de los remolinos turbulentos
por sus interacciones mutuas, y un proceso rápido como consecuencia de las interac-
ciones entre las fluctuaciones turbulentas y la deformación del flujo medio. La manera
más simple de representar el proceso lento considera que la transición hacia las condi-
ciones isótropas es proporcional al grado de anisotropía 𝑎𝑖𝑗 de los esfuerzos de Reynolds
dividido por un tiempo característico 𝑘/𝜀. Por el contrario, el efecto rápido únicamente
tiene en cuenta el proceso de producción causante de la generación de anisotropía.
𝜀 2 2
(︂ )︂ (︂ )︂
Π𝑖𝑗 = −𝐶1 𝑅𝑖𝑗 − 𝑘𝛿𝑖𝑗 − 𝐶2 𝑃𝑖𝑗 − 𝑃 𝛿𝑖𝑗 (3.24)
𝑘 3 3
Las seis ecuaciones para el transporte de los esfuerzos de Reynolds se resuelven
de manera conjunta, modelando además la tasa de disipación de la energía cinética
turbulenta de un modo prácticamente idéntico al modelo 𝑘 − 𝜀 (ecuación 3.25).
𝐷𝜀 𝜐𝑡 𝜀 𝜀2
(︂ )︂
= 𝑑𝑖𝑣 𝑔𝑟𝑎𝑑(𝜀) + 𝐶1𝜀 2𝜐𝑡 𝑆𝑖𝑗 · 𝑆𝑖𝑗 − 𝐶2𝜀 (3.25)
𝐷𝑡 𝜎𝜀 𝑘 𝑘
3.3.2. URANS
La aproximación Unsteady Reynolds-Averaged Navier Stokes (URANS) se carac-
teriza por ser una versión más compleja de RANS. En ella se obtiene una solución
temporal para problemas tridimensionales no estacionarios, con o sin tratamiento es-
pecial para inestabilidades del flujo [55].Las ecuaciones de URANS son similares a los
modelos RANS pero incluyendo un término no estacionario adicional. El aspecto más
importante es la selección correcta del modelo de turbulencia adecuado y el esquema
de discretización (que no sea excesivamente disipativo). URANS actúa bien para calles
de vórtices de von Kármán [56], cubos montados sobre superficies [57] y flujo alrededor
de coches [58] entre otras aplicaciones.
3.3.3. LES
El Large Eddy Simulation (LES) se trata de un paso intermedio a la hora de obtener
y calcular las interacciones de las fluctuaciones turbulentas. En este caso únicamente
se modelan los remolinos más pequeños, que presentan un comportamiento isótropo e
independiente del problema, mientras que el comportamiento de los remolinos grandes
sí que es calculado directamente. Éstos son los encargados de interactuar y extraer
energía del flujo medio, se caracterizan por una elevada anisotropía y una fuerte de-
pendencia tanto de la geometría del dominio del problema como de las condiciones de
contorno impuestas.
Esta separación entre los remolinos por sus tamaños se consigue definiendo una lon-
gitud de corte (cutoff width), que permite realizar un filtrado espacial de las ecuaciones
de Navier-Stokes no estacionarias. Para modelar los efectos de los remolinos pequeños
sobre el flujo resuelto se utilizan los modelos SGS (sub-grid scale).
Esta operación de filtrado espacial (ecuación 3.26) se lleva a cabo mediante una
función de filtro 𝐺(⃗𝑥, 𝑥⃗′ , Δ), con longitud de corte Δ.
∫︁∫︁∫︁ ∞
𝜑(⃗𝑥, 𝑡) = 𝐺(⃗𝑥, 𝑥⃗′ , Δ) 𝜑(𝑥⃗′ , 𝑡)𝑑𝑥′1 𝑑𝑥′2 𝑑𝑥′3 (3.26)
−∞
35
Memoria
⎨1/Δ3
⎧
𝑠𝑖 |⃗𝑥 − 𝑥⃗′ | ≤ Δ/2
𝐺(⃗𝑥, 𝑥⃗′ , Δ) = (3.27)
⎩ 0 𝑠𝑖 |⃗𝑥 − 𝑥⃗′ | ≥ Δ/2
donde 𝜏𝑖𝑗 = 𝜌𝑢𝑖⃗𝑢 − 𝜌𝑢𝑖⃗𝑢 = 𝜌𝑢𝑖 𝑢𝑗 − 𝜌𝑢𝑖 𝑢𝑗 . Se les llama esfuerzos SGS porque son
originados por el transporte convectivo e interacciones entre los remolinos modelados, es
decir, los que tienen una escala de longitud inferior al tamaño de malla. Además se ven
afectados por diversas contribuciones, cuya naturaleza se determina descomponiendo
las variables del flujo 𝜑(⃗𝑥, 𝑡) como suma de la función filtrada 𝜑(⃗𝑥, 𝑡) y las variaciones
espaciales pequeñas 𝜑′ (⃗𝑥, 𝑡).
36
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
1 1
(︃ )︃
𝜕𝑢𝑖 𝜕𝑢𝑗
𝑅𝑖𝑗 = −2𝜇𝑆𝐺𝑆 𝑆𝑖𝑗 + 𝑅𝑖𝑖 𝛿𝑖𝑗 = −𝜇𝑆𝐺𝑆 + + 𝑅𝑖𝑖 𝛿𝑖𝑗 (3.33)
3 𝜕𝑥𝑗 𝜕𝑥𝑖 3
Para saber cuánto se debe refinar una malla en la zona cercana a una pared se
utiliza el parámetro adimensional 𝑦 + , que se calcula mediante la fórmula 3.36 donde
𝑢𝜏 es la velocidad cortante o de fricción, 𝑦 la distancia a la pared y 𝜈 la viscosidad
cinemática. La ley de pared de von Karman, dicta que dependiendo del valor de este
37
Memoria
parámetro conviene utilizar una aproximación concreta para obtener la mayor precisión
posible en el cálculo de la velocidad (figura 3.3): la capa viscosa, la capa buffer y la
capa logarítmica.
𝑦𝑢𝜏
𝑦+ = (3.36)
𝜈
38
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
Por su parte, las condiciones de contorno a la salida del dominio son menos pro-
blemáticas, empleando la conocida condición de gradiente nulo para el flujo medio,
mientras que las propiedades fluctuantes son extrapoladas mediante:
𝜕𝜑 𝜕𝜑
+ 𝑢𝑛 =0 (3.37)
𝜕𝑡 𝜕𝑛
La demanda de recursos computacionales en términos de almacenamiento y volu-
men de operaciones es enorme, ya que se deben resolver las ecuaciones del flujo no
estacionarias. A pesar de ello, los avances en tecnología de computación están permi-
tiendo progresivamente la aplicación de estas técnicas LES para el estudio de problemas
con geometrías complejas. Se trata de una gran inversión de tiempo y recursos, pero
tiene en cuenta los detalles estadísticos relacionados con las fluctuaciones turbulentas
presentes en flujos con altas componentes de swirl.
3.3.4. DNS
En el Direct Numerical Simulation (DNS) se calcula el flujo medio y todas las
fluctuaciones turbulentas de la velocidad. Para ello, se resuelven las ecuaciones de
Navier-Stokes no estacionarias en mallas lo suficientemente finas (del orden de la escala
de Kolmogorov, donde tiene lugar la disipación de la energía) y con pasos temporales
extremadamente pequeños del orden del periodo de las fluctuaciones turbulentas más
rápidas. Estos cálculos son altamente costosos en términos de recursos computacionales,
motivo por el cual estos métodos todavía no se emplean en problemas industriales.
39
Memoria
donde las funciones 𝑁𝑖 (𝑥) son la base del desarrollo en serie de 𝑢, y se usa 𝑤(𝑥) =
𝑁𝑖 (𝑥), se obtiene el método de los elementos finitos. En él se usan funciones base
contínuas. En general, suelen ser polinomios (o polinomios definidos a trozos).
Cada 𝑁𝑖 es distinto de 0 sólo para cierta zona del dominio.
Métodos espectrales. Si el dominio en el que 𝑁 𝑖(𝑥) es distinto de 0 (soporte de
la base) es todo el intervalo, se tiene un método espectral. Éstos tienen conver-
gencia espectral: convergen de forma exponencial si la solución es suave, pero
lamentablemente, no capturan bien las discontinuidades.
Método de volúmenes finitos. Si se discretiza el dominio en 𝑛 volúmenes de control
y usamos 𝑛 funciones para la base tales que 𝑁𝑖 vale 1 sólo en el volumen de control
𝑖 y 0 fuera de él, se obtiene el método de volúmenes finitos. Este método es
extremadamente interesante para resolver las ecuaciones de Navier-Stokes (o una
aproximación de las mismas) y puede permitir capturar bien discontinuidades
40
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
En este apartado se tratarán por orden las tres etapas existentes y bien diferenciadas
que tienen lugar en la simulación CFD.
3.5.1. Pre-procesado
Para que las ecuaciones puedan ser aplicadas en el dominio computacional, debe
generarse un mallado, cuyas celdas o elementos representarán los volúmenes de control
que componen la geometría del problema. Es necesario poner especial atención en este
punto, pues de la calidad de la malla dependerá en gran medida la precisión de la
simulación. La solución del problema se define en los nodos dentro de cada celda,
con lo cual la cantidad de éstas influye en la validez de resultados, pero también en
el tiempo necesario para la resolución. A mayor número de celdas, tanto la precisión
como el coste computacional aumentan.
41
Memoria
En la teoría, una malla enorme con elementos infinitamente pequeños sería ideal
para unos resultados totalmente precisos, pero sería inabordable por su coste compu-
tacional. Por ello, independientemente del tipo de malla escogido, se debe asegurar que
la cantidad de celdas (así como el tamaño de las mismas) representa de forma sufi-
cientemente fidedigna la geometría del problema, a la vez que no supone un tiempo de
cálculo exagerado. En este punto se pone de manifiesto la importancia del análisis de
independencia de malla, cuyo objetivo es garantizar ambas cosas.
Este proceso es simple teóricamente: se estudia el mismo caso varias veces, modi-
ficando únicamente el refinamiento de la malla hasta que los resultados no sufren una
42
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
3.5.2. Solver
En el solver se resuelve la matriz o sistema de ecuaciones algebraicas que se ha
obtenido tras la discretización de las ecuaciones integrales y de la geometría y físicas
del problema. Esto se hace a través de un proceso iterativo, hasta la convergencia,
o sea, cuando la variación de los resultados entre dos iteraciones consecutivas queda
dentro de los límites en un criterio establecido por el usuario.
En el mercado hay infinidad de software, tanto de libre uso como de licencia comer-
cial, dedicado al cálculo CFD. En la sección 4.2 se exponen sendos programas utilizados
durante este trabajo para las simulaciones.
𝑢Δ𝑡
𝐶= (3.42)
Δ𝑥
Partiendo de la matriz con las ecuaciones de momento, se obtiene la presión al re-
solver la ecuación de momento (solve momentum). Entonces se aplica una corrección
en la presión (pressure corrector). Esta presión obtenida sirve para corregir las ecua-
ciones de momento y flujo (momentum and flux correctors), y se añaden a la matriz
los nuevos resultados. Este bucle (SIMPLE loop) puede ser repetido las veces que sea
necesario hasta alcanzar la convergencia. Una vez terminado se podrán resolver el resto
de ecuaciones de transporte en serie.
43
Memoria
44
3. MODELADO MEDIANTE MECÁNICA DE FLUIDOS COMPUTACIONAL
3.5.3. Post-procesado
Tras la simulación llega el último paso, donde se deben analizar los resultados. Para
ello es imprescindible conocer las magnitudes representativas del problema en cuestión,
y realizar gráficas que permitan obtener conclusiones sobre el estudio que se ha llevado
a cabo. Existe una gran variedad de posibilidades para el tratamiento de los resultados,
tanto desde el punto de vista gráfico (posibilidad de obtener resultados en secciones o
líneas, realización de animaciones, etc.) como a la hora de gestionar los resultados (estos
pueden ser exportados en múltiples formatos para su posterior tratamiento, se pueden
definir nuevas variables físicas a partir de las que vienen por defecto, guardado de las
trazas de las partículas de una fase líquida dentro de otra continua, etc.). En definitiva,
el ingeniero dispone de una amplia gama de utilidades dentro del post-procesado para
exprimir al máximo el valor de su simulación.
45
Memoria
Otra opción pasa por crear rutinas propias en algún lenguaje de programación para
tener mayor control sobre las magnitudes deseadas y la organización o almacenaje de las
mismas. Esto también permite obtener en gráficas con mayor grado de personalización
y automatización los datos generados por el solver.
46
Capítulo 4
Metodología
Para validar el código CFD se emplea la geometría del quemador KIAI mono-
inyector instalado en CORIA (COmplexe de Recherche Interprofessionel en Aérother-
mochimie) para estudios dedicados a la comprensión del fenómeno de ignición [68, 69].
Es una unidad mixta de investigación del CNRS, la Universidad de Rouen y el INSA
de Rouen, y forma parte del proyecto de investigación europeo KIAI (Knowledge for
Ignition, Acoustics and Instabilities) que aborda soluciones innovadoras para el desa-
rrollo de lean combustors fiables en motores de aviación. Los estudios presentados han
sido llevados a cabo con configuraciones académicas simplificadas, pero siendo sufi-
cientemente representativas para proporcionar resultados relevantes para aplicaciones
industriales.
47
Memoria
(a) Esquema descriptivo (b) Fotografía (c) Zoom del esquema descriptivo
48
4. METODOLOGÍA
49
Memoria
En este apartado se describen los principales programas de cálculo CFD que han sido
empleados en la realización del presente trabajo. También se ha utilizado una máquina
virtual debido a que cada uno necesitaba un sistema operativo distinto (Linux Ubuntu
16.04 y Microsoft Windows 10).
4.2.1. OpenFOAM
El software de CFD gratuito y de código abierto OpenFOAM se puede resumir como
un conjunto de módulos y bibliotecas en lenguaje C++ que permiten la resolución de
50
4. METODOLOGÍA
No obstante, esta comunidad desarrolladora sufre una gran fragmentación con nu-
merosos proyectos dando lugar a bifurcaciones del software que avanzan en paralelo.
La ausencia de interfaz gráfica y la pronunciada curva de aprendizaje necesaria para
agregar funcionalidades (User Applications) son los mayores inconvenientes y limitan-
tes del uso del software por parte del usuario medio. Además, algunas herramientas
innovadoras como el mallado adaptativo (AMR) todavía se encuentran en una fase
temprana con limitaciones en su uso.
En la figura 4.6 se muestra la estructura del directorio de un caso básico con el con-
junto de archivos mínimos para lanzar una simulación. Se compone de un directorio
/0 que contiene las condiciones iniciales y de contorno de la presión y velocidad, así
como de aquellas que requieran los modelos empleados. En el directorio /constant se
almacena la toda la información concerniente a la malla (polyMesh) y a las propiedades
físicas de los fenómenos acontecidos en la aplicación de estudio (modelos de turbulen-
cia, propiedades del fluido, etc.). Por último, en el directorio /system se encuentran
los parámetros de configuración relacionados con los esquemas de discretización em-
pleados para cada función de las ecuaciones (fvSchemes), los algoritmos numéricos con
sus respectivas tolerancias (fvSolution) y los parámetros temporales que controlan la
simulación (controlDict).
51
Memoria
4.2.2. CONVERGE
CONVERGE es un software comercial de CFD basado en el método de volúmenes
finitos que elimina el cuello de botella que supone la generación de la malla en las etapas
de pre-proceso de las simulaciones. Fue desarrollado inicialmente para aplicaciones con-
cernientes a los motores de combustión interna alternativos, pero durante los últimos
años su campo de aplicación se está ampliando e impulsando a otros ámbitos ingenie-
riles, especialmente al estudio de turbinas de gas. A diferencia del resto de programas
CFD, CONVERGE genera de manera automática una malla ortogonal, perfectamen-
te estructurada durante la propia simulación mediante unos sencillos parámetros de
control definidos por el usuario [72].
Para crear una geometría desde el propio software o importar un modelo 3D exis-
tente se puede utilizar CONVERGE Studio, herramienta complementaria y adaptada
a este solver (ya que ha sido desarrollado por la misma compañía) que permite el
pre-proceso y post-proceso del caso. A la hora de reconocer la geometría, el programa
necesita representar a través de triángulos las superficies, acción que lleva a cabo auto-
máticamente tras la importación del archivo. Sin embargo, el usuario tiene capacidad
de modificar dichos triángulos para reparar errores que se puedan encontrar, o incluso
52
4. METODOLOGÍA
Una vez terminado este paso, ya no es necesario repetirlo más, siempre y cuando
se vaya a usar la misma geometría. Será CONVERGE el que, interpretando la dis-
posición de los elementos, cree la malla hexaédrica estructurada (que presenta mayor
eficiencia computacional y simplifica los métodos numéricos del solver respecto a las no
estructuradas, como se ha mencionado en el apartado 3.5) siguiendo las indicaciones
del usuario sobre las características de la misma. Por estos motivos presenta enormes
ventajas en términos de eficiencia.
Grid Scaling, que refina o ensancha de manera global el dominio completo. Re-
sulta especialmente útil durante las primeras iteraciones de la simulación hasta
que el flujo se establece por completo.
53
Memoria
54
4. METODOLOGÍA
Cuando se trata de comparar gráficas, la manera más intuitiva consiste en hacer una
primera inspección visual. Sin embargo, este método es insuficiente para extraer toda
la información necesaria, sirve únicamente para obtener una primera aproximación.
La forma habitual de calcular errores (en este caso se entiende error como la dife-
rencia entre el valor de la simulación y el experimental) es mediante la fórmula del error
relativo (ecuación 4.1). Con este método, ampliamente utilizado en todos los campos
de la ciencia, se consigue un valor porcentual de la desviación.
|𝑚𝑎𝑔𝑛𝑖𝑡𝑢𝑑𝑟𝑒𝑎𝑙 − 𝑚𝑎𝑔𝑛𝑖𝑡𝑢𝑑𝑚𝑒𝑑𝑖𝑑𝑎 |
𝑒𝑟𝑟𝑜𝑟𝑟𝑒𝑙𝑎𝑡𝑖𝑣𝑜 = * 100 (4.1)
𝑚𝑎𝑔𝑛𝑖𝑡𝑢𝑑𝑟𝑒𝑎𝑙
Sin embargo, este método tiene algunas limitaciones, y es que para magnitudes
reales verdaderamente pequeñas, el error se dispara hasta límites excesivamente altos.
Esto es debido a que, al dividir entre números cada vez más pequeños, la ecuación
tiende a infinito. En el presente estudio, este problema ocurre a la hora de comparar
velocidades medias o desviaciones típicas, ya que conviven magnitudes de alto valor
absoluto (positivas y negativas) con magnitudes cercanas al cero.
|𝑚𝑎𝑔𝑛𝑖𝑡𝑢𝑑𝑟𝑒𝑎𝑙 − 𝑚𝑎𝑔𝑛𝑖𝑡𝑢𝑑𝑚𝑒𝑑𝑖𝑑𝑎 |
𝑅𝑃 𝐷 = * 100 (4.2)
(𝑚𝑎𝑔𝑛𝑖𝑡𝑢𝑑𝑟𝑒𝑎𝑙 + 𝑚𝑎𝑔𝑛𝑖𝑡𝑢𝑑𝑚𝑒𝑑𝑖𝑑𝑎 )/2
Entre los resultados de los casos experimentales que se pretende replicar en las
simulaciones se tomaron medidas de determinadas variables en ciertas localizaciones.
Sin embargo, cuando se trata de comparar de una manera eficiente unos casos con otros,
es necesario obtener un valor único de precisión representativo del caso, o al menos por
cada variable o conjunto de variables que se quiera estudiar de forma separada. Lo más
sencillo es realizar una media aritmética entre los porcentajes RPD de todas las medidas
de dicha variable, por lo que se ha elegido este método como manera de proceder. Por
ejemplo, si se quiere medir la velocidad axial en diferentes estaciones o posiciones de la
geometría, se calcula el RPD en cada punto en relación al resultado experimental, se
55
Memoria
promedia entre todos los puntos para cada estación, se promedian todas las estaciones
y se obtiene un escalar que, aunque no sea intuitivo, es comparable a otro obtenido de
igual manera.
Es importante destacar que este valor no es absoluto, es una magnitud relativa que
sirve para comparar en cierta variable unos casos con otros. Sin embargo no es correcto
comparar para un mismo caso distintas variables: si en un caso no premezclado, para la
velocidad media en las estaciones se obtiene un 60 % de error pero para la concentración
de especies se obtiene un 30 %, esto no significa que se haya conseguido un resultado
más fiel en la segunda que en la primera. Puede parecer que para una determinada
magnitud un RPD de (por ejemplo) 80 % sea algo desorbitado, pero este es resultado
de muchos promedios y no tiene por qué ser un error excesivo, cosa que se demuestra al
realizar gráficas de los perfiles y contornos de esas variables. También conviene recordar
que los resultados experimentales no están exentos de error, aunque en este trabajo se
tomen como referencia y se les suponga una precisión adecuada pueden alejarse de la
realidad.
56
4. METODOLOGÍA
Existen multitud de parámetros a la hora del mallado así como esquemas de discre-
tización y de turbulencia disponibles para el estudio del problema, por lo que al final
del capítulo se propone un estudio de independencia de malla como primer paso cuyo
objetivo es calibrar la precisión y duración de las simulaciones.
Las mallas se han creado con elementos tetraédricos, partiendo de un tamaño base
en la región del plenum y en la cámara, y refinando la región del swirler. Además,
tratando de capturar con mayor precisión la zona de recirculación, se ha refinado la
zona central de la cámara de combustión con esferas a diferentes distancias medidas
desde la entrada de la misma. Se puede ver un ejemplo de malla utilizada en la figura
4.11, con la estrategia mencionada.
57
Memoria
Además, en la tabla de la figura 4.12 están descritas todas las mallas utilizadas en
este apartado. La primera columna indica las distancias entre el eje de coordenadas y
el centro de cada esfera (excepto el tamaño base y el mallado concreto del swirler), la
primera fila es el número de celdas en millones de cada malla y la segunda fila indica
por columnas el tamaño de los elementos y el radio de las esferas de cada malla.
Para generar la malla estructurada por defecto, en primer lugar el programa sub-
divide el dominio en hexaedros teniendo en cuenta el tamaño base definido y el grid
scaling, que es utilizado para acelerar el establecimiento del flujo en el dominio, y así
también la convergencia del código. A continuación se obtiene la intersección entre es-
tas celdas y la triangulación de la geometría, consiguiendo una gran cantidad de celdas
cúbicas y un número menor de celdas irregulares (figura 4.14) correspondientes a los
bordes y superficies más complejos de la geometría. Después se aplica el parámetro de
fixed embedding en las zonas indicadas para refinar las celdas. Finalmente, el AMR se
tiene en cuenta durante la simulación, comparando los gradientes de las variables ele-
gidas con los límites establecidos para refinar en caso de sobrepasarlos. De esta manera
se consigue aumentar la resolución de la malla de manera local en secciones críticas
manteniendo elementos más gruesos en el resto del dominio.
58
4. METODOLOGÍA
59
Memoria
Aunque se han llevado a cabo una gran cantidad de simulaciones variando los valores
de intensidad turbulenta entre el 2 % y el 15 % en ambas entradas, y la escala de longitud
turbulenta entre un 5 % y un 12 % del diámetro, los mejores resultados se han obtenido
con la configuración de intensidad turbulenta del 2 % a la entrada del plenum y 10 %
en la línea de combustible.
𝜋𝑅𝑖
𝑡𝑟𝑒𝑠 = = 0.81 (4.3)
𝑢𝜃
En las simulaciones LES se debe refinar mucho más la malla, de modo que se consiga
obtener unos 𝑦 + lo suficientemente pequeños (menores de 4) como para que la ley de
pared se corresponda con la capa viscosa. Esto también permite capturar y calcular
con mayor precisión las estructuras del flujo más pequeñas, por eso a medida que se
refina la malla, mejores son los resultados.
60
4. METODOLOGÍA
Los parámetros más importantes que se han medido en los experimentos y simula-
ciones son la velocidad media (o mean velocity) y las medias cuadráticas o RMS (root
mean squared) en todas sus componentes: axial, radial y tangencial. En el caso no
premezclado también se ha medido la concentración del combustible.
Los valores RMS son una medida estadística de la magnitud de una cantidad va-
riable y representan la incertidumbre o el grado de dispersión entre unos datos y su
media. Se calculan como se indica en la ecuación 4.4 para la velocidad, mediante la
raíz cuadrada del promedio de todos los datos elevados al cuadrado (para eliminar los
signos). √︃
2 + 𝑣2 + 𝑣2
𝑣𝑎𝑥
𝑣𝑅𝑀 𝑆 = 𝑟𝑎𝑑 𝑡𝑎𝑛
(4.4)
3
61
Memoria
guardarse en ficheros pdf, tiff, u otras extensiones. Es importante mantener los casos
en un entorno de carpetas ordenado como se indica en los comentarios del código, y
asegurarse.
En las gráficas de las figuras 4.18 y 4.19 se muestran los porcentajes de RPD en
función del número de celdas de cada caso, tanto para la velocidad axial en la línea
63
Memoria
64
4. METODOLOGÍA
central como para la media de los perfiles de velocidades (axial, radial y tangencial) en
las estaciones.
Con estos datos, la malla más conveniente para realizar los estudios en OpenFOAM
del modelo 𝑘−𝜀 es la que se compone de entre 3 a 6 millones de celdas, si las limitaciones
de capacidad de cálculo son restrictivas. En caso contrario, se puede utilizar mallas de
17 millones de celdas o superior para el modelo RSM-LRR.
En las gráficas de las figuras 4.21 y 4.22 se muestran los porcentajes de RPD en
función del número de celdas de cada caso, igual que la sección anterior, tanto para
la velocidad axial en la línea central como para la media de los perfiles de velocidades
medias (axial, radial y tangencial) en las estaciones.
Como se puede observar en ellas, a medida que aumenta el número de celdas dis-
minuye la desviación de manera potencial. Además hay unas tendencias observables
según las características de la malla:
65
Memoria
66
4. METODOLOGÍA
67
Memoria
Ambas gráficas muestran una zona plana en general salvo pequeñas excepciones con
una desviación excesiva. Como la pendiente no es demasiado grande (excepto en las
mallas menos pesadas) la mejor malla en CONVERGE se obtiene entre 0.6-1.2 millones
de celdas.
68
Capítulo 5
Resultados
Las condiciones del caso premezclado han sido expuestas en el apartado 4.6 para
el estudio de independencia de malla, pero deben ser analizadas detenidamente y no
sólo mediante aproximación RANS, sino también mediante LES. La información que
se extraiga de este caso será de gran ayuda para analizar el caso siguiente.
La figura 5.1 indica que, para mallas de peso reducido, de todos los modelos de
turbulencia se obtienen resultados similares. En este caso destaca el modelo 𝑘 − 𝜀, que
tiene menor error que el resto. Sin embargo, a partir de los 3 millones de celdas se
estabiliza y pasa a ser el peor modelo con bastante diferencia sobre los demás. Los
otros modelos, a medida que aumenta el número de celdas, continúan disminuyendo
el error de forma significativa; no se alcanza una convergencia de malla porque son
modelos URANS/LES.
En la figura 5.2 se observa que todos los modelos tienen resultados similares. De
nuevo, como en la figura anterior, en mallas reducidas el modelo 𝑘−𝜀 obtiene resultados
ligeramente mejores, pero a partir de los 6 millones de celdas se diferencia del resto,
69
Memoria
70
5. RESULTADOS
Aunque las futuras simulaciones deban realizarse con las indicaciones que se han ido
dando en este apartado (por la capacidad computacional y tiempo de cálculo limitados),
para dibujar los perfiles y contornos se han escogido los casos más precisos de los ya
simulados, tratando de analizar las estructuras fluctuantes:
Entre las aproximaciones RANS, el número 18: con 17 millones de celdas y modelo
de turbulencia RSM-LRR (URANS).
Entre las aproximaciones LES, el número 17: con 17 millones de celdas y modelo
de turbulencia Dynamic Smagorinsky.
71
Memoria
72
5. RESULTADOS
y penetran en la misma con una longitud y diámetro determinados (figura 5.6). Con
ello se pretende capturar mejor las variables en esa zona.
La figura 5.7 refleja la comparativa entre los tres tipos. En cuanto a los dos primeros
no hay diferencia apreciable de desviación en línea central o en media de estaciones,
por lo que merece la pena un cilindro de longitud 1𝑐𝑚, refinando así la malla en un
lugar reducido y evitando que aumente el tiempo de simulación por la diferencia en
la cantidad de celdas. Una vez descartado el segundo, la comparativa entre el resto de
cilindros muestra muy poca diferencia, aunque favorable para el de mayor diámetro.
Sin embargo, no merece la pena duplicar el número de celdas y, por tanto, disparar el
tiempo de simulación por un cambio mínimo de la desviación, menor al 4 %. Además,
73
Memoria
74
5. RESULTADOS
en lo referente a las RMS sí que se reduce el error notablemente con el cilindro C1_3,
aunque sigue siendo alto.
A continuación, para el caso estudiado mediante LES también es necesario ver qué
modelo de turbulencia se ajusta más a los resultados experimentales. En la figura 5.8 se
ve una comparación entre los tres modelos: Dynamic Smagorinsky (DSGS), Dynamic
Structural (DSEM) y Smagorinsky-Lilly (SGS). En las velocidades de la línea central y
las estaciones, el modelo Dynamic Structural se desvía demasiado respecto de los otros
dos que dan resultados muy similares, por lo elegir uno u otro se debe decidir desde el
punto de vista del tiempo de simulación. Sin embargo, en cuanto a las RMS, todos los
modelos tienen unos resultados muy similares entre ellos, siendo ligeramente mejor el
modelo Dynamic Smagorinsky.
En esta ocasión, para las gráficas de perfiles y contornos se han seleccionado de nue-
vo dos casos, pero el primero sigue las indicaciones del apartado en cuanto a minimizar
el tiempo de simulación (con los modelos turbulentos y tamaños de malla adecuados)
mientras que el segundo se ha elegido por ser el más preciso para captar las estructuras
transitorias:
Entre los modelos RANS, el número 17: tamaño base de 3𝑚𝑚, 3 niveles de fixed
embedding y AMR, 1.2 millones de celdas y modelo de turbulencia 𝑘 − 𝜀 std.
Entre los modelos LES, el número 3: tamaño base de 2𝑚𝑚, 3 niveles de fixed
embedding y AMR, 17 millones de celdas y modelo de turbulencia Dynamic Sma-
gorinsky.
75
Memoria
25
Cordier2013EXP
URANS_18_OF
20
LES_17_OF
RANS_17_CG
Mean Axial Velocity [m/s]
15 LES_3_CG
10
-5
-10
-15
0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 0.05
Centerline Distance z [m]
La figura 5.9 muestra los resultados de velocidad media axial en la línea central. De
todos los casos, el más cercano a los datos experimentales es el LES de CONVERGE,
sin embargo hay una pequeña particularidad: todos ellos sobreestiman la velocidad en
las primeras distancias, especialmente en el origen. Esto puede deberse a que el modelo
de pared en el inyector no actúa correctamente (los 𝑦 + se encuentran en la capa buffer,
zona de mayor error al modelar), dando como resultado una mayor velocidad axial (en
torno a 20𝑚/𝑠) a la salida del mismo (en el centro de la entrada a la cámara). En
cualquier caso, los autores de las medidas experimentales reportaron dificultades para
realizar las mediciones en esa zona [70, 77, 78, 79].
En las figuras 5.10-5.12 se ven las tres componentes axial, radial y tangencial res-
pectivamente de las velocidades medias (a la izquierda) y las RMS (a la derecha) en
función de la distancia radial tomada desde la línea central de la cámara en cada
estación.
En primer lugar, la velocidad media axial es sobreestimada por todos los modelos en
las dos primeras estaciones, mientras que en el resto ocurre solamente con los modelos
76
5. RESULTADOS
Station z = 40mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 30mm
15
20
10
Mean Axial Velocity [m/s]
10
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 20mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 10mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 5mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
LES y cerca de la zona central. Esto indica que la zona de recirculación que se mostrará
en los contornos de simulación es más estrecha que la experimental, y ocurre acorde
a lo explicado en la figura anterior: existe una velocidad axial mayor de la esperada
en la entrada. En cuanto a las RMS, como ya se ha adelantado previamente en este
capítulo, el modelo RANS_17_CG da los peores resultados subestimando los valores.
En el resto de modelos llama la atención que alrededor de la zona 0.01𝑚 de distancia
radial exista un pico que los datos experimentales no reflejan.
77
Memoria
Station z = 40mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 30mm
15 10
10
Mean Radial Velocity [m/s]
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 20mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 10mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 5mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
en comparación al resto, aunque los resultados obtenidos por ellos tampoco logran
capturar correctamente las estructuras turbulentas cerca del centro.
78
5. RESULTADOS
Station z = 40mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 30mm
10
Mean Tangential Velocity [m/s]
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 20mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 10mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 5mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
79
Memoria
0.1 0.1
0.08 0.08
Axial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-5 0 5 2 4 6 8 10 12
0.1 0.1
Distance z [m]
0.08 0.08
Radial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-20 -10 0 10 20 2 4 6 8 10 12
0.1 0.1
Distance z [m]
Tangential
0.08 0.08
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
Distance x [m] Distance x [m]
80
5. RESULTADOS
del flujo. Además, esta estructura (descrita en 4.1) que se muestra en la componente
axial de velocidades medias, es de gran importancia por su función estabilizadora de
llama gracias al bloqueo aerodinámico formado. Su presencia dota al flujo del tiempo de
residencia, temperatura y turbulencia adecuadas para una atomización, evaporación,
mezclado y combustión eficientes.
0.09
120
0.08
100
0.07
Distance z [m]
0.06
80
0.05
60
0.04
0.03 40
0.02
20
0.01
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Distance x [m]
Otro método para identificar los vórtices es el Q-criterion [84], o segundo invariante
del tensor gradiente de velocidades. Separando las partes simétrica 𝑆 y antisimétrica
Ω del tensor, el Q-Criterion se establece de la siguiente manera:
1 1
∇𝑣 = 𝑆 + Ω donde 𝑆 = (∇𝑣 + ∇𝑣 𝑇 ), Ω = (∇𝑣 − ∇𝑣 𝑇 ) (5.1)
2 2
Para visualizar los vórtices se obtiene una isosuperficie de valores 𝑄 > 0 obtenidos de
la ecuación 5.2. Esto se observa en la figura 5.15 obtenida con el programa de pospro-
cesado EnSight, y ocurre porque las isosuperficies con 𝑄 positivas aislan áreas donde
81
Memoria
En el caso anterior (sección 5.1) se ha extraído información útil con la cual encarar
el resto de simulaciones del siguiente, como el tamaño de la malla debido a los niveles
de fixed embedding y AMR. De esta manera se puede continuar con el estudio del flujo
no premezclado en la cámara de combustión LDI.
En la gráfica de las velocidades (figura 5.16) se advierte que no hay grandes cambios
entre un modelo u otro, y es favorable para el dominio inicializado con aire y mezcla
83
Memoria
excepto para las RMS. Por otra parte, la gráfica de especies (figura 5.17) refleja unas
RMS similares y una diferencia muy importante en la media. Esta última razón hace que
sea imperioso elegir el dominio inicializado con aire en plenum y swirler, ya que además
no supone un perjuicio apreciable en otros aspectos como el tiempo de simulación.
40
Cordier2013EXP
35 Esclapez2014AVBP
URANS_15_CG
30
Mean Axial Velocity [m/s]
25
20
15
10
-5
-10
-15
0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 0.05
Centerline Distance z [m]
85
Memoria
La figura 5.20 muestra los resultados de velocidad media axial en la línea central.
A partir de los 15𝑚𝑚 se asemejan mucho los datos de la simulación con los de la
medición, pero en la zona más cercana existe una diferencia bastante grande. Esto se
debe, en parte, a los problemas que ocasionó la medición experimental en estos puntos,
que ya ocurrían en el caso premezclado pero se ven magnificados por la densidad del
combustible. En este sentido, el código Esclapez2014AVBP da unos resultados más
acordes con la velocidad teórica, y se parece más al caso experimental desde los 5𝑚𝑚.
En las figuras 5.21-5.23 se ven las tres componentes axial, radial y tangencial res-
pectivamente de las velocidades medias (a la izquierda) y las RMS (a la derecha) en
función de la distancia radial tomada desde la línea central de la cámara en cada
estación.
Station z = 40mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 30mm
15
20
10
Mean Axial Velocity [m/s]
10
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 20mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 10mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 5mm
15
20
10 10
0 5
-10 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
86
5. RESULTADOS
debido a su mayor densidad. A partir de 20𝑚𝑚 ya se asemejan más los resultados, pero
la zona de recirculación es más estrecha de lo que debería ser. Ocurre lo mismo con
los resultados de Esclapez2014AVBP, pero con un pico menos pronunciado y adaptán-
dose en una estación anterior (10𝑚𝑚). Las RMS del caso URANS_15_CG consiguen
acercarse de manera aceptable a la tendencia de los resultados experimentales, ex-
ceptuando la misma zona ya nombrada, y subestimando las medidas en las últimas
estaciones, mientras que el AVBP los sobreestima siempre. Aunque los resultados de
AVBP son algo más ajustados, los conseguidos en este proyecto no tienen nada que
envidiar teniendo en cuenta que hay una gran diferencia en el coste computacional.
Station z = 40mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 30mm
15 10
10
Mean Radial Velocity [m/s]
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 20mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 10mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 5mm
15 10
10
5 5
-5 0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
En segundo lugar, la velocidad media radial de ambos códigos CFD se ajusta casi
perfectamente en todas las estaciones. En cuanto a las RMS, de nuevo los resultados
del código AVBP son más precisos que los obtenidos en CONVERGE, los primeros
sobreestimando y los segundos subestimando los valores en relación a los resultados
87
Memoria
Station z = 40mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 30mm
10
Mean Tangential Velocity [m/s]
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 20mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 10mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Station z = 5mm
10
0
-10 5
-20
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
88
5. RESULTADOS
Station z = 22mm
0.6 0.2
0.4
0.1
0.2
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Mean CH4 Mass Fraction [-]
Station z = 14mm
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Station z = 10mm
0.6 0.2
0.4
0.1
0.2
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Station z = 6mm
0.6 0.2
0.4
0.1
0.2
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Station z = 22mm
0.2
0.15
0.1 0.1
0.05
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Mean CH4 Mass Fraction [-]
Station z = 14mm
RMS CH4 Mass Fraction [-]
0.2
0.15
0.1 0.1
0.05
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Station z = 10mm
0.2
0.15
0.1 0.1
0.05
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Station z = 6mm
0.2
0.15
0.1 0.1
0.05
0 0
-0.03 -0.02 -0.01 0 0.01 0.02 0.03
Figura 5.25: Figura 5.24 con zoom en la escala para mejor visuali-
zación.
89
Memoria
90
5. RESULTADOS
0.1 0.1
0.08 0.08
Axial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-5 0 5 2 4 6 8 10
0.1 0.1
Distance z [m]
0.08 0.08
Radial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-20 -10 0 10 20 2 4 6 8
0.1 0.1
Distance z [m]
Tangential
0.08 0.08
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
Distance x [m] Distance x [m]
91
Memoria
0.09
80
0.08
70
0.07
Distance z [m]
60
0.06
50
0.05
40
0.04
30
0.03
0.02 20
0.01 10
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Distance x [m]
92
Capítulo 6
Conclusiones y trabajos
futuros
93
Memoria
Las simulaciones del caso premezclado en CONVERGE han reportado muy buenos
resultados, haciendo evidente que el método de mallado adaptativo es una herramienta
potente capaz de reducir los tiempos de simulación manteniendo una precisión adecua-
da, a la altura de un código con mallado tradicional como OpenFOAM. Para el caso no
premezclado también se han obtenido grandes resultados con aproximación URANS,
en comparación a los que ofrecía AVBP (con mayor resolución y aproximación LES),
aunque se deja abierta la puerta para explorar una configuración que capte mejor las
velocidades axiales y, por ende, las concentraciones de combustible en la zona más
cercana a la entrada de la cámara.
Este TFM puede suponer, por tanto, el punto de partida de una investigación
con objetivos más ambiciosos, modelando la totalidad de los fenómenos inherentes
al proceso de inyección de combustible líquido: atomización, evaporación, mezclado e
interacción de las estructuras bifásicas, la combustión y las emisiones contaminantes.
94
Bibliografía
[1] ICAO. ICAO Environmental Report 2016. Cambridge University Press, page 250,
2016.
[2] Philippe Busquin, Pedro Argüelles, Manfred Bischoff, B.A.C. Droste, SirRichardH.
Evans, Walter Kröll, Jean-Luc Lagardère, Alberto Lina, John Lumsden, Denis
Ranque, Søren Rasmussen, Paul Reutlinger, SirRalph Robins, Helena Terho, and
Arne Wittlöw. European aeronautics: a vision for 2020. Air & Space Europe,
3(3-4):16–18, 2001.
[3] Joyce E Penner. Aviation and the Global Atmosphere: a special report of IPCC
Working Groups I and III. page 373, 1999.
[5] Scott C. Herndon, Joanne H. Shorter, Mark S. Zahniser, David D. Nelson, John
Jayne, Robert C. Brown, Richard C. Miake-Lye, Ian Waitz, Phillip Silva, Tho-
mas Lanni, Ken Demerjian, and Charles E. Kolb. NO and NO2 emission ratios
measured from in-use commercial aircraft during taxi and takeoff. Environmental
Science and Technology, 38(22):6078–6084, 2004.
[8] D. W. Bahr, F. Culick, M.V. Heitor, and J.H. Whitelaw. Unsteady combustion.
In Aircraft Turbine Engine NOx Emission Abatement, pages 243–264. Springer,
Dordrecht, 1996.
[10] A Seaton, W MacNee, K Donaldson, and D Godden. Particulate air pollution and
acute health effects. Lancet, 345(8943):176–178, 1995.
[11] Arthur H. Lefebvre and Dilip R. Ballal. Gas turbine combustion: alternative fuels
and emissions. Taylor & Francis, Boca Raton, FL, 3rd editio edition, 2010.
[13] K. K. Rink and A. H. Lefebvre. The influences of fuel composition and spray
characteristics on nitric oxide formation. Combustion Science and Technology,
68(1-3):1–14, 1989.
[14] G Leonard and S Correa. Second asme fossil fuel combustion symposium. AS-
ME/PD, 30:69–74, 1990.
[15] Chi-Ming Lee, Kathleen M. Tacina, and Changlie Wey. High pressure low NOx
emissions research: Recent progress at NASA Glenn Research Center. pages 1–8,
2007.
[16] S. Mosier and R. Pierce. Advanced combustion systems for stationary gas turbine
engines: Volume i. review and preliminary evaluation., 1980.
[18] D. W. Bahr. Technology for the design of high temperature rise combustors.
Journal of Propulsion and Power, 3(2):179–186, Mar 1987.
[20] Zhuohui J. He, Kristin Kopp-Vaughan, Changlie Wey, Clarence T. Chang, Al-
bert K. Cheung, Chi-Ming Lee, and Angela D. Surgenor. chapter Emission Cha-
racteristics of A P&W Axially Staged Sector Combustor. AIAA SciTech Forum.
American Institute of Aeronautics and Astronautics, Jan 2016. 0.
[23] J. Cai. Aerodynamics of lean direct injection combustor with multi-swirler arrays.
PhD thesis, University of Cincinnati, 2006.
[24] H. S. Alkabie, G. E. Andrews, and N. T. Ahmad. Lean low nox primary zones
using radial swirlers, 1988.
96
BIBLIOGRAFÍA
[25] H. S. Alkabie and G. E. Andrews. Ultra low nox emissions for gas and liquid fuels
using radial swirlers, 1989.
[26] H. S. Alkabie and G. E. Andrews. Radial swirlers with peripheral fuel injection
for ultra-low nox emissions, 1990.
[27] H. S. Alkabie, G. E. Andrews, and M. N. Kim. An ultra low nox pilot combustor
for staged low nox combustion. 11th ISABE-1993-7020, 1:215–225, 1993.
[28] Waldemar Lazik, Thomas Dörr, and Sebastian Bake. Low nox combustor deve-
lopment for the e3e core engine demonstrator, 09 2007.
[29] Daniel A. Nickolaus, D. Scott Crocker, David L. Black, and Clifford E. Smith.
Development of a lean direct fuel injector for low emission aero gas turbines, 2002.
[30] S. W. Shaffar and G. S. Samuelsen. A liquid fueled, lean burn, gas turbine com-
bustor injector. 139:41–57, 10 1998.
[31] T. Terasaki. Lean non-premixed combustion for low-nox gas turbine combustor.
Proceedings of the 1995 Yokohama International Gas Turbine Congress(II), pages
353–358, 1995.
[32] Robert Tacina. Low NOx Potential of Gas Turbine Engines. 28th Aerospace
Sciences Meeting, page 20, 1990.
[34] Robert Tacina, Chien-Pei Mao, and Changlie Wey. Experimental investigation
of a multiplex fuel injector module with discrete jet swirlers for low emission
combustors. 09 2004.
[35] R.R. Tacina, P. Lee, and C. Wey. A lean direct injection combustor using a 9
point swirl-venturi fuel injector, 2005.
[36] D. Dewanji and A. G. Rao. Spray combustion modeling in lean direct injec-
tion combustors, part II: Multi-point LDI. Combustion Science and Technology,
187(4):558–576, 2015.
[38] J.M. Beér and N.A. Chigier. Combustion aerodynamics. Fuel and energy science
series. Applied Science Publishers Ltd, 1972.
[39] Erol Kilik. The influence of swirler design parameters on the aerodynamics of
downstream recirculation region. PhD thesis, 1976.
97
BIBLIOGRAFÍA
[40] Robert Tacina, Jun Cai, and San-Mou Jeng. The Structure of a Swirl-Stabilized
Reacting Spray Issued From an Axial Swirler. 43rd AIAA Aerospace Sciences
Meeting & Exhibit, (January), 2005.
[43] S. Etemad and B. C. Forbes. Cfd modeling of a gas turbine combustor air swirler
effect, computational fluid dynamics in aero-propulsion, 1995.
[44] Ashraf A. Ibrahim and Milind A. Jog. Nonlinear breakup model for a liquid
sheet emanating from a pressure-swirl atomizer. Journal of Engineering for Gas
Turbines and Power, 129(4):945–953, Jan 2007.
[46] John D Anderson. Computational Fluid Dynamics: The Basics with Applications,
1995.
[47] Hernan Tinoco, Hans Lindqvist, Wiktor Frid, and Kraftgrupp Ab. Numerical
Simulation of Industrial Flows. 1990.
[49] H. Tennekes and J.L. Lumley. A First Course in Turbulence. MIT Press, 1972.
[51] B.E. Launder and B.I. Sharma. Application of the energy-dissipation model of
turbulence to the calculation of flow near a spinning disc. Letters in Heat and
Mass Transfer, 1(2):131 – 137, 1974.
[52] W.P. Jones and J.H. Whitelaw. Calculation methods for reacting turbulent flows:
A review. Combustion and Flame, 48:1 – 26, 1982.
98
BIBLIOGRAFÍA
[54] Wolfgang Rodi. Turbulent models and their application in hydraulics-a state of
art review. 12 1980.
[55] K. Hanjalic. Will rans survive les?: A view of perspectives. Journal of Fluids
Engineering, 127(5):831–839, Sep 2005.
[56] Stefan H Johansson, Lars Davidson, and Erik Olsson. Numerical simulation of vor-
tex shedding past triangular cylinders at high reynolds number using ak-𝜀 turbu-
lence model. International Journal for Numerical Methods in Fluids, 16(10):859–
878, 1993.
[57] Sven Perzon and Lars Davidson. On CFD and transient flow in vehicle aerodyna-
mics. SAE World Congress, (724):2000–01–0873, 2000.
[58] Gopal Patnaik, Jay P Boris, Theodore R Young, and Fernando F Grinstein. Lar-
ge scale urban contaminant transport simulations with miles. Journal of Fluids
Engineering, 129(12):1524–1532, 2007.
[59] Joel H. Ferziger. Large eddy numerical simulations of turbulent flows. AIAA
Journal, 15(9):1261–1267, Sep 1977.
[61] J.O. Hinze. Turbulence. McGraw-Hill classic textbook reissue series. McGraw-Hill,
1975.
[62] Massimo Germano, Ugo Piomelli, Parviz Moin, and William H. Cabot. A dyna-
mic subgrid?scale eddy viscosity model. Physics of Fluids A: Fluid Dynamics,
3(7):1760–1765, jul 1991.
[63] S.V Patankar and D.B Spalding. A calculation procedure for heat, mass and
momentum transfer in three-dimensional parabolic flows. International Journal
of Heat and Mass Transfer, 15(10):1787–1806, 1972.
[65] R.I Issa. Solution of the implicitly discretised fluid flow equations by operator-
splitting. Journal of Computational Physics, 62(1):40 – 65, 1986.
[66] Hrvoje Jasak. Error analysis and estimation for the finite volume method with
applications to fluid flows. 1996.
99
BIBLIOGRAFÍA
[70] Matthieu Cordier. Allumage et propagation de flamme dans les écoulements for-
tement swirlés : études expérimentales et numériques. PhD thesis, 2013.
[74] Convergent Science. Advanced converge cfd training 2.4 gas turbine combustion,
December 2017.
[75] NJDEP. Laboratory Quality Assurance and Quality Control Data Quality Assess-
ment and Data Usability Evaluation. (April), 2014.
[76] Lucas Esclapez. Numerical study of ignition and inter-sector flame propagation
in gas turbine. PhD thesis, Institut National Polytechnique de Toulouse (INP
Toulouse), 2015.
[77] David Barré, Matthias Kraushaar, Gabriel Staffelbach, Vincent Moureau, and
Laurent Y.M. Gicquel. Compressible and low Mach number LES of a swirl expe-
rimental burner. 3rd INCA Colloquim, 341(1-2):277–287, 2013.
[78] David Barre. Simulation Numerique De L’Allumage Dans Les Chambres De Com-
bustion Aeronautiques. PhD thesis, 2014.
[79] David Barré, Lucas Esclapez, Matthieu Cordier, Eleonore Riber, Bénédicte Cue-
not, Gabriel Staffelbach, Bruno Renou, Alexis Vandel, Laurent Y.M. Gicquel, and
Gilles Cabot. Flame propagation in aeronautical swirled multi-burners: Experi-
mental and numerical investigation. Combustion and Flame, 161(9):2387–2405,
2014.
100
BIBLIOGRAFÍA
[82] S Leibovich. The structure of vortex breakdown. Annual Review of Fluid Mecha-
nics, 10(1):221–246, 1978.
[83] Nicholas Syred. A review of oscillation mechanisms and the role of the precessing
vortex core (pvc) in swirl combustion systems. Progress in Energy and Combustion
Science, 32(2):93 – 161, 2006.
101
Apéndices
Apéndice A
Perfiles y contornos de
velocidades y energía cinética
turbulenta
0.09 30
0.08
25
0.07
Distance z [m]
0.06 20
0.05
15
0.04
0.03 10
0.02
5
0.01
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Distance x [m]
105
Apéndices
0.1 0.1
0.08 0.08
Axial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-5 0 5 2 4 6
0.1 0.1
Distance z [m]
0.08 0.08
Radial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-20 -10 0 10 20 2 4 6 8
0.1 0.1
Distance z [m]
Tangential
0.08 0.08
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
Distance x [m] Distance x [m]
106
A. PERFILES Y CONTORNOS DE VELOCIDADES Y ENERGÍA CINÉTICA
TURBULENTA
50
0.09
45
0.08
40
0.07
35
Distance z [m]
0.06
30
0.05
25
0.04
20
0.03
15
0.02 10
0.01 5
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Distance x [m]
107
Apéndices
0.1 0.1
0.08 0.08
Axial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-5 0 5 0 2 4 6
0.1 0.1
Distance z [m]
0.08 0.08
Radial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-20 -10 0 10 20 2 4 6
0.1 0.1
Distance z [m]
Tangential
0.08 0.08
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
Distance x [m] Distance x [m]
108
A. PERFILES Y CONTORNOS DE VELOCIDADES Y ENERGÍA CINÉTICA
TURBULENTA
0.09
25
0.08
0.07
20
Distance z [m]
0.06
0.05 15
0.04
10
0.03
0.02
5
0.01
0
-0.05 -0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05
Distance x [m]
109
Apéndices
0.1 0.1
0.08 0.08
Axial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-5 0 5 2 4 6 8
0.1 0.1
Distance z [m]
0.08 0.08
Radial
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
-20 -10 0 10 20 2 4 6
0.1 0.1
Distance z [m]
Tangential
0.08 0.08
0.06 0.06
0.04 0.04
0.02 0.02
0 0
-0.05 0 0.05 -0.05 0 0.05
Distance x [m] Distance x [m]
110
Apéndice B
B.1. LANZA
codigo/lanza.m
% C a r l o s Moreno Montagud 21/02/2018
%
% Folder s t r u c t u r e :
% \ ’ code ’ \ ’ geometry ’ \ ’ i n j e c t i o n ’ \ ’ f l u i d ’ \ ’ i n l e t T y p e ’ \ c a s e s \ ’
case ’ \ ’ f o l d e r ’ \
% \ mats \
% \ Literature \ " \ " \ " \ " \ cases \ ’
author ’ \
% \ mats \
% \ Results \ " \ " \ " \ " \’
solveMode ’ \ ’ code ’ \
%
% I m p o r t i n g CONVERGE:
% l a n z a .m
% − D e f i n e m i l i S e c (3 d i g i t s ) m i l i S e c ={’___’ } ;
% − Void mesh mesh ={ ’ ’};
% i m p o r t P r o p e r t i e s .m
% − Name o f c o l f i l e s c o l _ f i l e s ={’ v e l o c i t i e s ’ , ’
stress ’ , ’ species ’};
%
% I m p o r t i n g OpenFOAM:
% − T r a n s i e n t c a s e s ’ i n s t a n t f o l d e r must c o n t a i n a t l e a s t 3
digits after
111
Apéndices
%% S c r i p t
clearvars ;
clc ;
addpath ( genpath ( ’ f u n c t i o n s ’ ) ) ;
%% I n i t i a l i z e v a r i a b l e s
% Case v a r i a b l e s
c a s e s . code = { ’OpenFOAM ’ ; ’OpenFOAM ’ ; ’CONVERGE’ ; ’CONVERGE’ } ; %
’CONVERGE’ , ’OpenFOAM’ , ’ L i t e r a t u r e ’
c a s e s . geometry = { ’ KIAI ’ } ; % ’NASA’ , ’ KIAI ’
c a s e s . solveMode = { ’RANS ’ ; ’LES ’ ; ’RANS ’ ; ’LES ’ } ; % ’RANS’ , ’LES ’
c a s e s . i n j e c t i o n = { ’ premixed ’ } ; % ’ nonPremixed ’ , ’ premixed ’
c a s e s . f l u i d = { ’ gas ’ } ; % ’ gas ’ , ’ l i q u i d ’
c a s e s . i n l e t T y p e = { ’ nonReactive ’ } ; % ’ nonReactive ’ , ’ r e a c t i v e ’
c a s e s . numCase = { ’ 3 ’ ; ’ 2 ’ ; ’ 17 ’ ; ’ 3 ’ } ;
c a s e s . m i l i S e c = { ’ 119 ’ ; ’ 128 ’ ; ’ 433 ’ ; ’ 215 ’ } ; % ’ ’ , ’ 1 6 5 ’ , ’ 1 8 5 ’
c a s e s . mesh = { ’ 17M’ ; ’ 17M’ ; ’ ’ ; ’ ’ } ; % ’ ’ , ’6M’ , ’17M’
c a s e s . author = { ’ Cordier2013EXP ’ } ; % ’ JunCai ’ , ’ PatelMenon ’ , ’
Cordier2013EXP ’ , ’ Esclapez2014AVBP ’
% Post−p r o c e s s i n g f l a g s
f l a g s . plot_flag = true ;
% Import v a r i a b l e s
112
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
import = i m p o r t P r o p e r t i e s ( ) ;
% Plot variables
plotting = plotProperties () ;
% Path names
names = path_names ( c a s e s , p l o t t i n g . l e g e n d _ f l a g s ) ;
%% Load c a s e s
[ p l o t t i n g .EXP, p l o t t i n g . SIM ] = l o a d C a s e s ( names , import , f l a g s ) ;
%% Post−p r o c e s s i n g
i f flags . plot_flag
% Variables
p l o t t i n g . legend_auto_exp = names . legend_exp ;
p l o t t i n g . legend_auto_sim = names . legend_sim ;
p l o t t i n g . r e s u l t s _ f o l d e r = names . r e s u l t s _ f o l d e r ;
% Plot
num_comp = p l o t C a s e ( p l o t t i n g ) ;
end
restoredefaultpath ;
B.2. IMPORTPROPERTIES
codigo/importProperties.m
% 16/05/2018
% C a r l o s Moreno Montagud
% import . f i l e s r e f e r s t o t h e name o f t h e f i l e s . c o l e x p o r t e d
from CONVERGE.
% import . ∗ . v a r i a b l e s r e f e r s t o t h e v a r i a b l e s you can o b t a i n
from t h a t f i l e .
% import . ∗ . f i l e _ e n d r e f e r s t o t h e f i n a l p a r t o f t h e f i l e n a m e
where t h e d a t a
% w i l l be read when i m p o r t i n g EXP f i l e .
% import . ∗ . p l o t _ t y p e r e f e r s t o t h e way o f i m p o r t i n g and
showing t h e
113
Apéndices
% i n f o r m a t i o n : a s i n g l e l i n e , many s t a t i o n s a l o n g z a x i s or
contours .
% import . ∗ . ∗ . s t a t i s t i c i s i g n o r e d when i m p o r t i n g OpenFOAM c a s e
.
%
% The r e s t o f s u b s t r u c t u r e s are s e l f −e x p l a i n e d by name .
% Unused v a r i a b l e s w i l l be s k i p p e d so t h e s c r i p t won ’ t r e t u r n
errors , but
% can be commented t o s a v e some memory . Don ’ t d e l e t e them as
t h e y are good
% e x a m p les f o r new v a r i a b l e s and s t r u c t u r e s .
function import = i m p o r t P r o p e r t i e s ( )
%% General
% Geometry and meshing p r o p e r t i e s
import . b a s e _ s i z e = 0 . 0 0 3 ; % [m]
import .AMR = 3 ; % l e v e l o f r e f i n e m e n t
import . x l i m s = [ −0.05 0 . 0 5 ] ; % cc w i d t h [m]
import . y l i m s = [ −0.05 0 . 0 5 ] ; % cc h e i g h t [m]
% import . z l i m s = [ ] ;
%%CONVERGE f i l e s
import . c o l _ f i l e s = { ’ v e l o c i t i e s ’ , ’ s t r e s s ’ } ; % ’ v e l o c i t i e s ’ , ’
stress ’ , ’ species ’
import . o u t _ f i l e s = { ’ monitor ’ } ;
% Import v a r i a b l e s
import . v e l o c i t i e s _ f i l e . v a r i a b l e s = { ’ v e l o c i t i e s ’ , ’ t k e ’ } ;
import . s t r e s s _ f i l e . v a r i a b l e s = { ’ s t r e s s ’ } ;
import . s p e c i e s _ f i l e . v a r i a b l e s = { ’ s p e c i e s ’ } ;
import . m o n i t o r _ f i l e . v a r i a b l e s = { ’ p r e s s u r e ’ } ;
import . m o n i t o r _ f i l e . f s = 1 0 0 0 0 ; %
%import . m o n i t o r _ f i l e . p o i n t s = 2 ; % % % % % % %
%% OpenFOAM f i l e s
import . c r d _ f i l e s = { ’ ccx ’ , ’ ccy ’ , ’ c c z ’ } ;
import . v a r _ f i l e s = { ’UMean ’ , ’ UPrime2Mean ’ } ;
% Import v a r i a b l e s
114
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
import . UMean . v a r i a b l e s = { ’ v e l o c i t i e s ’ } ;
import . UMean . numCols = 3 ;
import . UMean . impCols = [ 3 1 2 ] ; % a x i a l , r a d i a l , t a n g e n t i a l (
in order )
import . UPrime2Mean . v a r i a b l e s = { ’ v e l o c i t i e s ’ , ’ t k e ’ } ;
import . UPrime2Mean . numCols = 6 ;
import . UPrime2Mean . impCols = [ 6 1 4 ] ; % a x i a l , r a d i a l ,
t a n g e n t i a l ( in order )
%%
% Velocities
import . v e l o c i t i e s . f i l e _ e n d = { ’ v e l o c i t y ’ } ;
import . v e l o c i t i e s . plot_type = { ’ c e n t e r l i n e ’ , ’ s t a t i o n s ’ , ’
contours ’ } ; % ’ c e n t e r l i n e ’ , ’ s t a t i o n s ’ , ’ contours ’
% Centerline
import . v e l o c i t i e s . c e n t e r l i n e . p l a n e = { ’ z ’ } ; %
import . v e l o c i t i e s . c e n t e r l i n e . s t a t i s t i c = { ’ mean ’ } ; % ’ mean
’ , ’RMS’
import . v e l o c i t i e s . c e n t e r l i n e . d i r e c t i o n = { ’ a x i a l ’ } ; % ’
axial ’ , ’ radial ’ , ’ tangential ’
import . v e l o c i t i e s . c e n t e r l i n e . l i m i t s = [ 0 0 . 0 5 ] ; % [ 0 0 . 0 5 ]
% Stations
import . v e l o c i t i e s . s t a t i o n s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
import . v e l o c i t i e s . s t a t i o n s . s t a t i s t i c = { ’ mean ’ , ’RMS ’ } ; % ’
mean ’ , ’RMS’
import . v e l o c i t i e s . s t a t i o n s . d i r e c t i o n = { ’ a x i a l ’ , ’ r a d i a l ’ , ’
tangential ’ }; % ’ axial ’ , ’ radial ’ , ’ tangential ’
import . v e l o c i t i e s . s t a t i o n s . s t a t i o n s = [ 5 10 20 30 4 0 ] ; %
[ 5 15 29 46 76 9 2 ] , [ 5 10 20 30 4 0 ]
% Contours
import . v e l o c i t i e s . c o n t o u r s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
import . v e l o c i t i e s . c o n t o u r s . s t a t i s t i c = { ’ mean ’ , ’RMS ’ } ; % ’
mean ’ , ’RMS’
import . v e l o c i t i e s . c o n t o u r s . d i r e c t i o n = { ’ a x i a l ’ , ’ r a d i a l ’ , ’
tangential ’ }; % ’ axial ’ , ’ radial ’ , ’ tangential ’
import . v e l o c i t i e s . c o n t o u r s . l i m i t s = [ 0 0 . 1 ] ; % [ −0.05 0 . 0 5
0 0.1]
% Tke
import . t k e . f i l e _ e n d = { ’ t u r b u l e n t ’ , ’ k i n e t i c ’ , ’ e n e r g y ’ } ;
115
Apéndices
import . t k e . plot_type = { ’ c o n t o u r s ’ } ; % ’ c e n t e r l i n e ’ , ’ s t a t i o n s
’ , ’ contours ’
%Contours
import . t k e . c o n t o u r s . p l a n e = { ’ xz ’ } ;
import . t k e . c o n t o u r s . s t a t i s t i c = { ’RMS ’ } ; % ’RMS’
import . t k e . c o n t o u r s . d i r e c t i o n = { ’ a x i a l ’ , ’ r a d i a l ’ , ’
tangential ’ }; % ’ axial ’ , ’ radial ’ , ’ tangential ’
import . t k e . c o n t o u r s . l i m i t s = [ 0 0 . 1 ] ; % [ −0.05 0 . 0 5 0 0 . 1 ]
% Stress
import . s t r e s s . f i l e _ e n d = { ’ s t r e s s ’ } ;
import . s t r e s s . plot_type = { ’ c o n t o u r s ’ } ; % ’ c e n t e r l i n e ’ , ’
s t a t i o n s ’ , ’ contours ’
% Contours
import . s t r e s s . c o n t o u r s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
import . s t r e s s . c o n t o u r s . d i r e c t i o n = { ’ s11 ’ , ’ s12 ’ , ’ s13 ’ , ’ s22
’ , ’ s23 ’ , ’ s33 ’ } ; % ’ s11 ’ , ’ s12 ’ , ’ s13 ’ . . .
import . s t r e s s . c o n t o u r s . l i m i t s = [ 0 . 0 0 1 0 . 1 ] ; % [ −0.05 0 . 0 5
0 0.1]
% Species
import . s p e c i e s . f i l e _ e n d = { ’ mass ’ , ’ f r a c t i o n ’ } ;
import . s p e c i e s . plot_type = { ’ s t a t i o n s ’ } ; % ’ c e n t e r l i n e ’ , ’
s t a t i o n s ’ , ’ contours ’
% Stations
import . s p e c i e s . s t a t i o n s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
import . s p e c i e s . s t a t i o n s . s t a t i s t i c = { ’ mean ’ , ’RMS ’ } ; % ’
mean ’ , ’RMS’
import . s p e c i e s . s t a t i o n s . d i r e c t i o n = { ’CH4 ’ } ; %
import . s p e c i e s . s t a t i o n s . s t a t i o n s = [ 6 10 14 2 2 ] ; % [ 6 10
14 2 2 ]
end
B.3. PLOTPROPERTIES
codigo/plotProperties.m
% 16/05/2018
% C a r l o s Moreno Montagud
116
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
function p l o t t i n g = p l o t P r o p e r t i e s ( )
%% F l a g s
plotting . lines_flag = true ;
plotting . stations_flag = true ;
plotting . contours_flag = true ;
p l o t t i n g . monitor_flag = f a l s e ;
p l o t t i n g . compare_flag = t r u e ;
%% Save
plotting . save_flag = true ;
p l o t t i n g . save_format = { ’ pdf ’ } ; % ’ pdf ’ ’ epsc ’ ’ png ’ t i f f ’
%% Legend
p l o t t i n g . legend_mode = ’ manual ’ ; % ’ auto ’ ’ manual ’
%%%%%%%%%%%%%%%%%%%%%%%%%
plotting . legend_flags = [1 0 0 0 1 1 1 0 ] ; % boolean values [
KR NP G NR X XXX CG 6M]
p l o t t i n g . legend_manual_exp = { ’ Cordier2013EXP ’ } ;
p l o t t i n g . legend_manual_sim = { ’URANS_18_OF ’ , ’LES_17_OF ’ , ’
RANS_17_CG ’ , ’LES_3_CG ’ } ;
%% C e n t e r l i n e
p l o t t i n g . c e n t e r l i n e . type = { ’ v e l o c i t i e s ’ } ;
p l o t t i n g . c e n t e r l i n e . v e l o c i t i e s . plane = { ’ z ’ };
p l o t t i n g . c e n t e r l i n e . v e l o c i t i e s . s t a t i s t i c = { ’ mean ’ } ;
plotting . centerline . v e l o c i t i e s . directions = { ’ axial ’ };
plotting . centerline . v e l o c i t i e s . file_end = { ’ velocity ’ };
p l o t t i n g . c e n t e r l i n e . v e l o c i t i e s . y_units = ’ [m/ s ] ’ ;
%p l o t t i n g . c e n t e r l i n e . v e l o c i t i e s . a s p e c t _ r a t i o = [ 5 1 1 ] ;
p l o t t i n g . c e n t e r l i n e . v e l o c i t i e s . l i m i t s . mean . a x i a l = [ 0 0 . 0 5 −15
3 5 ] ; % [ 0 0 . 0 5 −15 3 5 ] % % % % % %
%% S t a t i o n s
p l o t t i n g . s t a t i o n s . type = { ’ v e l o c i t i e s ’ } ; % ’ v e l o c i t i e s ’ , ’
species ’
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
117
Apéndices
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . s t a t i s t i c = { ’ mean ’ , ’RMS ’ } ; % ’
i n s t ’ , ’ mean ’ , ’RMS’
plotting . stations . velocities . directions = { ’ axial ’ , ’ radial ’ , ’
tangential ’ }; % ’ axial ’ , ’ radial ’ , ’ tangential ’
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . s t a t i o n s = [ 5 10 20 30 4 0 ] ; % [ 5
15 29 46 76 9 2 ] , [ 5 10 20 30 4 0 ]
plotting . stations . v e l o c i t i e s . file_end = { ’ velocity ’ };
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . y_units = ’ [m/ s ] ’ ;
%p l o t t i n g . s t a t i o n s . v e l o c i t i e s . a s p e c t _ r a t i o = [ 5 1 1 ] ;
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . l i m i t s . x l i m i t = [ −0.05 0 . 0 5 ] ;
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . l i m i t s . mean . a x i a l = [ −10 2 5 ] ;
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . l i m i t s . mean . r a d i a l = [−5 1 5 ] ;
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . l i m i t s . mean . t a n g e n t i a l = [ −25 5 ] ;
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . l i m i t s .RMS. a x i a l = [ 0 1 6 ] ;
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . l i m i t s .RMS. r a d i a l = [ 0 1 0 ] ;
p l o t t i n g . s t a t i o n s . v e l o c i t i e s . l i m i t s .RMS. t a n g e n t i a l = [ 0 1 0 ] ;
p l o t t i n g . s t a t i o n s . s p e c i e s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
p l o t t i n g . s t a t i o n s . s p e c i e s . s t a t i s t i c = { ’ mean ’ , ’RMS ’ } ; % ’ mean
’ , ’RMS’
p l o t t i n g . s t a t i o n s . s p e c i e s . d i r e c t i o n s = { ’CH4 ’ } ; %
p l o t t i n g . s t a t i o n s . s p e c i e s . s t a t i o n s = [ 6 10 14 2 2 ] ; % [ 6 10 14
22]
p l o t t i n g . s t a t i o n s . s p e c i e s . f i l e _ e n d = { ’ mass ’ , ’ f r a c t i o n ’ } ;
p l o t t i n g . s t a t i o n s . s p e c i e s . y_units = ’ [ −] ’ ;
%p l o t t i n g . s t a t i o n s . s p e c i e s . a s p e c t _ r a t i o = [ 5 1 1 ] ;
p l o t t i n g . s t a t i o n s . s p e c i e s . l i m i t s . x l i m i t = [ −0.03 0 . 0 3 ] ;
p l o t t i n g . s t a t i o n s . s p e c i e s . l i m i t s . mean . CH4 = [ 0 0 . 1 5 ] ;
p l o t t i n g . s t a t i o n s . s p e c i e s . l i m i t s .RMS. CH4 = [ 0 0 . 2 0 ] ;
%% Contours
p l o t t i n g . c o n t o u r s . type = { ’ v e l o c i t i e s ’ , ’ t k e ’ , ’ s t r e s s ’ } ; % ’
v e l o c i t i e s ’ , ’ tke ’ , ’ s t r e s s ’ , ’ yplus ’
p l o t t i n g . c o n t o u r s . v e l o c i t i e s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
p l o t t i n g . c o n t o u r s . v e l o c i t i e s . s t a t i s t i c = { ’ mean ’ , ’RMS ’ } ; % ’
mean ’ , ’RMS’
plotting . contours . v e l o c i t i e s . d i r e c t i o n s = { ’ axial ’ , ’ r a d i a l ’ , ’
tangential ’ }; % ’ axial ’ , ’ radial ’ , ’ tangential ’
plotting . contours . v e l o c i t i e s . file_end = ’ v e l o c i t y ’ ;
118
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
p l o t t i n g . c o n t o u r s . t k e . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
p l o t t i n g . contours . tke . file_end = ’ tke ’ ;
p l o t t i n g . c o n t o u r s . t k e . y_units = ’ [ −] ’ ;
p l o t t i n g . c o n t o u r s . t k e . t i t l e _ e n d = ’ Turbulent ␣ K i n e t i c ␣ Energy ’ ;
p l o t t i n g . c o n t o u r s . t k e . l i m i t s = [ −0.05 0 . 0 5 0 0 . 1 ] ; % [ −0.05
0.05 0 0 . 1 ]
p l o t t i n g . contours . tke . aspect_ratio = [1 1 1 ] ;
p l o t t i n g . c o n t o u r s . s t r e s s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
plotting . contours . s t r e s s . file_end = ’ s t r e s s ’ ;
p l o t t i n g . c o n t o u r s . s t r e s s . y_units = ’ [ −] ’ ;
p l o t t i n g . c o n t o u r s . s t r e s s . t i t l e _ e n d = ’ Reynolds ␣ S t r e s s ’ ;
%p l o t t i n g . c o n t o u r s . s t r e s s . l i m i t s = [ 0 0 . 1 ] ; % [ −0.05 0 . 0 5 0
0.1]
plotting . contours . s t r e s s . aspect_ratio = [1 1 1 ] ;
%p l o t t i n g . c o n t o u r s . y p l u s . p l a n e = { ’ xz ’ } ; % ’ xz ’ , ’ yz ’
%p l o t t i n g . c o n t o u r s . y p l u s . l i m i t s = [ −0.05 0 . 0 5 0 0 . 1 ] ; % [ −0.05
0.05 0 0 . 1 ]
%p l o t t i n g . c o n t o u r s . y p l u s .
end
B.4. PLOTCONFIG
codigo/plotConfig.m
% 14/05/2018
% C a r l o s Moreno Montagud
function c o n f i g = p l o t C o n f i g ( )
%% Line s t y l e
c o n f i g . sim_markers = [ ’ s ’ ’ ^ ’ ’ x ’ ’ d ’ ’ v ’ ’ . ’ ’> ’ ’ ∗ ’ ’< ’ ] ; %
Markers f o r s i m u l a t i o n s
119
Apéndices
c o n f i g . sim_marker_size = 6 ;
c o n f i g . sim_line_width = 1 . 5 ;
c o n f i g . l i n e _ s t y l e = [ ’ ␣− ’ ; ’−− ’ ; ’ −. ’ ; ’ ␣ : ’ ] ; % Line o r d e r
c o n f i g . c o l o r s = [ ’ b ’ ’ g ’ ’ r ’ ’ c ’ ’m’ ’ y ’ ’w ’ ] ; % Color o r d e r
c o n f i g . primes = [ 5 , 7 , 1 1 , 1 3 , 1 7 ] ; % S p a c i n g markers
%% Font s i z e s
config . li nes . legend_fontsize = 10;
%c o n f i g . l i n e s . t i t l e _ f o n t s i z e = 14;
config . lines . xlabel_fontsize = 14;
config . lines . ylabel_fontsize = 14;
%c o n f i g . c o n t o u r s . l e g e n d _ f o n t s i z e = 10;
config . contours . t i t l e _ f o n t s i z e = 12;
config . contours . xlabel1_fontsize = 10;
config . contours . ylabel1_fontsize = 10;
config . contours . xlabel2_fontsize = 14;
config . contours . ylabel2_fontsize = 14;
%% Contour
c o n f i g . cmap = j e t ; % j e t , g r a y
end
B.5. PATH_NAMES
codigo/functions/path_names.m
120
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
% 04/04/2018
% C a r l o s Moreno Montagud
%% Cases d a t a
l g e o = length ( c a s e s . geometry ) ;
l s o l = length ( c a s e s . solveMode ) ;
l i n j = length ( c a s e s . i n j e c t i o n ) ;
l f l u = length ( c a s e s . f l u i d ) ;
l i n l = length ( c a s e s . i n l e t T y p e ) ;
lnum = length ( c a s e s . numCase ) ;
l m i l = length ( c a s e s . m i l i S e c ) ;
l c o d = length ( c a s e s . code ) ;
lmes = length ( c a s e s . mesh) ;
i n j e c t i o n A b b r = c e l l f u n ( @abbr , i n j e c t i o n , ’ un ’ , 0 ) ;
inletTypeAbbr = c e l l f u n ( @abbr , i n l e t T y p e , ’ un ’ , 0 ) ;
f l u i d A b b r = c e l l f u n (@( x ) upper ( x ( 1 ) ) , f l u i d , ’ un ’ , 0 ) ;
codeAbbr = c e l l f u n ( @abbr , code , ’ un ’ , 0 ) ;
%% Names
% caseName NASA_RANS_NP_G_NR_X / JunCai
% matName NR_NP_G_NR_X_XXX_CG_(6M) / JunCai
g e o S o l v e = m a t 2 c e l l ( c e l l f u n (@( x ) x ( 1 ) , [ geometry solveMode ] ) ,
ones ( 1 , casesNum ) , 2 ) ;
121
Apéndices
i f strcmp ( code , ’ L i t e r a t u r e ’ )
f i l e = c a s e s . author ( 1 ) ;
names . caseName = f i l e ;
names . matName = f i l e ;
else
f o l d e r C e l l = m a t 2 c e l l ( [ geometry solveMode i n j e c t i o n A b b r
f l u i d A b b r inletTypeAbbr numCase ] , ones ( 1 , casesNum ) , 6 ) ;
f o l d e r = c e l l f u n (@( x ) s t r j o i n ( x , ’_ ’ ) , f o l d e r C e l l , ’ un ’ , 0 ) ;
names . caseName = f o l d e r ;
f i l e = strrep ( s t r c a t ( c e l l f u n (@( x ) s t r j o i n ( x , ’_ ’ ) , f i l e C e l l ,
’ un ’ , 0 ) , ’_ ’ ) , ’__ ’ , ’_ ’ ) ;
f i l e = c e l l f u n (@( x ) x ( 1 : end−1) , f i l e , ’ un ’ , 0 ) ;
names . matName = f i l e ;
end
% expName JunCai
names . expName = c a s e s . author ;
% legend_sim NR_NP_G_NR_X_XXX_CG_6M
names . legend_sim = strrep ( c e l l f u n (@( x ) s t r j o i n ( x ( l o g i c a l (
l e g e n d _ f l a g s ) ) , ’_ ’ ) , f i l e C e l l , ’ un ’ , 0 ) , ’__ ’ , ’_ ’ ) ;
% l e g en d _exp JunCai % % % % % % % % % % % % % %
names . legend_exp = c a s e s . author ;
122
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
i n j = unique ( i n j e c t i o n ) ; % nonPremixed
i f length ( i n j )==1
results_folder = strcat ( results_folder , inj , f i l e s e p ) ;
end
f l u = unique ( f l u i d ) ; % g a s
123
Apéndices
i f length ( f l u )==1
results_folder = strcat ( results_folder , flu , f i l e s e p ) ;
end
i n l = unique ( i n l e t T y p e ) ; % n o n R e a c t i v e
i f length ( i n l )==1
results_folder = strcat ( results_folder , inl , f i l e s e p ) ;
end
% m i l = u n i q u e ( m i l i S e c ) ; % XXX
% i f l e n g t h ( m i l )==1
% r e s u l t s F o l d e r = s t r c a t ( r e s u l t s F o l d e r , mil , ’ ms ’ , f i l e s e p ) ;
% end
names . r e s u l t s _ f o l d e r = r e s u l t s _ f o l d e r { 1 } ;
names . code = code ;
names . m i l i S e c = m i l i S e c ;
end
function s t r = abbr ( f s t r )
switch f s t r
124
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
c a s e ’CONVERGE’
s t r = ’CG ’ ;
c a s e ’OpenFOAM ’
s t r = ’OF ’ ;
c a s e ’ nonPremixed ’
s t r = ’NP ’ ;
c a s e ’ premixed ’
s t r = ’P ’ ;
c a s e ’ nonReactive ’
s t r = ’NR ’ ;
case ’ reactive ’
s t r = ’R ’ ;
otherwise
str = ’ error ’ ;
end
end
B.6. UPPERFIRST
codigo/functions/upperFirst.m
function s t r = u p p e r F i r s t ( data )
s t r = [ upper ( data ( 1 ) ) data ( 2 : end ) ] ;
end
B.7. LOADCASES
codigo/functions/import/loadCases.m
% 07/05/2018
% C a r l o s Moreno Montagud
function [EXP, SIM]= l o a d C a s e s ( names , import , f l a g s )
% Load SIM
for i =1: names . casesNum
125
Apéndices
i f e x i s t ( names . m a t F i l e { i } , ’ f i l e ’ )
% Load f i l e . mat
load ( names . m a t F i l e { i , 1 } , names . matName{ i , 1 } ) ;
SIM{ i , 1 } = eval ( names . matName{ i , 1 } ) ;
else
% Import c a s e and s a v e f i l e . mat
SIM{ i , 1 } = importCase ( names , import , i ) ;
eval ( s t r c a t ( names . matName{ i } , ’=SIM{ i , 1 } ; ’ ) ) ;
save ( names . m a t F i l e { i } , names . matName{ i }) ;
end
end
% Load EXP
i f f l a g s . p l o t _ f l a g | | f l a g s . compare_flag
for i =1:expNum
load ( names . e x p F i l e { i , 1 } , names . expName{ i , 1 } ) ;
EXP{ i , 1 } = eval ( names . expName{ i , 1 } ) ;
end
end
end
B.8. IMPORTCASE
codigo/functions/import/importCase.m
% 04/04/2018
% C a r l o s Moreno Montagud
c o l _ f i l e s = import . c o l _ f i l e s ;
o u t _ f i l e s = import . o u t _ f i l e s ;
%% Import
126
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
s w i t c h code
c a s e ’CONVERGE’
% Search f i l e s . c o l
c o l F o l d e r = s t r c a t ( im por tF ol der , ’ output ’ , f i l e s e p ) ;
m i l i S e c S t r = num2str ( s t r 2 d o u b l e ( m i l i S e c ) / 10 0 0 , ’ %1.5 e ’ )
;
m i l i S e c R e a d = 5− s t r 2 d o u b l e ( m i l i S e c S t r ( end ) ) ;
r e a d C o l F i l e S t r = s t r c a t ( c o l F o l d e r , c o l _ f i l e s , ’ ∗_+ ’ ,
m i l i S e c S t r ( 1 : m i l i S e c R e a d ) , ’ ∗ ’ , m i l i S e c S t r ( end−3:end )
, ’ . col ’ ) ;
r e a d C o l F i l e = c e l l f u n ( @dir , r e a d C o l F i l e S t r , ’ un ’ , 0 ) ;
readColFileMat = cell2mat ( readColFile ) ;
l a s t R e a d C o l F i l e = { r e a d C o l F i l e M a t ( end , : ) . name } ;
readColFilename = s t r c a t ( c o l F o l d e r , l a s t R e a d C o l F i l e ) ;
% Search f i l e s . o u t
r e a d O u t F i l e S t r = s t r c a t ( i mpo rt Fo lde r , o u t _ f i l e s , ’ ∗_ ’ , ’
point_ ’ , { ’ 1 ’ , ’ 2 ’ } , ’_ ’ , ’ volume ’ , ’ _avg . out ’ ) ;
r e a d O u t F i l e = c e l l f u n ( @dir , r e a d O u t F i l e S t r , ’ un ’ , 0 ) ;
readOutFileMat = c e l l 2 m a t ( r e a d O u t F i l e ) ;
a l l R e a d O u t F i l e = reshape ({ readOutFileMat . name } , s i z e (
readOutFileMat ) ) ;
readOutFilename = s t r c a t ( i mpo rt Fo ld er , a l l R e a d O u t F i l e ) ;
% siempre 1 columna
% Read v a r i a b l e s
delimiter = ’ ’ ;
startRow = 6 ;
% Import CONVERGE f i l e s
dataArray = importCVG ( readColFilename , readOutFilename ,
import , d e l i m i t e r , startRow ) ;
c a s e ’OpenFOAM ’
% Search 0 and l a s t ’ m i l i S e c ’ f o l d e r
readFile_temp = dir ( i m p o r t F o l d e r ) ;
readFile_temp = readFile_temp (~ isnan ( s t r 2 d o u b l e ({
readFile_temp . name }) ) ) ;
i f ~isempty ( m i l i S e c )
127
Apéndices
r e a d F i l e _ f i n d = s t r f i n d ({ readFile_temp . name } , [ ’ . ’ ,
miliSec ] ) ;
r e a d F i l e _ f i n d = c e l l f u n (@( x ) ~isempty ( x ) ,
readFile_find ) ;
readFile_find (1) = 1;
r e a d F i l e = readFile_temp ( r e a d F i l e _ f i n d ) ;
else
r e a d F i l e = readFile_temp ( [ 1 , end ] ) ;
end
l a s t R e a d F i l e = { r e a d F i l e . name } ;
r e a d F o l d e r s = s t r c a t ( i mp or tF old er , l a s t R e a d F i l e , f i l e s e p
);
% Read v a r i a b l e s
delimiter = ’␣ ’ ;
startRow = 2 3 ;
% Import OpenFOAM f i l e s
dataArray = importFOAM( r e a d F o l d e r s , import , d e l i m i t e r ,
startRow ) ; %% % % % % % % % % % % % %
case ’ Literature ’
% Read v a r i a b l e s
d e l i m i t e r = ’ ; ’ ; %% % % % % % % % % %e%
r r%o %
r %%%
startRow = 2 ;
% Import L i t e r a t u r e f i l e s
dataArray = importEXP ( i mpo rt Fo ld er , import , d e l i m i t e r ,
startRow ) ;
end
end
B.9. IMPORTCVG
codigo/functions/import/importCVG.m
% 01/02/2018
% Mario Belmar G i l
% 21/02/2018
128
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
% C a r l o s Moreno Montagud
b a s e _ s i z e = import . b a s e _ s i z e ;
AMR = import .AMR;
%% Col f i l e s
c o l _ f i l e s = s t r c a t ( import . c o l _ f i l e s , ’ _ f i l e ’ ) ;
n F i l e s = length ( c o l _ f i l e s ) ;
for f f =1: n F i l e s % ’ v e l o c i t i e s ’ , ’ s t r e s s ’ , ’ s p e c i e s ’
current_file = col_files { f f };
F=[];
% Import raw d a t a
CFD_raw = importdata ( readColFilename { f f }) ;
v a r i a b l e s _ n a m e s = a r r a y f u n (@( s t r ) strrep ( s t r , ’ . ’ , ’_ ’ ) ,
CFD_raw . c o l h e a d e r s ) ;
CFD_data = a r r a y 2 t a b l e (CFD_raw . data , ’ VariableNames ’ ,
variables_names ) ;
v a r i a b l e s = import . ( c u r r e n t _ f i l e ) . v a r i a b l e s ;
nVar = length ( v a r i a b l e s ) ;
for vv =1:nVar % ’ v e l o c i t i e s ’ , ’ t k e ’ , ’ s t r e s s ’ , ’ s p e c i e s ’
c u r r e n t _ v a r = v a r i a b l e s {vv } ;
129
Apéndices
switch current_plot
% CENTERLINE
case ’ c e n t e r l i n e ’
% Coordinates
coord_min = s t r u c t _ p l o t . l i m i t s ( 1 ) ;
coord_max = s t r u c t _ p l o t . l i m i t s ( 2 ) ;
d i s c r e t i z a t i o n = round ( ( coord_max−
coord_min ) ∗2^AMR/ b a s e _ s i z e ) ;
c o o r d s = linspace ( coord_min , coord_max
, discretization ) ’;
z e r o _ v e c t o r = zeros ( d i s c r e t i z a t i o n , 1 )
;
dataArray . c e n t e r l i n e . ( c u r r e n t _ p l a n e ) . (
current_var ) . coords = coords ;
s t a t i s t i c = struct_plot . s t a t i s t i c ;
n S t a t i s t i c = length ( s t a t i s t i c ) ;
for s s =1: n S t a t i s t i c % ’ mean ’
current_statistic = statistic {
ss };
%% % % % % p l%
a n e z a x i a l = W, e t c
head = [ h e a d e r s (
c u r r e n t _ s t a t i s t i c ) headers (
current_direction ) ] ;
try
% Interpolation
130
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
i f ~ i s f i e l d (F , head )
F . ( head ) =
scatteredInterpolant
(CFD_data . x ,
CFD_data . y ,
CFD_data . z ,
CFD_data . ( head ) ) ;
end
% Discretization
dataArray . c e n t e r l i n e . (
current_plane ) . (
current_var ) . (
current_statistic ) .(
current_direction ) .
c e n t e r l i n e = F . ( head ) (
zero_vector ,
zero_vector , coords ) ;
catch
disp ( ’ e r r o r ␣ cvg ␣ c e n t e r l i n e
’);
end
end
end
% STATIONS
case ’ stations ’
% Coordinates
coord_min = import . ( [ c u r r e n t _ p l a n e ( 1 ) ,
’ lims ’ ] ) (1) ;
coord_max = import . ( [ c u r r e n t _ p l a n e ( 1 ) ,
’ lims ’ ] ) (2) ;
d i s c r e t i z a t i o n = round ( ( coord_max−
coord_min ) ∗2^AMR/ b a s e _ s i z e ) ;
c o o r d s = linspace ( coord_min , coord_max
, discretization ) ’;
z e r o _ v e c t o r = zeros ( d i s c r e t i z a t i o n , 1 )
;
switch current_plane
c a s e ’ xz ’
131
Apéndices
c1 = c o o r d s ;
c2 = z e r o _ v e c t o r ;
c a s e ’ yz ’
c1 = z e r o _ v e c t o r ;
c2 = c o o r d s ;
end
dataArray . s t a t i o n s . ( c u r r e n t _ p l a n e ) . (
current_var ) . coords = coords ;
s t a t i s t i c = struct_plot . s t a t i s t i c ;
n S t a t i s t i c = length ( s t a t i s t i c ) ;
for s s =1: n S t a t i s t i c % ’ mean ’ , ’RMS
’
current_statistic = statistic {
ss };
stations = struct_plot .
stations ;
n S t a t i o n s = length ( s t a t i o n s ) ;
for mm=1: n S t a t i o n s % [ 5 10 20
30 4 0 ]
c u r r e n t _ s t a t i o n = num2str (
s t a t i o n s (mm) ) ;
s t a t i o n _ v e c t o r = repmat (
s t a t i o n s (mm) / 1 0 0 0 , [
discretization 1]) ;
head = [ h e a d e r s (
current_statistic )
headers (
current_direction ) ] ;
try
% Interpolation
i f ~ i s f i e l d (F , head )
132
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
F . ( head ) =
scatteredInterpolant
(CFD_data . x ,
CFD_data . y ,
CFD_data . z ,
CFD_data . ( head )
);
end
% Discretization
dataArray . s t a t i o n s . (
current_plane ) . (
current_var ) . (
current_statistic )
. ( current_direction
) . ( [ ’ station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm
’ ] ) = F . ( head ) ( c1 ,
c2 , s t a t i o n _ v e c t o r )
;
catch
disp ( ’ e r r o r ␣ cvg ␣
stations ’ ) ;
end
end
end
end
% CONTOURS
case ’ contours ’
% Coordinates
l i m i t s = [ import . ( [ c u r r e n t _ p l a n e ( 1 ) , ’
lims ’ ] ) struct_plot . l i m i t s ] ;
d i s c r e t i z a t i o n = round ( ( l i m i t s ( 2 )−
l i m i t s ( 1 ) ) ∗2^AMR/ b a s e _ s i z e ) ; % xy
dependant
v e r t e x 1 = linspace ( l i m i t s ( 1 ) , l i m i t s
( 2 ) , d i s c r e t i z a t i o n +1) ’ ; c e n t e r 1 =
v e r t e x 1 ( 1 : end−1)+d i f f ( v e r t e x 1 ) / 2 ;
133
Apéndices
v e r t e x 2 = linspace ( l i m i t s ( 3 ) , l i m i t s
( 4 ) , d i s c r e t i z a t i o n +1) ’ ; c e n t e r 2 =
v e r t e x 2 ( 1 : end−1)+d i f f ( v e r t e x 2 ) / 2 ;
coords = [ center1 , center2 ] ;
coord0 = zeros ( d i s c r e t i z a t i o n ) ;
[ coord1 , coord2 ] = meshgrid ( c e n t e r 1 ,
center2 ) ;
switch current_plane
c a s e ’ xz ’
c1 = coord1 ;
c2 = coord0 ;
c3 = coord2 ;
c a s e ’ yz ’
c1 = coord0 ;
c2 = coord1 ;
c3 = coord2 ;
end
dataArray . c o n t o u r s . ( c u r r e n t _ p l a n e ) . (
current_var ) . coords = coords ;
% Velocities
case ’ v e l o c i t i e s ’
for dd=1: n D i r e c t i o n % ’ a x i a l ’ ,
’ radial ’ , ’ tangential ’
current_direction =
d i r e c t i o n {dd } ;
s t a t i s t i c = struct_plot .
statistic ;
n S t a t i s t i c = length (
statistic ) ;
for s s =1: n S t a t i s t i c % ’
mean ’ , ’RMS’
current_statistic =
s t a t i s t i c { ss };
134
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
head = [ h e a d e r s (
current_statistic )
headers (
current_direction )
];
try
% Interpolation
i f ~ i s f i e l d (F , head
)
F . ( head ) =
scatteredInterpolant
(CFD_data . x
, CFD_data .
y , CFD_data
.z,
CFD_data . (
head ) ) ;
end
% Discretization
dataArray . c o n t o u r s
. ( current_plane
) . ( current_var )
.(
current_statistic
) .(
current_direction
) . contour = F . (
head ) ( c1 , c2 ,
c3 ) ;
catch
disp ( ’ e r r o r ␣ cvg ␣
contour ␣
velocities ’) ;
end
end
end
% Tke
case ’ tke ’
135
Apéndices
r e s u l t = zeros ( d i s c r e t i z a t i o n ,
discretization , nDirection ) ;
current_statistic =
struct_plot . s t a t i s t i c {1}; %
’RMS’
for dd=1: n D i r e c t i o n % ’ a x i a l ’ ,
’ radial ’ , ’ tangential ’
current_direction =
d i r e c t i o n {dd } ;
head = [ h e a d e r s (
current_statistic )
headers (
current_direction ) ] ;
try
% Interpolation
i f ~ i s f i e l d (F , head )
F . ( head ) =
scatteredInterpolant
(CFD_data . x ,
CFD_data . y ,
CFD_data . z ,
CFD_data . ( head )
);
end
% Discretization
v e l = F . ( head ) ( c1 , c2 ,
c3 ) ;
% Calculate
r e s u l t ( : , : , dd ) = v e l . ∗
vel ;
catch
disp ( ’ e r r o r ␣ cvg ␣
contour ␣ tke ’ ) ;
end
end
136
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
dataArray . c o n t o u r s . (
current_plane ) . ( current_var
) . contour = sum( r e s u l t , 3 )
/2;
% Stress
case ’ s t r e s s ’
t o t a l _ s t r e s s = zeros ( s i z e (
CFD_data . x ) ) ;
% Calculate
try
s t r e s s = CFD_data . (
head ) ;
total_stress =
t o t a l _ s t r e s s+s t r e s s
.^2;
% Symmetry
i f rem( s t r 2 d o u b l e (
current_direction
( 2 : end ) ) , 1 1 ) ~=0
total_stress =
t o t a l _ s t r e s s+
stress .^2;
end
catch
end
end
% Interpolation
%i f ~ i s f i e l d (F , head )
F. stress =
scatteredInterpolant (
CFD_data . x , CFD_data . y ,
137
Apéndices
CFD_data . z , sqrt (
total_stress ) ) ;
%end
% Discretization
dataArray . c o n t o u r s . (
current_plane ) . ( current_var
) . contour = F . s t r e s s ( c1 , c2
, c3 ) ;
% Otherwise
end
end
end
end
end
end
%% Out f i l e s
i f ~isempty ( readOutFilename )
o u t _ f i l e s = s t r c a t ( import . o u t _ f i l e s {1} , ’ _ f i l e ’ ) ; % monitor
only
f s = import . ( o u t _ f i l e s ) . f s ;
F = [];
o u t _ s u b f i l e s = readOutFilename ;
[ nSubfiles , nPoints ] = size ( o u t _ s u b f i l e s ) ;
for oo =1: n P o i n t s % ’ point_1 ’ , . . .
dataRead = [ ] ;
for f f =1: n S u b f i l e s % ’ monitor38 ’ , . . .
c u r r e n t _ s u b f i l e = o u t _ s u b f i l e s { f f , oo } ;
% Import
numCols = 4 ;
dataRead = [ dataRead ; c e l l 2 m a t ( i m p o r t T e x t F i l e (
c u r r e n t _ s u b f i l e , d e l i m i t e r , startRow , numCols , ’ %17
f ’ , false )) ];
%’ MultipleDelimsAsOne ’ , 1
end
138
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
F . ( c u r r e n t _ v a r ) . ( [ ’ point_ ’ , num2str ( oo ) ] ) =
dataRead ( : , vv+1) ;
i f oo==n P o i n t s
dataArray . psd . ( c u r r e n t _ v a r ) . cross_spectrum =
angle ( cpsd (F . ( c u r r e n t _ v a r ) . point_1 , F . (
c u r r e n t _ v a r ) . point_2 , [ ] , [ ] , [ ] , f s ) ) ;
%dataArray . psd . ( c u r r e n t _ v a r ) . msqrd_coherence =
mscohere (F . ( c u r r e n t _ v a r ) . point_1 , F . (
c u r r e n t _ v a r ) . point_2 , [ ] , [ ] , [ ] , f s ) ;
end
end
end
end
end
%% Secondary f u n c t i o n s
function s t r = h e a d e r s ( s t a t )
switch stat
% case ’ inst ’
% str = [ ] ;
c a s e ’ mean ’
s t r = ’BAR_ ’ ;
c a s e ’RMS ’
s t r = ’RMS_ ’ ;
c a s e ’ a x i a l ’ %% % % % % % % % % %yz % %p%
l a%
n e%
s t r = ’W’ ;
c a s e ’ r a d i a l ’ %% % % % % % % % % % % % % %
s t r = ’U ’ ;
c a s e ’ t a n g e n t i a l ’ %% % % % % % % % % % %
s t r = ’V ’ ;
otherwise
139
Apéndices
s t r = [ ’Y_ ’ s t a t ] ; %% % % % % % % % %
end
end
B.10. IMPORTEXP
codigo/functions/import/importEXP.m
% 21/02/2018
% C a r l o s Moreno Montagud
% Filename s t r u c t u r e :
% plane , s t a t i s t i c , d i r e c t i o n , f i l e _ e n d , p l o t _ t y p e
% A l l t o g e t h e r w i t h o u t s p a c e and upper c a s e f i r s t l e t t e r o f
the variable
% Example : x z M e a n A x i a l V e l o c i t y C e n t e r l i n e
%
% EXP f i l e s are s e p a r a t e d by ’ ; ’ d e l i m i t e r , f i r s t row are
t i t l e s and t h e
% number o f columns i s d o u b l e o f s t a t i o n s number (2 f o r
centerline ) .
%% Loops
f i l e s = import . c o l _ f i l e s ;
n F i l e s = length ( f i l e s ) ;
for f f =1: n F i l e s % ’ v e l o c i t i e s ’ , ’ s t r e s s ’ , ’ s p e c i e s ’
current_file = [ f i l e s { ff } , ’ _file ’ ] ;
v a r i a b l e s = import . ( c u r r e n t _ f i l e ) . v a r i a b l e s ;
nVar = length ( v a r i a b l e s ) ;
for vv =1:nVar % ’ v e l o c i t i e s ’ , ’ t k e ’ , ’ s t r e s s ’ , ’ s p e c i e s ’
c u r r e n t _ v a r = v a r i a b l e s {vv } ;
140
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
f i l e _ e n d = import . ( c u r r e n t _ v a r ) . f i l e _ e n d ;
s t r u c t _ p l o t = import . ( c u r r e n t _ v a r ) . ( c u r r e n t _ p l o t ) ;
% e r r o r , don ’ t import s t r e s s w i t h
L i t e r a t u r e ( s e e i m p o r t P r o p e r t i e s .m)
s t a t i s t i c = struct_plot . s t a t i s t i c ;
n S t a t i s t i c = length ( s t a t i s t i c ) ;
for s s =1: n S t a t i s t i c % ’ mean ’ , ’RMS’
current_statistic = s t a t i s t i c { ss };
f i l e C e l l = [{ current_statistic ,
current_direction } , file_end ,{
current_plot } ] ;
f i l e J o i n = s t r j o i n ( c e l l f u n ( @upperFirst
, f i l e C e l l , ’ un ’ , 0 ) , ’ ’ ) ;
readFilename = [ r e a d F o l d e r ,
current_plane , f i l e J o i n , ’ . csv ’ ] ;
switch current_plot
% CENTERLINE
case ’ c e n t e r l i n e ’
try
impData = i m p o r t T e x t F i l e (
readFilename , d e l i m i t e r ,
141
Apéndices
startRow , 2 , ’ %f ’ , f a l s e ) ;
catch
disp ( [ ’ e r r o r ␣ import ␣ exp ␣ ’ ,
current_var , ’ ␣
centerline ’ ]) ;
end
dataArray . c e n t e r l i n e . (
current_plane ) . ( current_var
) .( current_statistic ) .(
current_direction ) . (
c u r r e n t _ p l o t )=c e l l 2 m a t (
impData ) ;
% STATIONS
case ’ stations ’
stations = struct_plot .
stations ;
n S t a t i o n s = length ( s t a t i o n s ) ;
try
impData = i m p o r t T e x t F i l e (
readFilename , d e l i m i t e r ,
startRow , 2 ∗ n S t a t i o n s , ’ %
f ’ , false ) ;
catch
disp ( [ ’ e r r o r ␣ import ␣ exp ␣ ’ ,
current_var , ’ ␣ s t a t i o n s ’
]) ;
end
for mm=1: n S t a t i o n s % [ 5 10 20
30 4 0 ]
c u r r e n t _ s t a t i o n = num2str (
s t a t i o n s (mm) ) ;
dataArray . s t a t i o n s . (
current_plane ) . (
current_var ) . (
current_statistic ) .(
current_direction ) . ( [ ’
station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] )
= c e l l 2 m a t ( impData ( ( 2 ∗
142
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
mm−1) : ( 2 ∗mm) ) ) ;
end
end
end
end
end
end
end
end
end
B.11. IMPORTFOAM
codigo/functions/import/importFOAM.m
% 04/04/2018
% C a r l o s Moreno Montagud
b a s e _ s i z e = import . b a s e _ s i z e ;
AMR = import .AMR;
c o o r d _ f i l e s = import . c r d _ f i l e s ;
v a r _ f i l e s = import . v a r _ f i l e s ;
%% F o l d e r 0 ( c o o r d i n a t e s )
r e ad F i lename = s t r c a t ( r e a d F o l d e r s {1} , c o o r d _ f i l e s ) ;
n F i l e s = length ( c o o r d _ f i l e s ) ;
for f f =1: n F i l e s
current_file = coord_files { f f };
c u r r e n t _ f i l e n a m e = readFilename { f f } ;
numCols = 1 ;
CFD_data . ( c u r r e n t _ f i l e ) = c e l l 2 m a t ( i m p o r t T e x t F i l e (
c u r r e n t _ f i l e n a m e , d e l i m i t e r , startRow , numCols , ’ %f ’ , f a l s e )
);
143
Apéndices
end
%% F o l d e r end
r e ad F i l ename = s t r c a t ( r e a d F o l d e r s {2} , v a r _ f i l e s ) ;
n F i l e s = length ( v a r _ f i l e s ) ;
for f f =1: n F i l e s
current_file = var_files { f f };
c u r r e n t _ f i l e n a m e = readFilename { f f } ;
F=[];
i f strcmp ( c u r r e n t _ f i l e , ’ UPrime2Mean ’ )
CFD_data . ( c u r r e n t _ f i l e ) = r e a l ( sqrt (CFD_data . (
current_file ) ) ) ;
end
v a r i a b l e s = import . ( c u r r e n t _ f i l e ) . v a r i a b l e s ;
nVar = length ( v a r i a b l e s ) ;
for vv =1:nVar % ’ v e l o c i t i e s ’ , ’ t k e ’
c u r r e n t _ v a r = v a r i a b l e s {vv } ;
144
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
for dd=1: n D i r e c t i o n % ’ a x i a l ’ , ’ r a d i a l ’ , ’
tangential ’
c u r r e n t _ d i r e c t i o n = d i r e c t i o n {dd } ;
% mean , RMS
head = [ h e a d e r s ( c u r r e n t _ f i l e ) ’_ ’ h e a d e r s (
current_direction ) ] ;
c u r r e n t _ c o l = impCols ( dd ) ;
switch current_plot
% CENTERLINE
case ’ c e n t e r l i n e ’
coord_min = s t r u c t _ p l o t . l i m i t s ( 1 ) ;
coord_max = s t r u c t _ p l o t . l i m i t s ( 2 ) ;
d i s c r e t i z a t i o n = round ( ( coord_max−
coord_min ) ∗2^AMR/ b a s e _ s i z e ) ;
c o o r d s = linspace ( coord_min ,
coord_max , d i s c r e t i z a t i o n ) ’ ;
z e r o _ v e c t o r = zeros ( d i s c r e t i z a t i o n
, 1) ;
dataArray . c e n t e r l i n e . (
current_plane ) . ( current_var ) .
coords = coords ;
try
% Interpolation
i f ~ i s f i e l d (F , head )
F . ( head ) =
scatteredInterpolant (
CFD_data . ccx , CFD_data .
ccy , CFD_data . ccz ,
CFD_data . ( c u r r e n t _ f i l e )
( : , current_col ) ) ;
end
% Discretization
dataArray . c e n t e r l i n e . (
current_plane ) . ( current_var
) . ( headers ( c u r r e n t _ f i l e ) ) . (
current_direction ) .
c e n t e r l i n e = F . ( head ) (
145
Apéndices
zero_vector , zero_vector ,
coords ) ;
catch
disp ( ’ e r r o r ␣foam ␣ c e n t e r l i n e ’ ) ;
end
% STATIONS
case ’ stations ’
% Coordinates
coord_min = import . ( [ c u r r e n t _ p l a n e
(1) , ’ lims ’ ] ) (1) ;
coord_max = import . ( [ c u r r e n t _ p l a n e
(1) , ’ lims ’ ] ) (2) ;
d i s c r e t i z a t i o n = round ( ( coord_max−
coord_min ) ∗2^AMR/ b a s e _ s i z e ) ;
c o o r d s = linspace ( coord_min ,
coord_max , d i s c r e t i z a t i o n ) ’ ;
z e r o _ v e c t o r = zeros ( d i s c r e t i z a t i o n
, 1) ;
switch current_plane
c a s e ’ xz ’
c1 = c o o r d s ;
c2 = z e r o _ v e c t o r ;
c a s e ’ yz ’
c1 = z e r o _ v e c t o r ;
c2 = c o o r d s ;
end
dataArray . s t a t i o n s . ( c u r r e n t _ p l a n e )
. ( current_var ) . coords = coords ;
146
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
try
% Interpolation
i f ~ i s f i e l d (F , head )
F . ( head ) =
scatteredInterpolant
(CFD_data . ccx ,
CFD_data . ccy ,
CFD_data . ccz ,
CFD_data . (
current_file ) (: ,
current_col ) ) ;
end
% Discretization
dataArray . s t a t i o n s . (
current_plane ) . (
current_var ) . ( headers (
current_file ) ) .(
current_direction ) . ( [ ’
station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] )
= F . ( head ) ( c1 , c2 ,
station_vector ) ;
catch
disp ( ’ e r r o r ␣foam ␣ s t a t i o n s ’
);
end
end
% CONTOURS
case ’ contours ’
% Coordinates
l i m i t s = [ import . ( [ c u r r e n t _ p l a n e
(1) , ’ lims ’ ] ) struct_plot . l i m i t s
];
d i s c r e t i z a t i o n = round ( ( l i m i t s ( 2 )−
l i m i t s ( 1 ) ) ∗2^AMR/ b a s e _ s i z e ) ; %
xy dependant
147
Apéndices
v e r t e x 1 = linspace ( l i m i t s ( 1 ) ,
l i m i t s ( 2 ) , d i s c r e t i z a t i o n +1) ’ ;
c e n t e r 1 = v e r t e x 1 ( 1 : end−1)+d i f f
( vertex1 ) /2;
v e r t e x 2 = linspace ( l i m i t s ( 3 ) ,
l i m i t s ( 4 ) , d i s c r e t i z a t i o n +1) ’ ;
c e n t e r 2 = v e r t e x 2 ( 1 : end−1)+d i f f
( vertex2 ) /2;
coords = [ center1 , center2 ] ;
coord0 = zeros ( d i s c r e t i z a t i o n ) ;
[ coord1 , coord2 ] = meshgrid ( c e n t e r 1
, center2 ) ;
switch current_plane
c a s e ’ xz ’
c1 = coord1 ;
c2 = coord0 ;
c3 = coord2 ;
c a s e ’ yz ’
c1 = coord0 ;
c2 = coord1 ;
c3 = coord2 ;
end
dataArray . c o n t o u r s . ( c u r r e n t _ p l a n e )
. ( current_var ) . coords = coords ;
switch current_var
% Velocities
case ’ v e l o c i t i e s ’
try
% Interpolation
i f ~ i s f i e l d (F , head )
F . ( head ) =
scatteredInterpolant
(CFD_data . x ,
CFD_data . y ,
CFD_data . z ,
CFD_data . (
current_file )
( : , current_col )
148
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
);
end
% Discretization
dataArray . c o n t o u r s . (
current_plane ) . (
current_var ) . (
headers (
current_file ) ) .(
current_direction ) .
contour = F . ( head ) (
c1 , c2 , c3 ) ;
catch
disp ( ’ e r r o r ␣foam ␣
contour ␣ v e l o c i t i e s ’
);
end
% Tke
case ’ tke ’
r e s u l t = zeros (
discretization ,
discretization ,
nDirection ) ;
try
% Interpolation
i f ~ i s f i e l d (F , head )
F . ( head ) =
scatteredInterpolant
(CFD_data . x ,
CFD_data . y ,
CFD_data . z ,
CFD_data . (
current_file )
( : , current_col )
);
end
% Discretization
v e l = F . ( head ) ( c1 , c2 ,
c3 ) ;
% Calculate
149
Apéndices
r e s u l t ( : , : , dd ) = v e l . ∗
vel ;
catch
disp ( ’ e r r o r ␣foam ␣
contour ␣ tke ’ ) ;
end
dataArray . c o n t o u r s . (
current_plane ) . (
c u r r e n t _ v a r ) . contour =
sum( r e s u l t , 3 ) / 2 ;
end
end
end
end
end
end
CFD_data = r m f i e l d (CFD_data , c u r r e n t _ f i l e ) ;
end
end
%% Secondary f u n c t i o n s
function s t r = h e a d e r s ( s t a t )
switch stat
% case ’ inst ’
% str = [ ] ;
c a s e ’UMean ’
s t r = ’ mean ’ ;
c a s e ’ UPrime2Mean ’
s t r = ’RMS ’ ;
c a s e ’ a x i a l ’ %% % % % % % % % % %yz % %p%
l a%
n e%
s t r = ’W’ ;
c a s e ’ r a d i a l ’ %% % % % % % % % % % % % % %
s t r = ’U ’ ;
c a s e ’ t a n g e n t i a l ’ %% % % % % % % % % % %
s t r = ’V ’ ;
% otherwise
% s t r = [ ’Y_’ s t a t ] ; % % % % % % % % % %
150
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
end
end
B.12. IMPORTTEXTFILE
codigo/functions/import/importTextFile.m
% 21/02/2018
% C a r l o s Moreno Montagud
%% Format s t r i n g f o r each l i n e o f t e x t :
% A l l columns : d o u b l e ( %f )
% formatSpec = %f %f %f %f . . . % [ ^ \ n\ r ]
% For more i n f o r m a t i o n , s e e t h e TEXTSCAN documentation .
formatSpec = ’ %∗[^\n\ r ] ’ ;
for i =1: numCols
formatSpec = s t r c a t ( format , formatSpec ) ;
end
if parenthesis
formatSpec = s t r c a t ( ’ ( ’ , formatSpec ) ;
end
151
Apéndices
end
B.13. PLOTCASE
codigo/functions/plot/plotCase.m
% 21/02/2018
% C a r l o s Moreno Montagud
function num_comp = p l o t C a s e ( p l o t t i n g )
%% V a r i a b l e s
num_comp = [ ] ;
% Legends
i f strcmp ( p l o t t i n g . legend_mode , ’ auto ’ )
p l o t t i n g . legend_exp = p l o t t i n g . legend_auto_exp ;
p l o t t i n g . legend_sim = p l o t t i n g . legend_auto_sim ;
e l s e i f strcmp ( p l o t t i n g . legend_mode , ’ manual ’ )
p l o t t i n g . legend_exp = p l o t t i n g . legend_manual_exp ;
p l o t t i n g . legend_sim = p l o t t i n g . legend_manual_sim ;
end
%% P l o t l i n e s
if plotting . lines_flag
num_comp . c e n t e r l i n e = p l o t C e n t e r l i n e ( p l o t t i n g ) ;
end
%% P l o t s t a t i o n s
if plotting . stations_flag
num_comp . s t a t i o n s = p l o t S t a t i o n s S y m ( p l o t t i n g ) ;
end
%% P l o t c o n t o u r s
i f plotting . contours_flag
p l o t C o n t o u r s ( p l o t t i n g ) ; %% % % %
end
%% P l o t monitor
i f p l o t t i n g . monitor_flag
% plotPSD ( p l o t t i n g ) ; % % % % %
152
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
end
end
B.14. PLOTCENTERLINE
codigo/functions/plot/plotCenterline.m
% 15/05/2018
% C a r l o s Moreno Montagud
function c e n t e r l i n e = p l o t C e n t e r l i n e ( p l o t t i n g )
%% V a r i a b l e s
% Plotting
EXP = p l o t t i n g .EXP; numEXP = length (EXP) ;
SIM = p l o t t i n g . SIM ; numSIM = length (SIM) ;
save_flag = plotting . save_flag ;
% save_format = p l o t t i n g . save_format ;
legend_exp = p l o t t i n g . legend_exp ;
legend_sim = p l o t t i n g . legend_sim ;
centerline = [ ] ;
% Config
config = plotConfig () ;
exp_marker_size = c o n f i g . exp_marker_size ;
exp_line_width = c o n f i g . exp_line_width ;
exp_line_style = config . exp_line_style ;
sim_marker_size = c o n f i g . sim_marker_size ;
sim_line_width = c o n f i g . sim_line_width ;
153
Apéndices
%% P l o t
type = p l o t t i n g . c e n t e r l i n e . type ;
nType = length ( type ) ;
for i i =1:nType % ’ v e l o c i t i e s ’
c u r r e n t _ t y p e = type{ i i } ;
struct_type = p l o t t i n g . c e n t e r l i n e . ( current_type ) ;
file_end = struct_type . file_end ;
y_units = s t r u c t _ t y p e . y_units ;
i f a l l ( i s f i e l d ( s t r u c t _ t y p e , { ’ s t a t i s t i c ’ ’ d i r e c t i o n s ’ })
) %% % % % % % % % % %
s t a t i s t i c = struct_type . s t a t i s t i c ;
n S t a t i s t i c = length ( s t a t i s t i c ) ;
for kk =1: n S t a t i s t i c % ’ mean ’
c u r r e n t _ s t a t i s t i c = s t a t i s t i c {kk } ;
d i r e c t i o n s = struct_type . d i r e c t i o n s ;
n D i r e c t i o n s = length ( d i r e c t i o n s ) ;
for l l =1: n D i r e c t i o n s % ’ a x i a l ’
current_direction = directions { l l };
t i t C e l l = [ { c u r r e n t _ s t a t i s t i c } ,{
current_direction } , file_end ,{ ’
centerline ’ }];
t i t u l o = s t r j o i n ( c e l l f u n ( @upperFirst ,
t i t C e l l , ’ un ’ , 0 ) , ’ ␣ ’ ) ;
f i g = figure ( ’ Units ’ , ’ normalized ’ , ’
O u t e r P o s i t i o n ’ , [ 0 0 1 1 ] , ’Name ’ , t i t u l o ,
’ NumberTitle ’ , ’ o f f ’ , ’ PaperUnits ’ , ’
centimeters ’ ) ; % ’ PaperOrientation ’ , ’
landscape ’
154
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
save_fig = true ;
legends = {};
hold on
grid on
for mm=1:numEXP
try
e x p l i n e = EXP{mm} . c e n t e r l i n e . (
current_plane ) . ( current_type ) . (
current_statistic ) .(
current_direction ) . centerline ;
expx = e x p l i n e ( : , 1 ) ;
expy = e x p l i n e ( : , 2 ) ;
%marker_vector = round ( l i n s p a c e ( 1 ,
l e n g t h ( expy ) , 3 0 ) ) ;
plot ( expx , expy , . . .
’ LineStyle ’ , exp_line_style , . . .
’ Color ’ , ’ k ’ , . . .
’ marker ’ , exp_markers (rem(mm−1,
l e n e ) +1) , . . .
’ MarkerSize ’ , exp_marker_size
, . . . ’ MarkerIndices ’ ,
marker_vector , . . .
’ LineWidth ’ , exp_line_width ) ;
l e g e n d s {end+1} = legend_exp {mm} ;
catch
disp ( ’ e r r o r ␣ p l o t ␣ exp ␣ c e n t e r l i n e ’ ) ;
end
end
for mm=1:numSIM
%f i e l d x = s t r c a t ( ’ SIM{mm} . c e n t e r l i n e
. ’ , c u r r e n t _ t y p e , ’ . coords ’ ) ;
%f i e l d y = s t r c a t ( ’ SIM{mm} . c e n t e r l i n e
. ’ , current_type , ’ . ’ ,
current_statistic , ’. ’ ,
current_direction , ’ . line ’) ;
try
simx = SIM{mm} . c e n t e r l i n e . (
current_plane ) . ( current_type ) .
coords ;
155
Apéndices
simy = SIM{mm} . c e n t e r l i n e . (
current_plane ) . ( current_type ) . (
current_statistic ) .(
current_direction ) . centerline ;
marker_vector = round ( linspace ( 1 ,
length ( simy ) , primes (rem(mm−1,
l e n p ) +1) ) ) ;
plot ( simx , simy , . . .
’ L i n e S t y l e ’ , l i n e _ s t y l e (rem(mm
−1, l e n l ) + 1 , : ) , . . .
’ Color ’ , c o l o r s (rem(mm−1, l e n c )
+1) , . . .
’ marker ’ , sim_markers (rem(mm−1,
l e n s ) +1) , . . .
’ m a r k e r s i z e ’ , sim_marker_size
,...
’ M a r k e r I n d i c e s ’ , marker_vector
,...
’ LineWidth ’ , sim_line_width ) ;
l e g e n d s {end+1} = legend_sim {mm} ;
% Numeric comparison % % % % % % %
i f p l o t t i n g . compare_flag
for nn=1:numEXP
e x p l i n e = EXP{nn } .
centerline .(
current_plane ) . (
current_type ) . (
current_statistic ) .(
current_direction ) .
centerline ;
expx = e x p l i n e ( : , 1 ) ;
expy = e x p l i n e ( : , 2 ) ;
% Interpolation
simy2 = interp1 ( simx , simy ,
expx , ’ l i n e a r ’ , ’ e x t r a p ’ )
;
c e n t e r l i n e . ( current_type ) (
nn ,mm) = mean( abs ( expy−
156
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
% i f error save_fig = f a l s e
% Axis p r o p e r t i e s
%l i m i t s = e v a l ( [ ’ s t r u c t _ t y p e . l i m i t s . ’ ,
current_statistic , ’. ’ , current_direction
]) ;
axis ( s t r u c t _ t y p e . l i m i t s . ( c u r r e n t _ s t a t i s t i c
) . ( current_direction ) ) ;
%warning ( ’ ’ ) ; % C l e a r l a s t warning message
%t i t l e ( s t r c a t ( ’ \ b f \ f o n t s i z e { ’ , num2str (
s t a t i o n s _ f o n t s i z e ) , ’} Station z = ’ ,{ ’
’ } , c u r r e n t _ s t a t i o n , ’mm’ ) ) ;
%p b a s p e c t ( p l o t t i n g . v e l o c i t i e s _ a s p e c t _ r a t i o
);
% Legend
legend_str = s t r c a t ( ’ \ bf \ f o n t s i z e { ’ ,
num2str ( l e g e n d _ f o n t s i z e ) , ’ } ’ , strrep (
l e g e n d s , ’_ ’ , ’ \_ ’ ) ) ;
legend ( l e g e n d _ s t r ) ;
% Labels
x l a b = [ ’ \ b f \ f o n t s i z e { ’ , num2str (
x l a b e l _ f o n t s i z e ) , ’ }␣ C e n t e r l i n e ␣ D i s t a n c e
␣ z ␣ [m] ’ ] ;
xlabel ( x l a b ) ;
157
Apéndices
y l a b C e l l = [ { c u r r e n t _ s t a t i s t i c } ,{
c u r r e n t _ d i r e c t i o n } , f i l e _ e n d , { y_units } ] ;
y l a b = [ ’ \ b f \ f o n t s i z e { ’ , num2str (
y l a b e l _ f o n t s i z e ) , ’ }␣ ’ , s t r j o i n ( c e l l f u n (
@upperFirst , y l a b C e l l , ’ un ’ , 0 ) , ’ ␣ ’ ) ] ;
ylabel ( y l a b ) ;
% Save g r a p h i c s
i f s a v e _ f l a g && s a v e _ f i g
p l o t t i n g . p a p e r S i z e = [ 2 1 1 2 ] ; %21 1 4 . 8
p l o t t i n g . paperPos = [ 0 0 . 2 p l o t t i n g .
paperSize ] ;
plotting . results_file = strjoin (
t i t C e l l , ’_ ’ ) ;
end
B.15. PLOTCONTOURS
codigo/functions/plot/plotContours.m
% 16/05/2018
% C a r l o s Moreno Montagud
function p l o t C o n t o u r s ( p l o t t i n g )
%% V a r i a b l e s
% Plotting
% EXP = p l o t t i n g .EXP; numEXP = l e n g t h (EXP) ;
158
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
% Config
config = plotConfig () ;
% exp_marker_size = c o n f i g . exp_marker_size ;
% exp_line_width = c o n f i g . exp_line_width ;
% exp_line_style = config . exp_line_style ;
% sim_marker_size = c o n f i g . sim_marker_size ;
% sim_line_width = c o n f i g . sim_line_width ;
cmap = c o n f i g . cmap ;
%% P l o t
type = p l o t t i n g . c o n t o u r s . type ;
nType = length ( type ) ;
for i i =1:nType % ’ v e l o c i t i e s ’ , ’ t k e ’
c u r r e n t _ t y p e = type{ i i } ;
s t r u c t _ t y p e = eval ( [ ’ p l o t t i n g . c o n t o u r s . ’ , c u r r e n t _ t y p e ] ) ;
159
Apéndices
y_units = s t r u c t _ t y p e . y_units ;
t i t l e _ e n d = s t r u c t _ t y p e . t i t l e _ e n d ; %% % % % % % % % % % % % %
% Figure % % % % % % % % % % % % % % %
i f a l l ( i s f i e l d ( s t r u c t _ t y p e , { ’ s t a t i s t i c ’ ’ d i r e c t i o n s ’ })
) %% % % % % % % % % % % % % % % % % % %
s t a t i s t i c = struct_type . s t a t i s t i c ;
n S t a t i s t i c = length ( s t a t i s t i c ) ;
d i r e c t i o n s = struct_type . d i r e c t i o n s ;
n D i r e c t i o n s = length ( d i r e c t i o n s ) ;
for mm=1:numSIM % c a s e
% Figure
t i t u l o = [ u p p e r F i r s t ( c u r r e n t _ t y p e ) ’ ␣ Contours ␣
’ legend_sim {mm} ] ;
f i g = figure ( ’ Units ’ , ’ normalized ’ , ’
O u t e r P o s i t i o n ’ , [ 0 0 1 1 ] , ’Name ’ , t i t u l o , ’
NumberTitle ’ , ’ o f f ’ , ’ PaperUnits ’ , ’
centimeters ’ ) ; % ’ PaperOrientation ’ , ’
landscape ’
save_fig = true ;
legends = {};
for l l =1: n D i r e c t i o n s % ’ a x i a l ’ , ’ r a d i a l ’ ,
’ tangential ’
current_direction = directions { l l };
i n d e x = sub2ind ( [ n S t a t i s t i c
n D i r e c t i o n s ] , kk , l l ) ;
subAxPlot ( i n d e x ) = subplot ( n D i r e c t i o n s
, n S t a t i s t i c , index ) ;
160
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
%f i e l d x = s t r c a t ( ’ SIM{mm} . c o n t o u r s . ’ ,
current_plane , ’ . ’ , current_type , ’ .
coords ’ ) ;
%f i e l d y = s t r c a t ( ’ SIM{mm} . c o n t o u r s . ’ ,
current_plane , ’ . ’ , current_type , ’ . ’ ,
current_statistic , ’. ’ ,
c u r r e n t _ d i r e c t i o n , ’ . contour ’ ) ;
try
simx = SIM{mm} . c o n t o u r s . (
current_plane ) . ( current_type ) .
coords ;
simy = SIM{mm} . c o n t o u r s . (
current_plane ) . ( current_type ) . (
current_statistic ) .(
c u r r e n t _ d i r e c t i o n ) . contour ;
imagesc ( simx ( : , 1 ) , simx ( : , 2 ) , simy ) ;
axis xy
colormap ( cmap ) ;
colorbar ( ’ n o r t h o u t s i d e ’ ) ;
% Legend and l a b e l s
%t i t l e ( [ ’ \ b f \ f o n t s i z e { ’ , num2str (
t i t l e _ f o n t s i z e ) , ’} ’ , upperFirst
( current_statistic ) ’ ’
upperFirst ( current_direction ) ’
’ title_end ]) ;
%x l a b e l ( [ ’ \ b f \ f o n t s i z e { ’ , num2str (
x l a b e l 1 _ f o n t s i z e ) , ’} Distance
’ , c u r r e n t _ p l a n e ( 1 ) , ’ [m] ’ ] ) ;
%y l a b e l ( [ ’ \ b f \ f o n t s i z e { ’ , num2str (
y l a b e l 1 _ f o n t s i z e ) , ’} Distance
’ , c u r r e n t _ p l a n e ( 2 ) , ’ [m] ’ ] ) ;
pbaspect ( struct_type . aspect_ratio )
;
end
end
end
% Labels
161
Apéndices
p1 = get ( subAxPlot ( 1 ) , ’ p o s i t i o n ’ ) ; % p o s i t i o n
o f upper− l e f t a x e s
p2 = get ( subAxPlot ( end ) , ’ p o s i t i o n ’ ) ; %
p o s i t i o n o f lower−r i g h t a x e s
x t i t = s t r c a t ( ’ \ b f \ f o n t s i z e { ’ , num2str (
t i t l e _ f o n t s i z e ) , ’ }␣ ’ , c e l l f u n ( @upperFirst ,
s t a t i s t i c , ’ un ’ , 0 ) , { ’ ␣ ’ } , t i t l e _ e n d , { ’ ␣ ’ } ,
y_units ) ;
y t i t = s t r c a t ( ’ \ b f \ f o n t s i z e { ’ , num2str (
t i t l e _ f o n t s i z e ) , ’ }␣ ’ , c e l l f u n ( @upperFirst ,
d i r e c t i o n s , ’ un ’ , 0 ) ) ;
x l a b = [ ’ \ b f \ f o n t s i z e { ’ , num2str (
x l a b e l 1 _ f o n t s i z e ) , ’ }␣ D i s t a n c e ␣ ’ ,
c u r r e n t _ p l a n e ( 1 ) , ’ ␣ [m] ’ ] ;
y l a b = [ ’ \ b f \ f o n t s i z e { ’ , num2str (
y l a b e l 1 _ f o n t s i z e ) , ’ }␣ D i s t a n c e ␣ ’ ,
c u r r e n t _ p l a n e ( 2 ) , ’ ␣ [m] ’ ] ;
s u b p l o t L a b e l s ( n D i r e c t i o n s , n S t a t i s t i c , p1 , p2 ,
xlab , ylab , x t i t , y t i t ) ;
% Save g r a p h i c s f i g u r e % % % % %
i f s a v e _ f l a g && s a v e _ f i g
% Properties
subWidth = 9 . 9 ;
subHeight = 9 . 9 ;
o f f s e t = [0 −0.5];
p l o t t i n g . p a p e r S i z e = [ subWidth ∗ n S t a t i s t i c
subHeight ∗ n D i r e c t i o n s ] ;
p l o t t i n g . paperPos = [ o f f s e t subWidth ∗
n S t a t i s t i c subHeight ∗ n D i r e c t i o n s + 1 ] ;
p l o t t i n g . r e s u l t s _ f i l e = s t r c a t ( legend_sim {
mm} , ’_ ’ , ’ c o n t o u r ’ , ’_ ’ , f i l e _ e n d ) ;
162
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
t i t u l o = [ t i t l e _ e n d ’ ␣ ’ legend_sim {mm} ] ;
f i g = figure ( ’ Units ’ , ’ normalized ’ , ’
O u t e r P o s i t i o n ’ , [ 0 0 1 1 ] , ’Name ’ , t i t u l o , ’
NumberTitle ’ , ’ o f f ’ , ’ PaperUnits ’ , ’
centimeters ’ ) ; % ’ PaperOrientation ’ , ’
landscape ’
save_fig = true ;
legends = {};
f i e l d x = s t r c a t ( ’SIM{mm} . c o n t o u r s . ’ ,
c u r r e n t _ p l a n e , ’ . ’ , current_type , ’ . c o o r d s ’ ) ;
f i e l d y = s t r c a t ( ’SIM{mm} . c o n t o u r s . ’ ,
c u r r e n t _ p l a n e , ’ . ’ , current_type , ’ . c o n t o u r ’ ) ;
try
simx = eval ( f i e l d x ) ;
simy = eval ( f i e l d y ) ;
imagesc ( simx ( : , 1 ) , simx ( : , 2 ) , simy ) ;
axis xy
colormap ( cmap ) ;
colorbar
% Legend and l a b e l s
t i t l e ( [ ’ \ b f \ f o n t s i z e { ’ , num2str (
t i t l e _ f o n t s i z e ) , ’ }␣ ’ , t i t l e _ e n d ] ) ;
xlabel ( [ ’ \ b f \ f o n t s i z e { ’ , num2str (
x l a b e l 2 _ f o n t s i z e ) , ’ }␣ D i s t a n c e ␣ ’ ,
c u r r e n t _ p l a n e ( 1 ) , ’ ␣ [m] ’ ] ) ;
ylabel ( [ ’ \ b f \ f o n t s i z e { ’ , num2str (
y l a b e l 2 _ f o n t s i z e ) , ’ }␣ D i s t a n c e ␣ ’ ,
c u r r e n t _ p l a n e ( 2 ) , ’ ␣ [m] ’ ] ) ;
pbaspect ( struct_type . aspect_ratio ) ;
% Save g r a p h i c s f i g u r e
%%%%%%%%%%%%%%%%%%%%
i f s a v e _ f l a g && s a v e _ f i g
% Properties
subWidth = 2 1 ;
subHeight = 1 4 . 8 ;
offset = [0.4 0.4];
163
Apéndices
p l o t t i n g . p a p e r S i z e = [ subWidth
subHeight ] ;
p l o t t i n g . paperPos = [ o f f s e t subWidth
−0.5 subHeight − 0 . 5 ] ;
plotting . results_file = strcat (
legend_sim {mm} , ’_ ’ , ’ c o n t o u r ’ , ’_ ’ ,
file_end ) ;
end
B.16. PLOTSAVE
codigo/functions/plot/plotSave.m
% 31/05/2018
% C a r l o s Moreno Montagud
function p l o t S a v e ( f i g , p l o t t i n g )
164
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
% Export t o a l l f o r m a t s needed
for mm=save_format
format = c har (mm) ;
w r i t e F i l e = s t r c a t ( r e s u l t s _ f i l e n a m e , ’ . ’ , format ( 1 : 3 ) ) ;
w r i t e S e a r c h = s t r c a t ( r e s u l t s _ f i l e n a m e , ’ ∗ . ’ , format ( 1 : 3 ) ) ;
% Check i f a l r e a d y e x i s t
i f exist ( writeFile , ’ f i l e ’ )
f i l e s = dir ( w r i t e S e a r c h ) ;
f i l e n a m e = strrep ( f i l e s ( end ) . name , [ ’ . ’ format ( 1 : 3 ) ] , ’ ’
);
version = s t r 2 d o u b l e ( f i l e n a m e ( end ) ) ;
i f isnan ( version )
version = 0 ;
end
% V e r s i o n s from 1 −10. F u r t h e r v e r s i o n s o v e r w r i t e 10
r e s u l t s _ f i l e n a m e = s t r c a t ( r e s u l t s _ f i l e n a m e , num2str (
version +1) ) ;
end
% Save f i g u r e
format = s t r c a t ( ’−d ’ , format ) ;
print ( r e s u l t s _ f i l e n a m e , format ) ;
% p r i n t ( r e s u l t s _ f i l e n a m e , ’ − d t i f f ’ , ’ − r400 ’ ) % r e s o l u t i o n
end
end
B.17. PLOTSTATIONSSYM
codigo/functions/plot/plotStationsSym.m
% 15/05/2018
% C a r l o s Moreno Montagud
function s t a t i o n s = p l o t S t a t i o n s S y m ( p l o t t i n g )
%% V a r i a b l e s
% Plotting
165
Apéndices
stations = [ ] ;
% Config
config = plotConfig () ;
exp_marker_size = c o n f i g . exp_marker_size ;
exp_line_width = c o n f i g . exp_line_width ;
exp_line_style = config . exp_line_style ;
sim_marker_size = c o n f i g . sim_marker_size ;
sim_line_width = c o n f i g . sim_line_width ;
%% P l o t
type = p l o t t i n g . s t a t i o n s . type ;
nType = length ( type ) ;
for i i =1:nType % ’ v e l o c i t i e s ’ , ’ s p e c i e s ’
c u r r e n t _ t y p e = type{ i i } ;
struct_type = p l o t t i n g . s t a t i o n s . ( current_type ) ;
file_end = struct_type . file_end ;
y_units = s t r u c t _ t y p e . y_units ;
166
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
nPlane = length ( p l a n e ) ;
for j j =1: nPlane % ’ xz ’
current_plane = plane { j j };
i f a l l ( i s f i e l d ( s t r u c t _ t y p e , { ’ s t a t i s t i c ’ ’ d i r e c t i o n s ’ })
) %% % % % % % % % % % % % % % % % % % %
s t a t i s t i c = s t r u c t _ t y p e . s t a t i s t i c ; % ’ mean ’ , ’RMS’
n S t a t i s t i c = length ( s t a t i s t i c ) ;
d i r e c t i o n s = struct_type . d i r e c t i o n s ;
n D i r e c t i o n s = length ( d i r e c t i o n s ) ;
mean_d_temp = zeros (numEXP, numSIM , n D i r e c t i o n s ) ;
RMS_d_temp = zeros (numEXP, numSIM , n D i r e c t i o n s ) ;
for kk =1: n D i r e c t i o n s % ’ a x i a l ’ , ’ r a d i a l ’ , ’
t a n g e n t i a l ’ , ’CH4’
c u r r e n t _ d i r e c t i o n = d i r e c t i o n s {kk } ;
station_nums = s t r u c t _ t y p e . s t a t i o n s ;
n S t a t i o n s = length ( station_nums ) ;
mean_s_temp = zeros (numEXP, numSIM , n S t a t i o n s ) ;
RMS_s_temp = zeros (numEXP, numSIM , n S t a t i o n s ) ;
for l l =1: n S t a t i o n s % [ 5 10 20 30 4 0 ]
c u r r e n t _ s t a t i o n = num2str ( station_nums ( l l )
);
p l o t _ i n d e x = n S t a t i o n s+2− l l ;
subAxPlot ( p l o t _ i n d e x ) = subplot ( n S t a t i o n s
+1 ,1 , p l o t _ i n d e x ) ;
legends = {};
167
Apéndices
hold on
grid on
for mm=1:numEXP
try
meann = EXP{mm} . s t a t i o n s . (
current_plane ) . ( current_type ) . (
s t a t i s t i c {1}) . (
current_direction ) . ( [ ’ station_ ’
, c u r r e n t _ s t a t i o n , ’mm’ ] ) ;
mean_expx = −meann ( : , 1 ) ; %
%%%%%%%%%%
mean_expy = meann ( : , 2 ) ;
RMS = EXP{mm} . s t a t i o n s . (
current_plane ) . ( current_type ) . (
s t a t i s t i c {2}) . (
current_direction ) . ( [ ’ station_ ’
, c u r r e n t _ s t a t i o n , ’mm’ ] ) ;
RMS_expx = RMS( : , 1 ) ;
RMS_expy = RMS( : , 2 ) ;
% mean
yyaxis l e f t
%marker_vector = round ( l i n s p a c e ( 1 ,
l e n g t h ( mean_expy ) , 3 0 ) ) ;
plot ( mean_expx , mean_expy , . . .
’ LineStyle ’ , exp_line_style , . . .
’ Color ’ , ’ k ’ , . . .
’ marker ’ , exp_markers (rem(mm−1,
l e n e ) +1) , . . .
’ MarkerSize ’ , exp_marker_size
, . . . ’ MarkerIndices ’ ,
marker_vector , . . .
’ LineWidth ’ , exp_line_width ) ;
% RMS
yyaxis right
%marker_vector = round ( l i n s p a c e ( 1 ,
l e n g t h (RMS_expy) , 3 0 ) ) ;
plot (RMS_expx , RMS_expy , . . .
168
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
’ LineStyle ’ , exp_line_style , . . .
’ Color ’ , ’ k ’ , . . .
’ marker ’ , exp_markers (rem(mm−1,
l e n e ) +1) , . . .
’ MarkerSize ’ , exp_marker_size
, . . . ’ MarkerIndices ’ ,
marker_vector , . . .
’ LineWidth ’ , exp_line_width ) ;
l e g e n d s {end+1} = legend_exp {mm} ;
end
end
for mm=1:numSIM
try
simx = SIM{mm} . s t a t i o n s . (
current_plane ) . ( current_type ) .
coords ;
mean_simx = −simx ( simx >=0) ;
RMS_simx = simx ( simx >=0) ;
simy1 = SIM{mm} . s t a t i o n s . (
current_plane ) . ( current_type ) . (
s t a t i s t i c {1}) . (
current_direction ) . ( [ ’ station_ ’
, c u r r e n t _ s t a t i o n , ’mm’ ] ) ;
mean_simy = simy1 ( simx >=0) ;
simy2 = SIM{mm} . s t a t i o n s . (
current_plane ) . ( current_type ) . (
s t a t i s t i c {2}) . (
current_direction ) . ( [ ’ station_ ’
, c u r r e n t _ s t a t i o n , ’mm’ ] ) ;
RMS_simy = simy2 ( simx >=0) ;
% mean
yyaxis l e f t
marker_vector = round ( linspace ( 1 ,
length ( mean_simy ) , primes (rem(mm
−1, l e n p ) +1) ) ) ;
plot ( mean_simx , mean_simy , . . .
’ L i n e S t y l e ’ , l i n e _ s t y l e (rem(mm
−1, l e n l ) + 1 , : ) , . . .
169
Apéndices
’ Color ’ , c o l o r s (rem(mm−1, l e n c )
+1) , . . .
’ marker ’ , sim_markers (rem(mm−1,
l e n s ) +1) , . . .
’ m a r k e r s i z e ’ , sim_marker_size
,...
’ M a r k e r I n d i c e s ’ , marker_vector
,...
’ LineWidth ’ , sim_line_width ) ;
% RMS
yyaxis right
marker_vector = round ( linspace ( 1 ,
length (RMS_simy) , primes (rem(mm
−1, l e n p ) +1) ) ) ;
plot (RMS_simx , RMS_simy , . . .
’ L i n e S t y l e ’ , l i n e _ s t y l e (rem(mm
−1, l e n l ) + 1 , : ) , . . .
’ Color ’ , c o l o r s (rem(mm−1, l e n c )
+1) , . . .
’ marker ’ , sim_markers (rem(mm−1,
l e n s ) +1) , . . .
’ m a r k e r s i z e ’ , sim_marker_size
,...
’ M a r k e r I n d i c e s ’ , marker_vector
,...
’ LineWidth ’ , sim_line_width ) ;
l e g e n d s {end+1} = legend_sim {mm} ;
% Numeric comparison
i f p l o t t i n g . compare_flag
for nn=1:numEXP
meann = EXP{nn } . s t a t i o n s . (
current_plane ) . (
current_type ) . (
s t a t i s t i c {1}) . (
current_direction ) . ( [ ’
station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] ) ;
170
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
mean_expx = −meann ( : , 1 ) ; %
%%%%%%%%%%
mean_expy = meann ( : , 2 ) ;
RMS = EXP{nn } . s t a t i o n s . (
current_plane ) . (
current_type ) . (
s t a t i s t i c {2}) . (
current_direction ) . ( [ ’
station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] ) ;
RMS_expx = RMS( : , 1 ) ;
RMS_expy = RMS( : , 2 ) ;
% mean
mean_simy2 = interp1 (
mean_simx , mean_simy ,
mean_expx , ’ l i n e a r ’ , ’
extrap ’ ) ;
s t a t i o n s . ( current_type ) .
mean . ( c u r r e n t _ d i r e c t i o n
) . ( [ ’ station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] ) (
nn ,mm) = mean( abs (
mean_expy−mean_simy2 )
. / ( abs ( mean_expy )+abs (
mean_simy2 ) ) ) ; %% %
mean_s_temp ( nn ,mm, l l ) =
s t a t i o n s . ( current_type )
. mean . (
current_direction ) . ( [ ’
station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] ) (
nn ,mm) ;
% RMS
RMS_simy2 = interp1 (
RMS_simx , RMS_simy ,
RMS_expx , ’ l i n e a r ’ , ’
extrap ’ ) ;
171
Apéndices
s t a t i o n s . ( current_type ) .
RMS. ( c u r r e n t _ d i r e c t i o n )
. ( [ ’ station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] ) (
nn ,mm) = mean( abs (
RMS_expy−RMS_simy2) . / (
abs (RMS_expy)+abs (
RMS_simy2) ) ) ; %% %
RMS_s_temp( nn ,mm, l l ) =
s t a t i o n s . ( current_type )
.RMS. ( c u r r e n t _ d i r e c t i o n
) . ( [ ’ station_ ’ ,
c u r r e n t _ s t a t i o n , ’mm’ ] ) (
nn ,mm) ;
end
i f l l==n S t a t i o n s
s t a t i o n s . ( current_type ) .
mean . ( c u r r e n t _ d i r e c t i o n
) . t o t a l = mean(
mean_s_temp , 3 ) ;
mean_d_temp ( : , : , kk ) =
s t a t i o n s . ( current_type )
. mean . (
current_direction ) .
total ;
s t a t i o n s . ( current_type ) .
RMS. ( c u r r e n t _ d i r e c t i o n )
. t o t a l = mean(
RMS_s_temp , 3 ) ;
RMS_d_temp ( : , : , kk ) =
s t a t i o n s . ( current_type )
.RMS. ( c u r r e n t _ d i r e c t i o n
) . total ;
% d i s p ( [ ’ s t a t i o n s [ mean
RMS] ’ c u r r e n t _ t y p e ’ ’ c u r r e n t _ d i r e c t i o n ] ) ;
% disp ([ stations .(
c u r r e n t _ t y p e ) . mean . ( c u r r e n t _ d i r e c t i o n ) . t o t a l , s t a t i o n s . (
c u r r e n t _ t y p e ) .RMS. ( c u r r e n t _ d i r e c t i o n ) . t o t a l ] ) ;
172
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
i f kk==n D i r e c t i o n s
s t a t i o n s . ( current_type
) . mean . t o t a l = mean
( mean_d_temp , 3 ) ;
s t a t i o n s . ( current_type
) .RMS. t o t a l = mean(
RMS_d_temp, 3 ) ;
disp ( [ ’ s t a t i o n s ␣ ’
current_type ’ ␣ [
mean␣RMS] ’ ] ) ;
disp ( [ s t a t i o n s . (
c u r r e n t _ t y p e ) . mean .
total , stations . (
c u r r e n t _ t y p e ) .RMS.
total ]) ;
end
end
% d i s p ( [ ’ s t a t i o n s [ mean RMS] ’
c u r r e n t _ t y p e ’ ’ c u r r e n t _ d i r e c t i o n ’ ’ c u r r e n t _ s t a t i o n ’mm
’]) ;
% disp ( [ s t a t i o n s . ( current_type
) . mean . ( c u r r e n t _ d i r e c t i o n ) . ( [ ’ s t a t i o n _ ’ , c u r r e n t _ s t a t i o n , ’mm
’ ] ) s t a t i o n s . ( c u r r e n t _ t y p e ) .RMS. ( c u r r e n t _ d i r e c t i o n ) . ( [ ’
s t a t i o n _ ’ , c u r r e n t _ s t a t i o n , ’mm’ ] ) ] ) ;
end
end
end
% i f error save_fig = f a l s e
% Axis p r o p e r t i e s
%e v a l ( s t r c a t ( ’ l i m i t s = p l o t v a r . a x i s . ’ ,
velDir , ’ ; ’ ) ) ;
xlim ( s t r u c t _ t y p e . l i m i t s . x l i m i t ) ;
yyaxis l e f t
ylim ( s t r u c t _ t y p e . l i m i t s . mean . (
current_direction ) ) ;
173
Apéndices
yyaxis right
ylim ( s t r u c t _ t y p e . l i m i t s .RMS. (
current_direction ) ) ;
%warning ( ’ ’ ) ; % C l e a r l a s t warning message
t i t l e ( s t r c a t ( ’ \ b f \ f o n t s i z e { ’ , num2str (
s t a t i o n s _ f o n t s i z e ) , ’ } S t a t i o n ␣ z ␣=␣ ’ ,{ ’ ␣ ’
} , c u r r e n t _ s t a t i o n , ’mm’ ) ) ;
%p b a s p e c t ( p l o t t i n g . v e l o c i t i e s _ a s p e c t _ r a t i o
);
end
% Legend
l e g e n d s = s t r c a t ( ’ \ b f \ f o n t s i z e { ’ , num2str (
l e g e n d _ f o n t s i z e ) , ’ } ’ , strrep ( l e g e n d s , ’_ ’ , ’ \_
’));
l g = legend ( l e g e n d s , ’ O r i e n t a t i o n ’ , ’ h o r i z o n t a l ’
) ; % create legend
l p o s = get ( l g , ’ p o s i t i o n ’ ) ; % g e t l e g e n d
position
subAxPlot ( 1 ) = subplot ( n S t a t i o n s +1 ,1 ,1) ; %
create subplot
set ( subAxPlot ( 1 ) , ’ V i s i b l e ’ , ’ o f f ’ ) ; % h i d e
subplot
pos = get ( subAxPlot ( 1 ) , ’ p o s i t i o n ’ ) ; % g e t
subplot position
set ( l g , ’ p o s i t i o n ’ , [ pos ( 1 ) +(pos ( 3 )−l p o s ( 3 ) ) /2
pos ( 2 )+pos ( 4 ) /4− l p o s ( 4 ) l p o s ( 3 ) l p o s ( 4 ) ∗ 2 ] )
;
% Labels
p1 = get ( subAxPlot ( 2 ) , ’ p o s i t i o n ’ ) ; % p o s i t i o n
o f upper− l e f t a x e s
p2 = get ( subAxPlot ( end ) , ’ p o s i t i o n ’ ) ; %
p o s i t i o n o f lower−r i g h t a x e s
posOuter = [ p1 ( 1 ) p2 ( 2 ) p2 ( 1 )+p2 ( 3 )−p1 ( 1 ) p1
( 2 )+p1 ( 4 )−p2 ( 2 ) ] ;
hAxOuter = axes ( ’ p o s i t i o n ’ , posOuter , ’ c o l o r ’ , ’
none ’ , ’ v i s i b l e ’ , ’ o f f ’ ) ; % an o u t e r a x i s f o r
global p
174
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
x l a b = [ ’ \ b f \ f o n t s i z e { ’ , num2str (
x l a b e l _ f o n t s i z e ) , ’ }␣ R a d i a l ␣ D i s t a n c e ␣ ’ ,
c u r r e n t _ p l a n e ( 1 ) , ’ ␣ [m] ’ ] ;
hX=text ( 0 . 5 , − 0 . 0 5 , xlab , ’ r o t a t i o n ’ , 0 , ’
horizontalalignment ’ , ’ center ’ , ’
v e r t i c a l a l i g n m e n t ’ , ’ top ’ ) ;
y l a b C e l l = s t r j o i n ({ t i t u l o , y_units } , ’ ␣ ’ ) ;
y l a b F i r s t = c e l l f u n ( @upperFirst , s t a t i s t i c , ’ un ’
,0) ;
y l a b = s t r c a t ( ’ \ b f \ f o n t s i z e { ’ , num2str (
y l a b e l _ f o n t s i z e ) , ’ } ’ ,{ ’ ␣ ’ } , y l a b F i r s t , { ’ ␣ ’ } ,
ylabCell ) ;
hY1=text ( − 0 . 0 5 , 0 . 5 , y l a b {1} , ’ r o t a t i o n ’ , 9 0 , ’
horizontalalignment ’ , ’ center ’ , ’
v e r t i c a l a l i g n m e n t ’ , ’ bottom ’ ) ;
hY2=text ( 1 . 0 5 , 0 . 5 , y l a b {2} , ’ r o t a t i o n ’ , −90 , ’
horizontalalignment ’ , ’ center ’ , ’
v e r t i c a l a l i g n m e n t ’ , ’ bottom ’ ) ;
% Save g r a p h i c s f i g u r e
i f s a v e _ f l a g && s a v e _ f i g
% Properties
subHeight = 6 ;
subWidth = 2 5 ;
o f f s e t = [ −0.5 0 ] ;
p l o t t i n g . p a p e r S i z e = [ subWidth subHeight ∗ (
n S t a t i o n s −1) ] ;
p l o t t i n g . paperPos = [ o f f s e t subWidth
subHeight ∗ ( n S t a t i o n s −3/4) ] ;
p l o t t i n g . r e s u l t s _ f i l e = s t r j o i n ( t i t C e l l , ’_
’);
175
Apéndices
end
end
end
%% Secondary f u n c t i o n s
function s t r = u p p e r F i r s t ( data )
s t r = [ upper ( data ( 1 ) ) data ( 2 : end ) ] ;
end
B.18. SUBPLOTLABELS
codigo/functions/plot/subplotLabels.m
% 24/05/2018
% C a r l o s Moreno Montagud
176
B. CÓDIGO MATLAB PARA POSPROCESAR RESULTADOS
end
177
Documento II
PLIEGO DE
CONDICIONES
Índice de documento
1 Introducción
3 Seguridad estructural
4 Superficie y cubicación
6 Disposiciones generales
7 Iluminación de emergencia
11 Electricidad estática
181
Capítulo 1
Introducción
183
Capítulo 2
Obligaciones y derechos de
los trabajadores
186
Capítulo 3
Seguridad estructural
187
Capítulo 4
Superficie y cubicación
Los locales de trabajo reunirán las siguientes condiciones mínimas: Tres metros de
altura desde el piso al techo. Dos metros cuadrados de superficie por cada trabajador.
Diez metros cúbicos para cada trabajador. No obstante, en los establecimientos comer-
ciales, de servicios y locales destinados a oficinas y despachos, la altura a la que se
refiere el primer punto puede quedar reducida hasta 2:5 m, pero respetando la cubi-
cación que establece el tercer punto, y siempre que el aire se renueve suficientemente.
Para el cálculo de la superficie y el volumen, no se tendrán en cuenta los espacios
ocupados por máquinas, aparatos, instalaciones y materiales.
189
Capítulo 5
191
Capítulo 6
Disposiciones generales
193
Capítulo 7
Iluminación de emergencia
195
Capítulo 8
Ventilación, temperatura y
humedad
En los lugares de trabajo y sus anexos se mantendrá, por medios naturales o ar-
tificiales, unas condiciones atmosféricas adecuadas, evitando el aire viciado, exceso de
calor o de frío, humedad o sequía y los olores desagradables. En ningún caso, el anhídri-
do carbónico ambiental, podrá sobrepasar la proporción de 50/10000, y el monóxido de
carbono la de 1/10000. En los locales de trabajo cerrados, el suministro de aire fresco
y limpio por hora y por trabajador, será de al menos 30𝑚3 , salvo que se efectúe una
renovación total del aire varias veces por hora, no inferior a seis veces para trabajos
sedentarios, ni a diez veces para trabajos que exijan un esfuerzo físico superior al nor-
mal. En el otro extremo, la circulación de aire en locales cerrados se acondicionará de
modo que los trabajadores no estén expuestos a corrientes molestas y que la velocidad
del aire no exceda de 15𝑚/𝑚𝑖𝑛 con temperatura normal, ni de 45𝑚/𝑚𝑖𝑛 en ambientes
extremadamente calurosos. En los centros de trabajo expuestos a altas y bajas tem-
peraturas, serán evitadas las variaciones bruscas por el medio que se considere más
eficaz. Cuando la temperatura sea extremadamente distinta entre los lugares de tra-
bajo, deberán existir locales de paso para que los operarios se adapten gradualmente
de unas condiciones a las otras. De acuerdo con todo lo anterior, se fijan como lími-
tes de temperatura y humedad en locales y para los distintos trabajos, siempre que el
procedimiento de fabricación lo permita, los siguientes:
debajo del 50 %. En aquellos trabajos en los que por exigencias del proceso los locales
estén sometidos a calor o frío extremo, se eliminará la permanencia de los operarios,
estableciendo los turnos adecuados en cada caso.
198
Capítulo 9
Ruidos, vibraciones y
trepidaciones
199
Capítulo 10
En las instalaciones y equipos eléctricos, para la protección de las personas contra los
contactos con partes habitualmente en tensión, se adoptaran algunas de las siguientes
prevenciones:
Se alejarán las partes activas de la instalación a distancia suficiente del lugar don-
de las personas habitualmente se encuentran o circulan, para evitar un contacto
fortuito o por la manipulación de objetos conductores, cuando éstos puedan ser
utilizados cerca de estas partes activas de la instalación.
Se recubrirán las partes activas con el aislamiento apropiado, que permita con-
servar indefinidamente las propiedades del conductor y que limiten la corriente
de contacto a un valor inocuo para las personas.
Se interpondrán obstáculos que impidan todo contacto accidental con las partes
activas de la instalación. Los obstáculos de protección deben estar fijados en
forma segura y ser capaces de resistir los esfuerzos mecánicos usuales.
Para la protección contra los riesgos de contacto con las masas de las instalaciones que
puedan qudar accidentalmente con tensión, se adoptarán, en corriente alterna, uno o
varios de los siguientes dispositivos de seguridad.
Puesta a tierra de las masas. Las masas deben estar unidas eléctricamente a una
toma de tierra o a un conjunto de tomas de tierra interconectadas, que tengan
una resistencia apropiada. Las instalaciones, tanto con neutro aislado como con
neutro a tierra, deben estar permanentemente controladas por un dispositivo
que indique automáticamente la existencia de cualquier defecto de aislamiento, o
que separe automáticamente la instalación o la parte defectuosa de la fuente de
energía de la que se alimenta.
201
Pliego de condiciones
202
Capítulo 11
Electricidad estática
La humedad relativa del aire se mantendrá siempre con un valor por debajo del
50 %.
Las cargas de electricidad estática que puedan acumularse en los cuerpos metá-
licos, serán neutralizadas por medio de la conexión de conductores a tierra. La
forma de realizar estas conexiones puede variar dependiendo del tipo de máquina.
203
Capítulo 12
Recomendaciones sobre
materias inflamables
205
Pliego de condiciones
206
Capítulo 13
Prevención y extinción de
incendios
En los centros de trabajo que ofrezcan peligro de incendios, con o sin explosión,
se adoptarán las prevenciones que se indican a continuación, combinando su empleo
con la protección general más próxima que puedan prestar los servicios públicos contra
incendios:
207
Pliego de condiciones
208
Documento III
PRESUPUESTO
Índice de documento
1 Introducción
4 Costes generales
5 Presupuesto total
211
Capítulo 1
Introducción
Con todo, teniendo en cuenta que el Trabajo Final de Máster únicamente ha cons-
tado de un estudio computacional, se realizará un presupuesto parcial de la misma y un
presupuesto total. Adicionalmente, se aplican los precios correspondientes a las tarifas
legales vigentes en el caso de ser conocidas o aplicado estimaciones de no serlo, siempre
manteniendo una coherencia en la misma. El presupuesto total se define como la suma
de los presupuestos parciales, con un 50 % adicional en concepto de los medios auxilia-
res y/o costes imprevistos. Finalmente, se le aplica el IVA vigente. A continuación se
muestran cada uno de los costes descritos anteriormente.
213
Capítulo 2
En esta fase han participado como recursos humanos un total de tres profesionales:
Además de los recursos humanos, se deben incluir en este apartado los costes de
amortización del equipo informático utilizado como herramienta de trabajo. Esta etapa
ha abarcado un período de seis meses. También se incluye en este apartado el período
de amortización del software utilizado.
215
Presupuesto
Una vez obtenidos los datos de tasa horaria recogidas en la tabla 2.1, estimando
las horas de dedicación de cada uno de los integrantes del trabajo se puede obtener la
cuantía del coste de cada uno de ellos, ver tabla 2.2:
216
Capítulo 3
El coste de amortización por año del equipo informático se calcula como la ecuación
3.1.
𝑉𝐶 − 𝑉𝑅
𝑎= (3.1)
𝑛
Donde:
𝑉𝐶 es el valor de compra del equipo. Para el caso del portátil particular es de 670
EUR, y para la estación de trabajo 5000 EUR.
𝑉𝑅 es el valor residual del equipo. Para el caso de un portátil y una estación de
trabajo se estima en un 20 % del valor de compra.
𝑛 es el periodo de amortización. Para el caso de un portátil particular se ha
estimado en 5 años, y para la estación de trabajo 4.
Para calcular el coste de las horas de cálculo se estima el precio por núcleo y hora de
0.01 EUR. Las simulaciones se han realizado en aproximadamente 5 meses, utilizando
Altamira, 2 cuentas de RIGEL con 64 núcleos cada una y la estación de trabajo con
32 núcleos. En la tabla 3.1 se recoge el desglose de núcleos y horas.
Por otra parte se incluye el coste de electricidad de los equipos en este apartado por
la brevedad del mismo. Se ha calculado teniendo en cuenta la aproximación de horas de
dedicación en la tabla 2.2 y que el precio de la electricidad ronda los 0,15 EUR/kWh.
En la tabla 3.2 se recogen todos y cada uno de los costes de licencia y equipos
necesarios en el desarrollo del presente trabajo. Los precios anuales se han multiplicado
por 0.67 debido a que el TFM ha durado 8 meses.
218
Capítulo 4
Costes generales
En este apartado se tienen en cuenta los costes relacionados con los imprevistos y
medios auxiliares. Se estiman como el 5 % de la suma de los costes de personal y de
equipos, según lo recogido en la tabla 4:
219
Capítulo 5
Presupuesto total
El presupuesto total queda definido como la suma de los costes de recursos humanos,
costes de equipo y costes generales, aplicándoles el IVA correspondiente. Pueden verse
en la tabla 5:
Por lo tanto, puede decirse que el presupuesto total estimado para la realización del
Trabajo Final de Máster es de VEINTICINCO MIL QUINIENTOS OCHO EUROS Y
TRENTA Y NUEVE CÉNTIMOS.
221