AUTOMATIZACIÓN DEL DISEÑO CFD TRIDIMENSIONAL DE PEQUEÑOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR

AUTOMATIZACIÓN DEL DISEÑO CFD TRIDIMENSIONAL DE PEQUEÑOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR FEBRERO 2025 Antonio Orti Hortelano DIRECTOR DEL TRABAJO FIN DE GRADO : Andrés Sebastián Herrera Rubén Abbas Cámara TRABAJO FIN DE GRADO PARA LA OBTENCIÓN DEL TÍTULO DE GRADUADO EN INGENIERÍA EN TECNOLOGÍAS INDUSTRIALES AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Antonio Orti Hortelano Trabajo de Fin de Grado para la obtenci´on del t´ıtulo Grado en Ingenier´ıa en Tecnolog´ıas Industriales Directores: Andr´es Sebasti´an Herrera Rub´en Abbas C´amara UNIVERSIDAD POLIT´ ECNICA DE MADRID Escuela T´ecnica Superior de Ingenieros Industriales 2025 “Muchas veces habr´a de proseguir la fuerza lo que hubiera comenzado y no consiguiese la prudencia” Hern´an Cort´es AGRADECIMIENTOS A Dios y a mi familia. A Andr´es y a Rub´en, por estos meses de intenso trabajo, sin su ayuda y gu´ıa, el proyecto no habr´ıa sido posible. i RESUMEN EJECUTIVO RESUMEN EJECUTIVO El fin de este Trabajo de Fin de Grado titulado “Automatizaci´on del dise˜no CFD tridimensional de peque˜nos turbocompresores radiales para bomba de calor” consiste en desarrollar una metodolog´ıa de automatizaci´on del proceso de simulaci´on num´erica de micro-turbocompresores centr´ıfugos. El proyecto se enmarca dentro de los trabajos de investigaci´on del grupo “Thermal Energy for Sustainability” de la Universidad Polit´ecnica de Madrid, en la l´ınea de desarrollo y dise˜no de micro-turbocompresores centr´ıfugos para bomba de calor. Se trata de compresores abiertos con difusor sin ´alabes. Seg´un la EEA (European Environment Agency), el 50 % de la energ´ıa consumida en la Uni´on Europea es destinada a generaci´on de calor y fr´ıo industrial y residencial, y el 70 % de esta energ´ıa se obtiene de combustibles f´osiles. La bomba de calor se presenta como una de las principales tecnolog´ıas para acelerar la descarbonizaci´on y reducir el impacto en los procesos de generaci´on de calor y refrigeraci´on en el medio ambiente. Concretamente, los micro-turbocompresores para bombas de calor a peque˜na escala bajo estudio emplean fluidos no ideales u org´anicos con el objetivo de reducir las emisiones contaminantes debidas al actual empleo de refrigerantes con alto potencial de calentamiento global (GWP). Una de las principales ventajas de la bomba de calor es su alto rendimiento medido por el COP (Coefficient of Performance), que puede ser mayor que uno, pues es un sistema que transfiere calor de un foco fr´ıo a un foco caliente, siendo el COP la relaci´on entre la potencia que se transmite en forma de calor y la potencia que se aporta al compresor. Estos sistemas pueden, por tanto, ofrecer eficiencias mucho mayores que los sistemas convencionales de generaci´on de calor. Por ello, se hace esencial el desarrollo de compresores eficientes que reduzcan la cantidad de trabajo aportado y mejoren la eficiencia de esta tecnolog´ıa. El proyecto consiste en el desarrollo de un c´odigo en MATLAB que se divide en dos partes principales. Por un lado, el desarrollo de un c´odigo que defina y dise˜ne geometr´ıas tridimensionales de compresores centr´ıfugos, definiendo todas las variables correspondientes. Posteriormente, el c´odigo de automatizaci´on que realice la conexi´on MATLAB-ANSYS Workbench (WB) y permita la automatizaci´on de la simulaci´on de las geometr´ıas. En este proyecto, la parte del c´odigo de dise˜no y definici´on de geometr´ıas tomar´a como punto de partida variables unidimensionales obtenidas a partir de c´odigos de dise˜no l´ınea media. Para definir la geometr´ıa, se definir´an principalmente el n´umero de ´alabes, las curvas que definen el canal meridional, las curvas que definen las distribuciones de ´angulos del ´alabe, la distribuci´on de espesores del ´alabe y la posici´on del ´alabe intermedio (si lo hay). Las diferentes curvas se definen mediante curvas de B´ezier, que a su vez vienen definidas por puntos de control que son los que se calcular´an en el c´odigo de dise˜no de geometr´ıas de MATLAB. Mediante la modificaci´on de estos puntos de control y las variables wrap yrake (Figura 0.2), se pueden obtener diferentes geometr´ıas que parten de las mismas variables unidimensionales obtenidas del dise˜no de l´ınea media, que permanecen constantes con la idea de obtener el dise˜no tridimensional ´optimo completando el dise˜no de l´ınea media. Para poder automatizar la simulaci´on, se parte de un proyecto previamente creado (proyecto base) en Ansys WB, donde existe una topolog´ıa (donde se aplica la malla) que se ii Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR aplicar´a a las nuevas geometr´ıas que se simulen y que se podr´a ir modificando, as´ı como unas condiciones de contorno que tambi´en se podr´an modificar seg´un el problema objeto de estudio y un sistema de obtenci´on de resultados configurado. Ansys WB cuenta con los archivos journal, se trata de archivos de texto que contienen scripts con secuencia de comandos escritos en API y Python que permitir´an conectar MATLAB con ANSYS WB. En estos archivos se escriben las instrucciones oportunas para realizar las siguientes acciones principales: Abrir el proyecto base. Importar la geomtr´ıa dise˜nada en MATLAB. Moficicar los par´ametros pertinentes. Actualizar el proyecto (comando update). Guardar los resultados como otro proyecto, quedando el proyecto base intacto y listo para una nueva geomtr´ıa. Para la validaci´on del c´odigo de dise˜no de geometr´ıas y c´odigo de automatizaci´on, se ha seleccionado el compresor centr´ıfugo dise˜nado por Schiffmann et al [1] con valores experimentales. Se modificar´an el canal meridional, la distribuci´on de ´angulos y los valores del wrap yrake angles (principales variables tridimensionales) con el fin de analizar los efectos que estas variaciones tienen sobre la evoluci´on de las variables termodin´amicas, la relaci´on de compresi´on y la eficiencia; simulando un total de diecis´eis geometr´ıas de manera autom´atica. Se tomar´an como punto de partida, por tanto, valores unidimensionales aportados por el estudio de l´ınea media del profesor Schiffmann y los valores de wrap yrake angles ser´an definidos. En la Figura 0.1 se muestra el ejemplo de un rodete dise˜nado con el c´odigo de dise˜no de geometr´ıas de MATLAB a partir de estos valores unidimensionales. Las variables de entrada y un esquema del proceso se muestran en la Figura 0.2. (a) (b) Figura 0.1: Ejemplo de rodete tridimensional dise˜nado c´odigo de dise˜no de MATLAB. Antonio Orti Hortelano iii RESUMEN EJECUTIVO Figura 0.2: Variables de entrada y proceso caso de estudio Para analizar los efectos del canal meridional, se simulan 4 geometr´ıas con distintos canales meridionales (Figura 0.3) aumentando la curvatura del canal, pero con igual distribuci´on de ´angulo y mismos wrap yrake angles. Para analizar el efecto de la distribuci´on de ´angulo, se modifican las distribuciones de ´angulo disminuyendo la curvatura del ´alabe (Figura 0.4); se utiliza un mismo canal meridional as´ı como un mismo wrap yrake angle. Figura 0.3: Distintos canales meridionales. Figura 0.4: Distintas distribuciones de ´angulo. Se ha observado que a medida que se curva el canal meridional, disminuye la entalp´ıa a la salida del rodete, lo que disminuye la relaci´on de compresi´on. Asimismo, la entrop´ıa aumenta con esta curvatura, dando lugar a una reducci´on de la eficiencia a medida que se curva el canal meridional. En la Figura 0.5 se muestra la evoluci´on de la presi´on de la geometr´ıa con el canal menos curvado y del m´as curvado, donde se aprecia la disminuci´on de la presi´on a la salida del rodete. iv Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR (a) (b) Figura 0.5: Evoluci´on de la presi´on en el canal meridional. Las diferentes distribuciones de ´angulo tienen un impacto crucial sobre la relaci´on de compresi´on (Figura 0.7). A medida que disminuye la curvatura del ´alabe, se producen menos p´erdidas por rozamiento del fluido con las paredes, siendo m´as eficiente la transferencia de energ´ıa mec´anica al fluido en forma de entalp´ıa y vi´endose aumentada la relaci´on de compresi´on. En cambio, la eficiencia sigue un comportamiento menos predecible (Figura 0.6). Se produce un aumento de entrop´ıa a la salida del rodete en las ´ultimas geometr´ıas frente a un aumento de entalp´ıa. En la ´ultima geometr´ıa, que es la que mayor relaci´on de compresi´on presenta, el rendimiento cae por el efecto m´as significativo del aumento de entrop´ıa. Figura 0.6: Eficiencias total a total. Figura 0.7: Relaci´on de compresi´on total a total. Para analizar el efecto del wrap yrake angles, se han simulado diferentes geometr´ıas variando el valor del wrap angle de mayor a menor con un rake, distribuci´on de ´angulos y canal meridional constantes. Para analizar el efecto del rake angle, se simulan tres geometr´ıas con el mismo wrap angle, distribuci´on de ´angulos y canal meridional. El rake de la primera geometr´ıa tiene un valor de 30º, el segundo de 0ºy el tercero de -30º. El aumento del wrap angle conlleva un aumento de la longitud del canal de paso del fluido a trav´es del compresor, lo que implica mayor rozamiento con las paredes. De los resultados de las simulaciones, se observa c´omo disminuye la entrop´ıa a la salida del compresor y el aumento de entalp´ıa. Asimismo, se aprecia un aumento del n´umero de Mach en los ´angulos wrap m´as grandes, lo que supone mayores p´erdidas por exceso de aceleraci´on del Antonio Orti Hortelano v RESUMEN EJECUTIVO flujo. En la Figura 0.8 se muestra la evoluci´on del n´umero de Mach, y c´omo el flujo se ve acelerado en la geometr´ıa con un wrap de 80o. (a) (b) Figura 0.8: Evoluci´on del n´umero de Mach geometr´ıa con (a) Wrap 80oy (b) Wrap 60o. En cuanto a la influencia del rake angle, los valores negativos de este ´angulo reducen los flujos secundarios, lo que se traduce en una menor generaci´on de entrop´ıa en comparaci´on con ´angulos positivos, aumentando as´ı la eficiencia. Sin embargo, estos rodetes experimentan mayores tensiones mec´anicas, lo que genera gradientes de presi´on m´as elevados. Como resultado, se produce un incremento tanto en la entalp´ıa como en la relaci´on de compresi´on. En la Figura 0.9 se muestra la generaci´on de entrop´ıa en una geometr´ıa con rake positivo y otra con rake negativo. (a) (b) Figura 0.9: Evoluci´on de la entrop´ıa en geometr´ıa con rake (a) 30oy (b) −30o. Los resultados obtenidos destacan la influencia significativa de los par´ametros tridimensionales en la relaci´on de compresi´on y la eficiencia, resaltando la importancia de su optimizaci´on. “Plantilla en LaTeX acorde con la Normativa para la elaboraci´on de informes de TFT de la ETSII (UPM)” by Javier Soto P´erez-Olivares is licensed under a Creative Commons Attribution 4.0 International License. vi Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR 5.2. DiagramadeGantt............................... 89 5.3. Objetivos de desarrollo sostenible . . . . . . . . . . . . . . . . . . . . . . . 92 Antonio Orti Hortelano xiii AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR 1. Introducci´on 1.1. Motivaci´on Debido a la creciente concienciaci´on por parte de la sociedad y de las ´elites pol´ıticas sobre los efectos de la actividad humana en el cambio clim´atico y la sostenibilidad del planeta, se han llevado a cabo grandes esfuerzos en los ´ultimos a˜nos para reducir el impacto de las actividades industriales en el medio ambiente. Dentro de estos esfuerzos destacan las investigaciones para reducir el efecto de las emisiones contribuyentes al cambio clim´atico en la generaci´on de calor y la refrigeraci´on. La bomba de calor es una de las tecnolog´ıas que se est´an potenciando para reducir el impacto medioambiental en la generaci´on de calor y fr´ıo, pues estos sistemas est´an considerados como elemento clave para la descarbonizaci´on, ya que aprovecha fuentes de energ´ıa renovables y reduce significativamente el consumo de combustibles f´osiles al transferir calor en lugar de generarlo directamente. Su alta eficiencia energ´etica permite obtener m´as calor por cada unidad de electricidad utilizada, disminuyendo las emisiones de CO2en sectores como la calefacci´on residencial, la industria y el transporte. Seg´un la EEA (European Environment Agency [2]), actualmente la demanda de calor y fr´ıo, tanto industrial como residencial, supone un 50 % del total de la energ´ıa demandada en Europa, de la cual m´as del 70 % proviene de combustibles f´osiles (gas natural). El principal factor nocivo asociado con los sistemas de bomba de calor es el uso de refrigerantes que tienen un alto potencial de calentamiento global (GWP) y contribuyen al efecto invernadero si no se manejan adecuadamente. Para mitigar estos impactos, se ha implementado el uso de refrigerantes de bajo GWP (R-32 y R-290), as´ı como la mejora en la eficiencia energ´etica de los sistemas para reducir el consumo de electricidad. Una bomba de calor toma calor de un foco fr´ıo y lo transfiere a otro m´as caliente. Para ello, sigue un ciclo de Rankine inverso compuesto de dos intercambiadores de calor, una v´alvula de laminaci´on y un compresor (Figura 1.1). Los intercambiadores de calor son un condensador y un evaporador. En el evaporador se obtiene la potencia calor´ıfica, mientras que el condensador la transmite al foco caliente. El fluido pasa por el evaporador en fase l´ıquida, y al absorber calor, pasa a fase gaseosa. El compresor eleva la presi´on hasta la presi´on del condensador, donde se extrae la energ´ıa calor´ıfica del fluido condensado que posteriormente se lamina para volver a pasar por el evaporador. Figura 1.1: Ciclo bomba de calor [3]. Puesto que se trata de un sistema que transfiere calor, su eficiencia se mide por el COP Antonio Orti Hortelano 1 1. Introducci´on (Coefficient of Performance), que es la relaci´on entre la potencia calor´ıfica que transmite y la potencia que se aporta para su funcionamiento, en este caso, al compresor. Puede tener rendimientos mayores que la unidad, lo que hace que tenga m´as rendimiento que los sistemas convencionales de generaci´on de calor y fr´ıo. Actualmente, los compresores que se utilizan para impulsar las bombas de calor son m´aquinas volum´etricas (como compresores Scroll, tornillo o alternativas). El uso de compresores centr´ıfugos en bombas de calor supondr´ıa una gran reducci´on de peso del sistema, evitar el uso de aceites,un aumento del caudal de trabajo, y por tanto, un aumento de la potencia de la bomba y del COP. Este proyecto se nutre del trabajo “Dise˜no de turbocompresores radiales con fluido no ideal para aplicaciones de bomba de calor” [4], donde se realiza un estudio de l´ınea media del compresor, utilizado como herramienta de predise˜no, y se obtienen par´ametros geom´etricos unidimensionales optimizados, que ser´an el punto de partida para definir las variables tridimensionales. Debido a la gran variedad de variables geom´etricas que hay que definir en un compresor centr´ıfugo, el an´alisis sigue una estrategia que consta de dos etapas. En primer lugar, se realiza el predise˜no mediante c´odigos de l´ınea media; posteriormente, se define la geometr´ıa tridimensional y se analiza mediante simulaci´on num´erica. Esta segunda etapa es la que se implementa en el proyecto. El fin principal es desarrollar un c´odigo de dise˜no y automatizaci´on de peque˜nos compresores centr´ıfugos para bomba de calor, de tal manera que se puedan dise˜nar diferentes geometr´ıas tridimensionales a partir de variables unidimensionales optimizadas, que, junto con el dise˜no de l´ınea media, completa el proceso de an´alisis y dise˜no de compresores centr´ıfugos, acelerando dicho proceso de dise˜no y simulaci´on. 1.2. Estado del arte En este apartado, se analizar´a la evoluci´on hist´orica de los compresores centr´ıfugos as´ı como las distintas metodolog´ıas utilizadas para su dise˜no. El desarrollo de compresores centr´ıfugos se remonta a finales de siglo XIX. En sus inicios, los compresores se utilizaban para comprimir aire en grandes instalaciones industriales y con el fin de conseguir altas relaciones de compresi´on. A comienzos de este siglo un nuevo enfoque se ha llevado a cabo, consiste en el desarrollo de peque˜nos turbocompresores centr´ıfugos para aplicaciones a peque˜na escala con fluidos no ideales u org´anicos. Con el objetivo de reducir emisiones contaminantes y contribuir a la aceleraci´on de la descarbonizaci´on, destacan l´ıneas de trabajo que optimizan el dise˜no de estas turbom´aquinas para aplicaciones de bomba de calor. Entre estas investigaciones, se encuentra la llevada a cabo por los profesores Schiffmann y Favrat, en el ´ Ecole Polytechnique F´ed´erale de Lausanne, desarrollando compresores centr´ıfugos con 20 mil´ımetros de di´ametro de salida del rodete, sin aceite, girando a velocidades de hasta 210.000 rpm y consiguiendo relaciones de compresi´on de hasta 3.3 y eficiencias superiores al 78 % (Figura 1.2). 2 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Figura 1.2: Rodete de Schiffmann [1] Otras l´ıneas de trabajo como las del profesor Casey, investigan el dise˜no de peque˜nos compresores centr´ıfugos con alto r´egimen de giro con el objetivo de sustituir compresores m´as grandes o con menor r´egimen de giro en aplicaciones como bomba de calor, utilizando refrigerantes org´anicos [5]. Asimismo, Wright et al [6] investigan en el Laboratorio Nacional de Sandia (Sandia National Laboratories) el desarrollo de ciclos de CO2supercr´ıtico para bomba de calor. Tambi´en destaca la labor del grupo de investigaci´on “Thermal Energy for Sustainability” de la Universidad Polit´ecnica de Madrid, que desarrolla, entre otros trabajos, proyectos para el dise˜no de turbomaquinaria a peque˜na escala y con fluidos no ideales desde una perspectiva anal´ıtica, de simulaci´on num´erica y experimental, y en cuyo marco se encuadra el presente proyecto. Algunos ejemplos de estos proyectos son la influencia del n´umero de Reynolds en micro-turbocompresores en condiciones fuera de dise˜no [7] y el estudio de los efectos aerot´ermicos en micro-turbocompresores [8]. Existen diferentes m´etodos para el dise˜no de turbom´aquinas, entre los que destaca el dise˜no mediante c´odigos de l´ınea media (anal´ıtico), simulaci´on num´erica y mediante resoluci´on del plano meridional (anal´ıtico). A continuaci´on se presenta brevemente cada una de ellas. C´odigos de l´ınea media Uno de los grandes problemas en el dise˜no tridimensional de los compresores centr´ıfugos es la gran cantidad de variables geom´etricas por determinar para definir y optimizar la geometr´ıa. El dise˜no mediante c´odigos de l´ınea realiza un an´alisis unidimensional, simplificando el n´umero de inc´ognitas y aproxima el valor de las variables termodin´amicas al valor en la l´ınea media que caracteriza la trayectoria que sigue el fluido por el compresor, todo ello mediante modelos te´oricos y correlaciones emp´ıricas de p´erdidas. Se trata de un an´alisis que permite obtener resultados de manera mucho m´as r´apida y menos costosa, por lo que se suele utilizar como herramienta de predise˜no de compresores centr´ıfugos. Permite analizar las variables fundamentales y caracterizar el movimiento del fluido de forma m´as sencilla y con menor costo computacional que el an´alisis mediante simulaci´on num´erica. Destacan modelos de l´ınea media como los de Meroni et al [9], que validan su modelo con Antonio Orti Hortelano 3 1. Introducci´on compresores de aire, CO2y R134a en el punto de dise˜no; profesor Casey [10], la validaci´on de su modelo se basa en an´alisis con par´ametros adimensionales y compara con mapas de rendimiento de compresores con difusor sin ´alabes, y Li et al [11], donde se utiliza como modelo de validaci´on un compresor con difusor sin ´alabes y con difusor con ´alabes, empleando aire como fluido en condiciones fuera de dise˜no. Resoluci´on del plano meridional Este m´etodo, tambi´en conocido como Throughflow Code Method, consiste en el desarrollo de c´odigos que modelan o aproximan el comportamiento fluido en un plano (soluci´on bidimensional). Calculan, por tanto, la evoluci´on de las variables termodin´amicas en el plano meridional del compresor. El fluido a trav´es del compresor centr´ıfugo es desviado en la direcci´on circunferencial por acci´on de la fuerza centr´ıfuga y el empuje de los ´alabes. Se trata de un m´etodo que permite el control de la curva del canal meridional as´ı como del ´alabe sin necesidad de entrar en la precisi´on y costo computacional de la simulaci´on num´erica. Este m´etodo emplea puntos de c´alculo fijos (q-o’s) (Figura 1.3) , situados de manera perpendicular a la direcci´on del flujo y una malla meridional que no es fija, sino que va variando durante la simulaci´on computacional para adaptarse a las l´ıneas de corriente y obtener as´ı un plano meridional que sea una imagen fiel del flujo tridimensional que atraviesa el compresor. El principal objetivo es optimizar el gradiente de velocidades en el canal meridional para tener transiciones suaves. Figura 1.3: L´ınea de corriente y q-o’s en el canal meridional de un compresor centr´ıfugo con difusor sin ´alabes [12] Destacan trabajos como los desarrollados por Yang et al [13], donde se analiza el rendimiento de compresores centr´ıfugos en condiciones fuera de dise˜no mediante este m´etodo, y el del profesor Casey [14], donde se desarrolla una metodolog´ıa para predecir el estrangulamiento del flujo y la distribuci´on de presi´on en los ´alabes. Simulaci´on num´erica (CFD) 4 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR La simulaci´on num´erica mediante t´ecnicas de fluidodin´amica computacional es el procedimiento de an´alisis y caracterizaci´on que m´as coste temporal y computacional tiene, pero tambi´en es el m´as preciso a la hora de evaluar las variables termodin´amicas del fluido en una geometr´ıa tridimensional. El proceso consiste en la definici´on de una malla computacional compuesta por vol´umenes finitos que act´uan como elementos diferenciales donde ser´an resueltas las ecuaciones de gobierno del movimiento o ecuaciones de Navier-Stokes. Cuanto m´as densa sea la malla, mayor precisi´on se obtiene en la soluci´on, pero bajo un mayor coste computacional. Se trata, por tanto, de elegir una malla adecuada que sea capaz de capturar bien las capas l´ımite y que suponga un tiempo de simulaci´on razonable. Una vez creada la malla, se imponen las condiciones de contorno del problema: Definici´on de superficies de contorno, temperatura y presi´on de entrada, caudal m´asico... as´ı como el modelo de turbulencia adecuado. Existen diferentes t´ecnicas de simulaci´on, adem´as de modelos de turbulencia. En CFD, el m´as utilizado como soluci´on de compromiso entre precisi´on y coste computacional son los modelos RANS (Reynolds-Averaged Navier-Stokes). Finalmente, se resuelven las ecuaciones de Navier-Stokes de forma iterativa hasta que se alcanzan unos criterios de convergencia y error de los residuales. Destacan investigaciones como las llevadas a cabo por Xu et al [15], donde se lleva a cabo la optimizaci´on del canal meridional y la distribuci´on de ´angulos de los ´alabes del compresor centr´ıfugo para optimizar la relaci´on de compresi´on y la eficiencia. Pakle y Jiang [16], llevan a cabo el dise˜no mediante t´ecnicas de simulaci´on num´erica de peque˜nos compresores con el objetivo de optimizar el caudal m´asico que puede pasar a trav´es del compresor y ampliar el rango de operaci´on. 1.3. Objetivos El principal objetivo del presente proyecto es automatizar el proceso de simulaci´on y dise˜no en Ansys CFX de peque˜nos turbocompresores centr´ıfugos para bomba de calor. Para ello, ser´a necesaria la consecuci´on de los siguientes objetivos secundarios: 1. Desarrollar un c´odigo en MATLAB para crear geometr´ıas en lote. Para ello se define la geometr´ıa 3D del compresor centr´ıfugo a partir de valores unidimensionales optimizados mediante dise˜no de l´ınea media. 2. Implementar la conexi´on Ansys WB-MATLAB, de manera que a partir del c´odigo de automatizaci´on en MATLAB podamos simular diferentes geometr´ıas a trav´es de un bucle de manera autom´atica, sin necesidad de abrir la interfaz gr´afica de Ansys WB. 3. Ejecuci´on del bucle integrado de simulaci´on partiendo de valores unidimensionales obtenidos de dise˜no de l´ınea media. 4. Evaluaci´on de los resultados mediante t´ecnicas de post-proceso. Antonio Orti Hortelano 5 1. Introducci´on 1.4. Estructura de la memoria La estructura de esta memoria est´a compuesta por cinco secciones. Una primera secci´on introductoria, donde se explica la motivaci´on y estado del arte de este proyecto, as´ı como los objetivos que se pretenden conseguir. En la segunda secci´on se explica la metodolog´ıa del dise˜no que se ha seguido para definir la geometr´ıa 3D de los compresores centr´ıfugos objeto de estudio y la automatizaci´on de su simulaci´on en Ansys. En particular, se define esta geometr´ıa a partir de variables unidimensionales optimizadas mediante algoritmos gen´eticos y otras funciones del trabajo de fin de m´aster “Dise˜no de turbocompresores radiales para bomba de calor”, realizado por Jaime L´azaro Abrisqueta [4]. En la tercera secci´on, se realiza un an´alisis detallado de los resultados obtenidos mediante la automatizaci´on de las simulaciones. Estas simulaciones se basan en las geometr´ıas derivadas de las variables de l´ınea media de compresores experimentales previamente validados. Estos compresores sirven como referencia para validar tanto nuestro c´odigo desarrollado en MATLAB como el proceso completo de automatizaci´on implementado. La cuarta secci´on es el apartado de conclusiones, donde se analizar´an los diferentes impactos que ha tenido el proyecto, fruto de sus resultados. Finalmente, la ´ultima secci´on donde se detalla el presupuesto, planificaci´on temporal y el impacto econ´omico, social y ambiental del proyecto. 6 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR 2. Metodolog´ıa En esta secci´on se explican los fundamentos te´oricos del an´alisis y dise˜no de l´ınea media, se explica la metodolog´ıa seguida para llevar a cabo el dise˜no tridimensional del c´odigo de dise˜no de geometr´ıas tridimensionales implementado y los fundamentos de din´amica de fluidos computacional que se han aplicado en este proyecto. Posteriormente, se detalla la parte de c´odigo de automatizaci´on integrado en MATLAB que se ha llevado a cabo para la automatizaci´on del proceso de simulaci´on y los casos de estudio implementados. 2.1. Dise˜no de l´ınea media de compresor centr´ıfugo El dise˜no de l´ınea media en compresores centr´ıfugos es fundamental para el desarrollo y optimizaci´on de estas turbom´aquinas, ya que define la trayectoria general del flujo a trav´es del sistema. Este enfoque simplifica la geometr´ıa y la din´amica del flujo, permitiendo analizar y optimizar variables clave como las caracter´ısticas geom´etricas principales y la eficiencia. La l´ınea media act´ua como una referencia esencial para relacionar el comportamiento del flujo con componentes como el rodete y el difusor. La integraci´on de correlaciones emp´ıricas y simulaciones ayuda a predecir el comportamiento del flujo, mientras que las t´ecnicas de optimizaci´on geom´etrica garantizan un dise˜no que equilibra eficiencia, estabilidad y viabilidad. Los compresores centr´ıfugos estudiados en este proyecto constan de un solo rodete y un difusor sin ´alabes. Se trata, por tanto, de compresores centr´ıfugos de un solo escalonamiento. El rodete es abierto, girando sobre un eje (rotor) y el difusor es el estator. Las velocidades que caracterizan el flujo son las siguientes: Velocidad absoluta (c): Velocidad del flujo respecto a un sistema de ejes fijos o inerciales. Velocidad relativa (w): Velocidad del flujo referida a unos ejes m´oviles no inerciales. Velocidad perif´erica (u): Velocidad tangencial con la gira el ´alabe. La velocidad perif´erica relaciona la velocidad absoluta con la relativa: c =w +u (2.1) Las velocidades absolutas y relativas se pueden descomponer en una componente tangencial, componente radial y componente axial. Los tri´angulos de velocidades de entrada y salida del rodete se muestran en las Figuras 2.1 y 2.2 : Antonio Orti Hortelano 7 2. Metodolog´ıa Figura 2.1: Tri´angulo de velocidades rodete [17]. Figura 2.2: Detalle de tri´angulo de velocidades a la salida [17]. La entrada del fluido se considerar´a axial por no existir conductos con formas ni ´alabes gu´ıa. La entrada ser´a axial en lo que se refiere a la velocidad absoluta, pero al estar el rodete girando, la velocidad relativa no ser´a axial, lo que implica tener que curvar el ´alabe a la entrada. Este ´angulo ser´a el β1. En cuanto a la salida del flujo, existen tres posibles curvaturas que se le puede dar al ´alabe definidas por el ´angulo β′ 2(´angulo geom´etrico). Se considera que β′ 2>0 si los ´alabes est´an curvados hacia atr´as, β′ 2= 0 si la velocidad relativa tiene la misma direcci´on y sentido que la meridiana y β′ 2<0 si los ´alabes est´an curvados hacia delante. Como consecuencia del giro del rodete, la masa del fluido en el interior del mismo sufre un movimiento relativo de giro respecto de los ´alabes. La superposici´on de este movimiento a la del movimiento principal del fluido lleva a una disminuci´on de la velocidad del fluido en la cara de presi´on y a un aumento en la cara de succi´on. Esto ocasiona la desviaci´on de la velocidad de salida. Por tanto, la velocidad relativa real de salida del rodete ser´a la suma de la velocidad tangencial de giro que sufre la masa fluida m´as la velocidad relativa. Es por ello que existe un ´angulo β2que es el ´angulo que caracteriza la velocidad relativa del flujo de salida. Figura 2.3: Tri´angulo de velocidades a la salida del rodete [12]. Se define el factor de deslizamiento (σ) como consecuencia del giro relativo del fluido dentro del compresor como: σ=cu2 cu2∞ (2.2) 8 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR La ecuaci´on general de una curva de B´ezier de grado nest´a dada por: B(t) = n X i=0 n i(1 −t)n−itiPi,(2.9) donde Pison los puntos de control, t∈[0,1] es el par´ametro que recorre la curva, y n i representa el coeficiente binomial. [22] Debido a la gran flexibilidad que proporcionan estas curvas, s´olo con modificar los puntos intermedios de control se podr´an generar diferentes canales meridionales en cada iteraci´on. Del vector optimizado, se conocen el valor de los di´ametros de hub,shroud, di´ametro de salida y di´ametro del difusor, as´ı como la longitud axial del rodete y el ancho de salida del rodete. Es necesario definir la regi´on del inlet al rodete, que estimamos en un 30 % de la longitud axial para darle suficiente espacio y que al simular, el solver pueda obtener resultados en un n´umero considerable de elementos de malla. Mediante la modificaci´on de los puntos de control de B´ezier respecto de las posiciones de los puntos de control inicial y final de cada curva, que son fijos pues representan los radios D1s, D1h y D2 (variables que provienen del dise˜no de l´ınea media y no se modifican) se pueden generar distintas configuraciones de canal meridional. En la Figura 2.8 se muestra un ejemplo del canal meridional del rodete definido en MATLAB y los puntos de B´ezier que definen las curvas. Figura 2.8: Ejemplo canal meridional del rodete. Distribuci´on de ´angulos En el dise˜no de l´ınea media, se considera la l´ınea que corresponde al radio que divide el flujo en dos secciones de ´area iguales, utilizando ´unicamente los ´angulos de entrada y salida del rodete en esta capa media. Sin embargo, para el dise˜no tridimensional, la curvatura Antonio Orti Hortelano 15 2. Metodolog´ıa de los ´alabes se define como una distribuci´on de ´angulos a lo largo de varias capas. Estas capas representan fracciones de la altura del ´alabe medidas desde el hub hasta el shroud, donde la capa 0 corresponde al hub y la capa 1 al shroud. La definici´on de la distribuci´on de ´angulos es fundamental en el dise˜no de la geometr´ıa. Concretamente, la distribuci´on de ´angulos en el hub es la que m´as influencia tiene sobre las variables fundamentales del compresor [23]. En cuanto a los grados de libertad, existen una gran cantidad de variables necesarias para definir los ´angulos de un ´alabe tridimensional, es por ello que se utilizan distribuciones de ´angulos en cada capa del ´alabe para definir la geometr´ıa. No ser´a necesario definir expl´ıcitamente la distribuci´on de ´angulos en todas las capas, ya que BladeGen 1cuenta con la funcionalidad General Spanwise Distribution [24] que permite generar la geometr´ıa del ´alabe definiendo ´unicamente la distribuci´on de ´angulos en el hub y el shroud (aunque se podr´ıan definir m´aximo cinco capas). La curvatura de los ´alabes se puede definir o bien con el ´angulo en coordenadas cil´ındricas θo bien con la distribuci´on β, que se relacionan como sigue: tan(β) = dθ m′(2.10) siendo m′la coordenada meridional adimensionalizada, definida como: m′=Z√dZ2+dR2 R(2.11) Las restricciones clave en el dise˜no de la curvatura de los ´alabes las imponen las variables tridimensionales ∆θwrap (Figura 2.10) , σTE yσLE (Figuras 2.11,2.9) , pues determinan cu´anto se extiende el ´alabe a lo largo del rodete, la diferencia de ´angulos entre el hub y shroud en el borde de salida y la inclinaci´on que tiene el ´alabe en el borde de ataque. El ´angulo de inclinaci´on en el borde de ataque del ´alabe (lean angle) se define como sigue: σLE =r1,tip sin(∆θ) qr2 1,tip +r2 1,hub −2r1,tipr1,hub cos(∆θ) (2.12) Dado que se est´a trabajando con microcompresores que operan a altas velocidades, la velocidad perif´erica genera tensiones muy elevadas en esta zona del ´alabe. Estas tensiones se ven incrementadas debido al ´angulo de inclinaci´on del ´alabe, lo que representa un desaf´ıo adicional. En los microcompresores, el reducido espesor del ´alabe no proporciona la rigidez necesaria para soportar dichas tensiones, a diferencia de los rodetes de mayor tama˜no, que s´ı pueden resistirlas debido a su estructura m´as robusta, por lo que en este proyecto se trabajar´a siempre con σLE = 0 [12]. 1esta herramienta se explicar´a en secciones posteriores (secci´on 2.3.2). 16 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Figura 2.9: Lean y Rake Angle[25] Figura 2.10: Wrap angle [12]. Figura 2.11: Rake angle [16]. Como punto de partida, se adoptar´an las distribuciones de ´angulos propuestas por Aungier [26], las cuales ofrecen distribuciones t´ıpicas de ´angulos y una aproximaci´on v´alida para el dise˜no inicial del ´alabe. βS=β1S+ (β2−β1S)(3ξ2−2ξ3) (2.13) βH=β1H+Aξ +Bξ2+Cξ3(2.14) βH= 90K+ (1 −K)(β2+β1H)/2 (2.15) A=−4(β2−2βH+β1H) (2.16) B= 11β2−16βH+ 5β1H(2.17) Antonio Orti Hortelano 17 2. Metodolog´ıa C=−6β2+ 8βH−2β1H(2.18) Donde βHes la distribuci´on de los ´angulos βen el hub yβSes la distribuci´on de los ´angulos βen el shroud. Ambos dependen del par´ametro adimensional ξque va de 0 a 1 y representa la longitud meridional del ´alabe en tanto por uno. La distribuci´on βen el hub depende del par´ametro βHque cambia seg´un los valores de A, B y C, que a su vez dependen de los valores de βen el hub yshroud en el borde de ataque (β1Hyβ1S respectivamente) y del valor del ´angulo en el borde de salida β2. En la Figura 2.12 se muestra c´omo es la forma de estas distribuciones t´ıpicas: Figura 2.12: Distribuciones de ´angulos t´ıpicas [26] Mediante estas ecuaciones, se representan distribuciones para el hub yshroud a partir de los ´angulos β1h,β1sβ2. Con estas distribuciones se consigue una transici´on suave entre ´angulos, evitando inconsistencias y discontinuidades que provoquen separaciones en el flujo. Para establecer los ´angulos β1hyβ1sa partir del ´angulo β1mid obtenido en el dise˜no de l´ınea media, se aplicar´a el principio de Free Vortex en el borde de ataque, mientras que el borde de salida se considerar´a constante con un ´angulo β2. De este modo, se podr´a definir la distribuci´on de ´angulos en las capas intermedias. El principio de Free Vortex [27] se utiliza en el dise˜no de compresores centr´ıfugos para establecer una distribuci´on de velocidades tangenciales y ´angulos que permita un flujo uniforme y eficiente. En este enfoque, se fuerza a que las velocidades perif´ericas en el borde de ataque (u) sean iguales a lo largo de todo el radio del ´alabe. Esto significa que en el hub (uh), en el radio medio (RMS), y en el shroud (us), la velocidad perif´erica se define como: u=ωr, donde ωes la velocidad angular de rotaci´on y rel radio correspondiente (hub,rms o shroud). Esta condici´on asegura que el momento angular del flujo sea uniforme a lo largo de todas las capas. 18 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Como consecuencia de esta distribuci´on uniforme de velocidad perif´erica, el dise˜no impone que la relaci´on tan(β)/r sea constante, donde βrepresenta el ´angulo de flujo en cada capa, pudiendo obtener as´ı β1hyβ1s. Una vez obtenidas las distribuciones de beta en hub yshroud, se modificar´an para ajustar los ´angulos ∆θwrap yσT E. De la Figura 2.10: ∆θwrap =θhub,TE −θhub,LE (2.19) mientras que el ´angulo rake en el borde de estela σT E se define como F. Merli y P. Gaetani [25]: σTE =2r2sin ∆θ 2 b2,ref ,(2.20) donde ∆θ=θhub,TE −θtip,TE. Para ajustar las distribuciones beta en hub yshroud y obtener los valores de wrap yrake buscados, primero se ajusta el valor del par´ametro K que aparece en la ecuaci´on 2.15. Al aumentar este par´ametro, la convexidad de la curva βHaumentar´a, disminuyendo el valor de θHy pudiendo controlar ∆θwrap puesto que los valores de θse han impuesto nulos en el borde de ataque (θh,LE = 0, θs,LE = 0). El valor de θh,TE se obtiene de la Ecuaci´on 2.10. Una vez controlado el wrap, se puede controlar el rake variando la distribuci´on βen el shroud atendiendo a la Ecuaci´on 2.20, de tal manera que se obtienen los puntos de B´ezier de la curva, y se ajusta su posici´on para que el valor de θs,T E tome el valor de manera que la diferencia ∆θtome el valor al que le corresponde un σT E deseado. Al igual que para controlar el wrap, se obtiene el valor de θs,T E integrando la Ecuaci´on 2.10 y sabiendo que θs,LE es cero. Las curvas βresultantes dan lugar a la primera aproximaci´on de las curvas, a partir de estas se obtienen los puntos de B´ezier de las curvas θy ser´an estas curvas las que se modifiquen de manera que los ´angulos wrap yrake est´en siempre controlados. Esto justifica el orden mostrado en el bucle de optimizaci´on del c´odigo (Figura 2.7). En la Figura 2.13, se muestra un ejemplo de las curvas θ, definidas por curvas de B´ezier, as´ı como sus puntos de control. Mediante la modificaci´on de los puntos de B´ezier de las curvas θ, se puede controlar la curvatura de la distribuci´on βy conseguir as´ı ´alabes m´as o menos curvados, estando siempre controlados el wrap yrake ´angulos, que es el motivo por el que no se modifican directamente las distribuciones β. Antonio Orti Hortelano 19 2. Metodolog´ıa Figura 2.13: Forma de la distribuci´on θen hub yshroud. Distribuci´on del espesor Para definir el espesor de los ´alabes se toman las variables t1H,t1tip yt2. La distribuci´on del espesor (Figura 2.14) se define mediante curvas de B´ezier con 3 puntos de control, siendo el punto de control intermedio el que define c´omo evoluciona el espesor a lo largo del ´alabe. Se considerar´a una distribuci´on en el shroud constante, mientras que en el hub se definir´a tomando una forma c´oncava, que representa distribuciones t´ıpicas de espesor. Se ha tomado como criterio de dise˜no la posici´on del punto intermedio de control un 30 % superior a t1H,TE. Estas formas reducen las tensiones mec´anicas en el rodete [16]. 20 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Figura 2.14: Ejemplo distribuci´on de espesor. Posici´on del ´alabe intermedio La posici´on del ´alabe intermedio es un par´ametro fundamental para la eficiencia del compresor. Se emplea como un elemento complementario a los ´alabes principales, con el objetivo de dirigir el flujo de manera eficiente sin disminuir la secci´on efectiva de entrada y reducir el deslizamiento. Esto permite evitar fen´omenos como el bloqueo s´onico y, adem´as, reduce la cantidad de material necesario para su fabricaci´on, lo que contribuye a incrementar la eficiencia global del rodete . La posici´on del ´alabe intermedio ser´a determinada a trav´es de la optimizaci´on realizada en el dise˜no de l´ınea media, representada por la variable fL,splitter, que define la fracci´on de la longitud axial del rodete donde comienza el ´alabe intermedio. 2.3. Simulaci´on num´erica basada en t´ecnicas de din´amica de fluidos computacional En esta secci´on se har´a una breve explicaci´on de los fundamentos te´oricos en los que se basan las t´ecnicas de din´amica de fluidos computacional, as´ı como de su implementaci´on para este proyecto en el software Ansys. Finalmente, se presentar´a el modelo de referencia que se ha llevado a cabo y sus configuraciones. 2.3.1. Fundamentos La simulaci´on CFD (Din´amica de Fluidos Computacional) es una herramienta clave para analizar y optimizar el comportamiento de turbom´aquinas, permitiendo estudiar el flujo Antonio Orti Hortelano 21 2. Metodolog´ıa de fluidos de forma detallada. Este enfoque se basa en la resoluci´on de las ecuaciones que gobiernan el movimiento del fluido, empleando m´etodos num´ericos y modelos adecuados para representar fen´omenos complejos como la turbulencia. Ecuaciones de Navier-Stokes El movimiento del fluido en turbom´aquinas se describe mediante las ecuaciones de NavierStokes, las cuales representan la conservaci´on de la masa, la cantidad de movimiento y la energ´ıa. Estas ecuaciones en su forma diferencial son [12]: Conservaci´on de la masa ∂ρ ∂t +∂(ρui) ∂xi =Smass i,(2.21) donde: - ρ: Densidad del fluido. - ui: Componentes de la velocidad del fluido. Conservaci´on de la cantidad de movimiento ∂(ρui) ∂t +∂(ρuiuj) ∂xj =−∂p ∂xi +∂τij ∂xj +Smom i(2.22) donde: - p: Presi´on est´atica del fluido. - τij: Tensor de tensiones viscosas. - 3 ecuaciones, una por cada componente. Conservaci´on de la energ´ıa total ∂(ρht) ∂t +∂(ρhtui) ∂xi =∂qj ∂xj +∂ ∂xj (τijui) + Senergy i(2.23) donde: - ht: Entalp´ıa total del fluido. - qj: Flujo de calor en la direcci´on j. Los t´erminos ˙ Sison t´erminos fuente en cada ecuaci´on. Las inc´ognitas son las tres componentes de la velocidad, la presi´on est´atica y la entalp´ıa total. Estas ecuaciones son no lineales y acopladas, lo que hace que su soluci´on anal´ıtica sea imposible en la mayor´ıa de los casos, especialmente en geometr´ıas complejas como las presentes en turbom´aquinas. Soluci´on num´erica de las ecuaciones de Navier-Stokes Para resolver las ecuaciones de Navier-Stokes en el flujo que atraviesa las turbom´aquinas, se emplean m´etodos num´ericos como el m´etodo de vol´umenes finitos, el m´etodo de diferencias finitas o el m´etodo de elementos finitos. En el contexto de turbom´aquinas, la simulaci´on CFD implica los siguientes pasos clave: Discretizaci´on de la geometr´ıa: La geometr´ıa del rodete, difusor y dem´as componentes se representa mediante una malla computacional, la cual puede ser estructurada o no estructurada dependiendo de la complejidad. Definici´on de condiciones de contorno: Se especifican valores en las fronteras, como el caudal de entrada, la presi´on de entrada y las condiciones peri´odicas entre ´alabes. 22 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Resoluci´on iterativa: Las ecuaciones se resuelven iterativamente hasta alcanzar una soluci´on convergente que satisfaga las ecuaciones de conservaci´on. El m´etodo de vol´umenes finitos se ha convertido en el enfoque m´as utilizado en la simulaci´on num´erica de din´amica de fluidos, principalmente debido a su capacidad para garantizar soluciones f´ısicamente consistentes. Esto se debe a que este m´etodo se basa en la formulaci´on integral de las ecuaciones de Navier-Stokes, lo que permite que las propiedades conservadas se calculen directamente sobre vol´umenes discretos, preservando los principios fundamentales de la f´ısica. Asimismo, la conservaci´on de las propiedades f´ısicas es independiente de la resoluci´on de la malla [28] [12]. M´etodo de vol´umenes finitos En el m´etodo de vol´umenes finitos se discretiza el espacio por vol´umenes de control colocados sobre una malla que sirve para delimitar los vol´umenes de control. Las ecuaciones de Navier-Stokes se integran de tal manera que se obtiene una soluci´on para cada diferencial del volumen. Existen dos tipos de mallados: Malla estruturada: Se trata de mallados que siguen una conectividad, puediendo ser identificado cada nodo por su posici´on i,j,k . En la Figura 2.15 se muestra un ejemplo de malla estructurada. Malla no estructurada: En los mallados no estructurados, los elementos de la malla no tienen un orden particular y no siguen un patr´on definido. En 2D est´an formados por elementos cuadrangulares y triangulares, mientras que en 3D est´an formados por elementos prism´aticos y piramidales. En la Figura 2.16 se muestra un ejemplo de malla no estructurada. Figura 2.15: Malla Estructurada [29] Figura 2.16: Malla No Estructurada[29]. Para convertir las ecuaciones que gobiernan al fluido en forma diferencial a forma integral se emplea el teorema de la divergencia. Las ecuaciones de Navier-Stokes en forma integral se pueden escribir como [28]: ∂ ∂t ZΩ  Q dV +IΓ F·ndS =IΓ SdV (2.24) Donde  Qes el vector de variables conservativas,  Fel vector de flujos y S, los t´erminos fuente. Ω representa el volumen de control y Γ la superficie de dicho elemento. Si, como Antonio Orti Hortelano 23 2. Metodolog´ıa se ha mencionado, se discretiza el espacio en vol´umenes de control m´as peque˜nos, se tiene que Ω=PΩi. Resolviendo la ecuaci´on 2.24 para cada volumen de control i, se aproxima el valor de las magnitudes fluidas al valor en el centro de gravedad de cada volumen de control: ∂ Q ∂t =−1 Ωi [ZΓ  F·ndS +ZΩi SdV ] (2.25) Uno de los mayores problemas en la soluci´on por vol´umenes finitos es evaluar las variables en las caras de los vol´umenes de control, cuyo valor depende de los vol´umenes que coinciden en dicha cara. Se trata de un problema num´erico que se puede resolver fundamentalmente con dos algoritmos: algoritmos upwind y algoritmos centrados [28]. Adem´as de la discretizaci´on espacial, puesto que el problema es no estacionario, se ha de llevar a cabo una discretizaci´on temporal del problema. Se discretiza el tiempo quedando la parte derecha de la ecuaci´on 2.25 Qn+1−Qn ∆tsiendo ∆tel paso temporal. Existen dos procedimientos de integraci´on [30]: M´etodo expl´ıcito: Las variables Qien el tiempo tn+1 se calculan a partir de la variable Qien el tiempo tn. M´etodo impl´ıcito: Donde las variables Qn+1 aparecen a ambos lados de la ecuaci´on y habr´a por tanto que resolver un sistema de ecuaciones. El m´etodo expl´ıcito es m´as sencillo de aplicar, aunque por problemas de estabilidad debe tener un paso de tiempo peque˜no. El m´etodo impl´ıcito es tedioso al tener que resolver un sistema de ecuaciones y tener el problema con tantas inc´ognitas; no obstante, ambas formas son eficientes. Condiciones de contorno Las condiciones de contorno en las simulaciones de turbom´aquinas determinan c´omo interact´ua el flujo dentro del dominio computacional con su entorno. Estas condiciones establecen los valores de entrada y salida necesarios para resolver las ecuaciones de flujo, permitiendo analizar el comportamiento de la turbom´aquina en distintas condiciones operativas. Generalmente, se establecen como condiciones de control la temperatura y presi´on totales, as´ı como el ´angulo de entrada del flujo. Adem´as, hay que definir el tipo de contorno de las superficies de los canales por donde pasa el flujo [12]. Estos suelen ser: Superficie de entrada: Es por donde el flujo entra a la turbom´aquina. Superficie de salida: Es la superficie por la que el fluido abandona la turbom´aquina. Pared: Superficie de los ´alabes donde se impone la condici´on de velocidad nula del fluido. 24 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR ∆y=L∆y+s80 R1/4 e ReL (2.26) donde: L: Es la longitud de la cuerda del ´alabe (blade chord). ∆y+: Es el valor objetivo especificado de y+. Rex: Es el n´umero de Reynolds basado en la distancia a lo largo de la cuerda (medido desde el borde de ataque). ReL: Es el n´umero de Reynolds basado en la longitud de la cuerda (L). Figura 2.19: ∆y[33] Five-Edge Vertex Mesh Size Reduction Este factor trata de reducir la falta de homogeneidad en los nodos de la topolog´ıa donde se unen 5 bordes de celdas. En esta opci´on, para valores menores que 1, las celdas que se unen en ese nodo se expanden para intentar mejorar la ortogonalidad de la malla, tal y como se puede ver en las Figuras 2.20 y 2.21. Figura 2.20: Factor 1 [34] Figura 2.21: Factor 0.5 [34]. Passage Tab En esta secci´on se define la distribuci´on de los elementos de la malla desde la ra´ız (hub) hasta la punta (shroud). Existen tres formas de distribuir los elementos: Element count and Size: Este m´etodo permite definir la cantidad de elementos que se distribuyen a lo largo de la superficie y su tama˜no. Antonio Orti Hortelano 31 2. Metodolog´ıa Boundary layer: Se define la distribuci´on de los elementos en tres zonas espec´ıficas: la capa l´ımite en el hub, la capa l´ımite en el shroud, y la secci´on intermedia entre ambas. Expansion rate: A trav´es de este m´etodo se puede definir el tama˜no de los elementos as´ı como la expansi´on o crecimiento que van a tener cerca de las paredes. Uniform: Se impone que todos los elementos sean del mismo tama˜no. Proportional: TurboGrid calcula autom´aticamente la cantidad de elementos, as´ı como su distribuci´on desde la ra´ız a la punta, tratando de cumplir con la relaci´on de aspecto establecida mediante el Aspect Ratio. Mediante los ratios de expansi´on (expansion rates) en la superficie de entrada como la de salida, se puede controlar c´omo se expanden las celdas de la malla. Habr´a que lograr un compromiso entre calidad de malla y n´umero de elementos ajustando su transici´on y tama˜no. Mesh statistics Para evaluar la calidad de la malla, TurboGrid proporciona estad´ısticas basadas en m´etricas espec´ıficas. Estas estad´ısticas permiten determinar el porcentaje de celdas con errores, ofreciendo as´ı una referencia clara sobre la calidad general de la malla. Las m´etricas definen l´ımites aceptables para las variables de an´alisis de la malla, ayudando a identificar problemas. Incluye las siguientes m´etricas: Element Volume Ratio: Esta m´etrica compara el volumen m´as grande de un elemento en contacto con un nodo con el volumen m´as peque˜no. Sirve como indicador de la uniformidad en la expansi´on local de los elementos. Minimum Volume: Este valor asegura que no se generen vol´umenes negativos en la malla. Representa el volumen m´as bajo de cualquier elemento que toque un nodo. Edge Length Ratio: Mide la relaci´on entre la arista m´as larga y la m´as corta de una cara. Se utiliza para evaluar el aspecto de los elementos y detectar posibles distorsiones. Maximum Face Angle:´ Angulo m´as grande entre los bordes de una cara, ´util para evaluar la inclinaci´on. Minimum Face Angle:´ Angulo m´as peque˜no en las caras, indicando posibles distorsiones. Connectivity Number: Cuenta cu´antos elementos est´an conectados a un nodo. Un valor alto puede afectar a la soluci´on computacional por falta de ortogonalidad. CFX-Pre En este m´odulo se establecen las condiciones de contorno para resolver el problema. En primer lugar, se debe definir la configuraci´on b´asica, que se muestra en la Tabla 2.4: 32 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Par´ametro Configuraci´on Tipo de m´aquina Compresor centr´ıfugo Coordenadas del marco [0,0,0] Eje de rotaci´on Z Tipo de an´alisis Steady State Tabla 2.4: Configuraci´on b´asica Posteriormente, se definen las variables f´ısicas del problema (Tabla 2.5): Fluido Tipo de fluido Datos del modelo Presi´on de referencia Transferencia de calor Turbulencia Condiciones de frontera de flujo Entrada: Presion Total Salida: Caudal m´asico Salida Caudal m´asico Entrada Presi´on total Temperatura total Direcci´on del flujo Par´ametros del solucionador Esquema de advecci´on Control de convergencia Opci´on de escala de tiempo Tabla 2.5: Par´ametros f´ısicos del modelo A continuaci´on, se definen las caracter´ısticas de cada superficie del canal, que se divide en tres dominios: Inlet,rotor R1 y difusor S1. Para ello, habr´a que definir: Inlet:En esta regi´on, se especifican las condiciones de presi´on y temperatura totales del fluido al ingresar al dominio. La superficie se considera estacionaria. Outlet:Aqu´ı se aplica como condici´on de contorno el flujo m´asico que atraviesa la salida del canal. Se trata de una superficie estacionaria. Paredes fijas: Estas superficies se representan en el modelo como rotantes pero con una velocidad opuesta al rotor para simular un estado estacionario. Adem´as, se fija velocidad del fluido en la superficie nula para tener en cuenta los efectos de la capa l´ımite. Se define un acabado con baja rugosidad as´ı como la condici´on de superficies adiab´aticas. Paredes en movimiento: Estas superficies giran a la misma velocidad que el rotor. Antonio Orti Hortelano 33 2. Metodolog´ıa Al igual que en las paredes fijas, velocidad inicial igual a cero en la superficie, un acabado con baja rugosidad as´ı como la condici´on de superficies adiab´aticas Superficies peri´odicas: En estas superficies, se asegura la conservaci´on de propiedades como la masa, el momento, la turbulencia y la transferencia de calor para el flujo que atraviesa las paredes laterales. Interfaces entre los extremos de los ´alabes y la carcasa: Estas interfaces se definen de manera similar a las laterales, con las mismas propiedades de conservaci´on. Superfice entre dominios: Esta regi´on conecta el rotor con la entrada del difusor, as´ı como el inlet con el rotor. Adem´as, se pueden definir los criterios que utilizar´a el Solver para obtener el resultado de la simulaci´on. Se podr´a definir: N´umero m´aximo de iteraciones. Error de los residuales. Criterio de convergencia. Interrupciones. Los residuales en ANSYS CFX Solver representan el error num´erico entre iteraciones sucesivas de las ecuaciones resueltas en una simulaci´on. Reflejan qu´e tan cerca est´a la soluci´on actual de cumplir con las ecuaciones de conservaci´on en cada punto de la malla. Al definir un valor m´aximo para los residuales, garantizamos que la soluci´on se mantenga dentro de un rango de error especificado. Una vez que los residuales est´an por debajo de este l´ımite, la simulaci´on se detiene. Sin embargo, puede suceder que la soluci´on oscile y no est´e realmente dentro del rango de error establecido. Para evitarlo, se implementa un criterio de convergencia que asegura la estabilidad de la soluci´on, definiendo un rango adecuado para este prop´osito. Asimismo, se pueden definir interrupciones, de manera que la simulaci´on se pare cuando estas se cumplan. Por ejemplo, se pueden introducir interrupciones de tiempo de iteraci´on o de valores m´ınimos de variables cuando queramos detener una simulaci´on cuando los valores est´en fuera de rango. Solver El Solver es el m´odulo de Ansys CFX que simula el paso del flujo por la turbom´aquina y resuelve las ecuaciones de gobierno 2.21 , 2.22 y 2.23 iterativamente en cada volumen de control de la malla. En la interfaz, se muestran diferentes gr´aficas que evolucionan seg´un se resuelven dichas ecuaciones. Se pueden ver, entre otras variables que se pueden definir, el valor de los residuales, transferencia de calor y turbulencia, como las que se muestran a continuaci´on a modo de ejemplo. Los residuales se pueden imponer en CFX Pre en el rms o como 34 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR MAX, calculando el valor global. A parte de calcular las diferentes variables, el Solver proporciona una tabla donde se muestra la calidad de la malla (Tabla 2.6): Domain Name Orthog. Angle Exp. Factor Aspect Ratio Inlet 80.8 OK 2 OK 235 OK R1 30.4 ok 8 ok 599 OK S1 57.5 OK 5 OK 836 OK Global 30.4 ok 8 ok 836 OK Percentage %! %ok %OK %! %ok %OK %! %ok %OK Inlet 0 0 100 0 0 100 0 0 100 R1 0 3 97 0 <1 100 0 0 100 S1 0 0 100 0 0 100 0 0 100 Global 0 3 97 0 <1 100 0 0 100 Tabla 2.6: Mesh Statistics Los diagn´osticos de calidad de malla eval´uan par´ametros como la ortogonalidad, la expansi´on y la relaci´on de aspecto. Estos par´ametros se clasifican como buenos, aceptables o deficientes, etiquetados como OK, ok o !, respectivamente. Se presenta el valor m´ınimo o m´aximo relevante y el porcentaje de celdas en cada categor´ıa. ´ Angulo ortogonal (Orthogonal Angle): Se refiere a los ´angulos entre las caras de las celdas. Cualquier desviaci´on significativa que implique ´angulos <20 grados puede afectar negativamente la precisi´on de la soluci´on (Figura 2.22). Factor de expansi´on (Expansion Factor): Mide las transiciones entre las celdas, asegurando que sean suaves y no existan cambios abruptos en las superficies de las celdas que puedan ocasionar alteraciones repentinas en la soluci´on (Figura 2.23). Relaci´on de aspecto (Aspect Ratio): Se refiere a la relaci´on entre el lado m´as largo y el m´as corto de la superficie de una celda. En el caso ideal, esta relaci´on es igual a uno (Figura 2.24). Antonio Orti Hortelano 35 2. Metodolog´ıa Figura 2.22: Ortogonalidad [35]. Figura 2.23: Factor de expansi´on. Figura 2.24: Relaci´on de aspecto [35]. Los criterios de calidad de la malla se representan en la Tabla 2.7: Maximum aspect ratio OK <10000 ok 10000 <100000 !>100000 Maximum expansion factor OK <5 ok 5 <20 !>20 Minimum orthogonal angle OK <50 ok 50 >20 !<20 Tabla 2.7: Criterios de calidad de la malla [36]. CFX-Post En CFX-Post se pueden tratar los datos obtenidos en el Solver, pudiendo crear figuras a partir de variables definidas en diferentes secciones o superficies del canal. Asimismo, se pueden crear tablas donde guardar valores definidos por expresiones. 36 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR 2.3.3. Modelo CFD de referencia A continuaci´on, se explica el modelo que se ha implementado en Ansys Workbench, para la simulaci´on de las geometr´ıas. Como se explic´o en el apartado 2.2, la geometr´ıa se define en el c´odigo de Matlab y las variables se escriben en un archivo .bgi de BladeGen. En la Figura 2.25 se muestra un ejemplo de una geometr´ıa dise˜nada en el modelo de referencia. (a) (b) Figura 2.25: Ejemplo de rodete tridimensional del modelo de referencia. En el modelo de referencia se va a utilizar la geometr´ıa de Schiffmann et al [37]. Se trata de una geometr´ıa con valores unidimensionales experimentales validados que se ha empleado para el dise˜no de microturbocompresores centr´ıfugos para bomba de calor. Un ejemplo de un rodete dise˜nado con el c´odigo de dise˜no de geometr´ıas de MATLAB a partir de estos valores unidimensionales se muestra en la Figura 2.25. Las condiciones de dise˜no y variables geom´etricas se muestran en la Tabla 2.8 y son los datos que se utilizar´an para generar geometr´ıas en 3D y automatizar su simulaci´on. Antonio Orti Hortelano 37 2. Metodolog´ıa Ref. Unidades Schiffmann [4] Fluido R134a T01 K 265 P01 bar 1.65 Nkrpm 180 Geometr´ıa Rodete D1,hub m 0.004 D1,tip m 0.0112 D2m 0.02 b2m 0.002 β1g,tip °-45 β1g,hub °-56 β2g°-45 L2m 0.0077 Nbf - 9 ft,split - 0.5 t1,hub mm 0.1 t1,tip mm 0.1 t2mm 0.1 Ra0.01 Geometr´ıa Difusor D3/D21.65 b3/b20.75 Tabla 2.8: Geometr´ıa Schiffmann [37]. Mallado de la geometr´ıa Para realizar el mallado, se emplear´a el m´etodo ATM, pues genera topolog´ıas de gran calidad de manera automatizada. Para definir el tama˜no de la malla, se utilizar´a el global size factor, pues permite controlar el tama˜no de la malla e ir aumentando el n´umero de elementos para conseguir mallas m´as finas o m´as gruesas. Se utilizar´a un valor objetivo de expansi´on (target expansion rate) para tratar de controlar la relaci´on de aspecto y que no haya cambios abruptos en la malla. El espesor de la capa l´ımite ser´a controlado con el m´etodo first element offset que permite establecer un valor objetivo para la altura de los elementos de la primera fila de la malla estableciendo un y+objetivo con el que TurboGrid trata de cumplir. Asimismo, para controlar el tama˜no de los elementos cerca de las paredes, se utilizar´a el m´etodo y+, que es funci´on del n´umero de Reynolds y permite controlar el espaciado ∆y entre las celdas. El control de la altura de los elementos en la direcci´on perpendicular al flujo (spanwise direction, Figura 2.26) se llevar´a a cabo mediante el factor Proportional, pues mediante un ´unico par´ametro permite hacer m´as densa o m´as gruesa la malla en esta regi´on. 38 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Figura 2.26: Spanwise direction En las regiones de inlet youtlet se controlar´a el ratio de expansi´on en la entrada y salida. Hay que prestar especial atenci´on en el outlet, donde las l´ıneas de flujo se separan y las celdas de la malla pueden ser muy alargadas en comparaci´on con su anchura, dando errores en las estad´ısticas de malla. Se trata, por tanto, de combinar estos par´ametros para obtener una malla de calidad, que capture bien las superficies l´ımite, pero sin crear una malla muy densa, pues supondr´a un costo computacional y temporal que se ver´a multiplicado ya que se est´an simulando geometr´ıas en lote. Los par´ametros que se controlar´an desde el c´odigo en MATLAB son los mostrados en la Figura 2.27. Figura 2.27: Par´ametros de la malla a controlar. Se muestran en la Tabla 2.9 los valores del n´umero de elementos obtenidos en cada dominio. Se trata de una malla de densidad media que representa valores en consonancia con la literatura [38]. Antonio Orti Hortelano 39 2. Metodolog´ıa Dominio N´umero de elementos Inlet 28.800 R1 (Rodete) 1.434.117 S1 (Difusor) 405.600 Total 1.868.517 Tabla 2.9: N´umero de elementos de la malla. Set-up A continuaci´on se definen las condiciones de contorno y variables para resolver el problema. Se definen las interfaces de la geometr´ıa explicadas en el apartado 2.3.2. El canal se divide en tres cuerpos o dominios: Inlet, R1 (rotor) y S1 (difusor). Dentro de cada dominio se definen las interfaces que se muestran en la Figura 2.28: Figura 2.28: Interfaces CFX-PRE Posteriormente, se concretan las variables f´ısicas del problema (Tabla 2.10): 40 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR turboMeshComponent1.Update(AllDependencies=True) setupComponent1.Update(AllDependencies=True) solutionComponent1.Update(AllDependencies=True) resultsComponent1.Update(AllDependencies=True) Cada l´ınea actualiza un componente diferente: turboMeshComponent1: Componente de mallado. setupComponent1: Configuraci´on del solucionador. solutionComponent1: Soluci´on del modelo. resultsComponent1: Resultados. Comando Save: El comando Save guarda el proyecto en una ubicaci´on espec´ıfica del sistema de archivos. El comando guarda la simulaci´on como un nuevo proyecto, lo que nos permite ver todos los resultados de la simulaci´on. Es decir, los resultados que se guardar´ıan si hici´eramos la simulaci´on en la interfaz gr´afica, se guardar´an en la nueva direcci´on de memoria sin alterar el proyecto base. Si no se guarda de esta forma, los resultados de la simulaci´on se reescribir´ıan. Save(FilePath="C:/SharedFolderWorkstation/...") Puesto que el archivo journal es un archivo de texto, solo bastar´a escribir una vez los comandos mencionados e ir sustituyendo las direcciones de memoria de los archivos as´ı como los valores num´ericos de los par´ametros que se deseen modificar para cada geometr´ıa dise˜nada en el c´odigo de MATLAB. Se crear´a, por tanto, un archivo de texto journal para cada geometr´ıa que se vaya a simular, pues al menos las direcciones de memoria de cada archivo deben modificarse. Finalmente, se ejecuta el archivo journal en el ejecutor de Ansys Workbench RunWB2.exe mediante el comando system() de MATLAB. C:\Program Files\ANSYS Inc\version\Framework\bin\Win64\RunWB2.exe -B -R script.wbjn Con esta instrucci´on, Workbench recibe la informaci´on en modo lote sin abrir su interfaz gr´afica, lo que facilita la automatizaci´on. Por tanto, aplicando este flujo de automatizaci´on en un bucle de MATLAB, se podr´a simular lotes de diferentes geometr´ıas. Los resultados de la simulaci´on podr´an ser posteriormente estudiados, viendo cu´al es la geometr´ıa ´optima o la influencia de un par´ametro deseado. En la Figura 2.32, se muestra un esquema de c´omo el c´odigo en MATLAB interact´ua con los m´odulos de Ansys WB. Antonio Orti Hortelano 47 2. Metodolog´ıa Figura 2.32: Esquema de control del proceso. 2.5. Casos de estudio Para finalizar el proyecto, se simular´an las geometr´ıas tridimensionales obtenidas a partir del an´alisis de l´ınea media realizado por el profesor Schiffmann [37] de compresor centr´ıfugo para bomba de calor y validar el c´odigo de dise˜no y automatizaci´on. Se pondr´a el foco en el estudio del efecto sobre la presi´on y la eficiencia que tienen diferentes canales meridionales, distribuciones de ´angulos as´ı como diferentes wrap yrake angles, que sumados a las variables unidimensionales definen la geometr´ıa tridimensional. Todo ello analizando la evoluci´on de las diferentes propiedades (entalp´ıa, entrop´ıa, n´umero de Mach y presi´on) a lo largo del proceso de compresi´on aportadas por el software de simulaci´on Ansys. Las variables geom´etricas y condiciones de control se encuentran definidas en la Tabla 2.8, que se present´o en el apartado 2.3.3 para justificar la configuraci´on realizada en el flujo de trabajo de Ansys Workbench. El proceso ser´a, por tanto, como se muestra en la Figura 2.33: 48 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Figura 2.33: Entradas caso de estudio. Con el objetivo de analizar el efecto que tienen el canal meridional, la curvatura del ´alabe (definida por las distribuciones de ´angulo) y los ´angulos wrap yrake, se han realizado simulaciones variando cada una de las variables y dejando el resto constante. Se han dise˜nado cuatro geometr´ıas con cuatro canales meridionales diferentes, pero con la misma distribuci´on de ´angulos y mismos wrap yrake angles. La evoluci´on del canal meridional de cada geometr´ıa es fruto de la modificaci´on de los puntos de control de las curvas de B´ezier que definen el canal meridional, y se modificar´an de manera que el canal aumente su curvatura al disminuir la altura de los puntos de control respecto de los radios que definen la geometr´ıa del rodete . A partir de ahora, estas geometr´ıas ser´an mencionadas como muestra la Tabla 2.14: Geometr´ıa Canal meridional GCM1 Canal meridional 1 GCM2 Canal meridional 2 GCM3 Canal meridional 3 GCM4 Canal meridional 4 Tabla 2.14: Nomenclatura geometr´ıas con diferentes canales meridionales. La segunda parte comprende la simulaci´on de geometr´ıas dise˜nadas con diferentes distribuciones de ´angulo en el hub. Para ello, se han modificado los puntos de control de la distribuci´on del ´angulo θdel ´alabe, y mediante su variaci´on se controlar´a la distribuci´on del ´angulo β. Todas ellas tienen el mismo canal meridional as´ı como iguales wrap yrake angles. A partir de ahora, se mencionar´an estas geometr´ıas como muestra la Tabla 2.15: Antonio Orti Hortelano 49 2. Metodolog´ıa Geometr´ıa Distribuci´on de ´angulo GDA1 Distribuci´on hub 1 GDA2 Distribuci´on hub 2 GDA3 Distribuci´on hub 3 GDA4 Distribuci´on hub 4 GDA5 Distribuci´on hub 5 Tabla 2.15: Nomenclatura geometr´ıas con diferentes distribuciones de ´angulo. La tercera parte consiste en la simulaci´on de diferentes geometr´ıas variando el wrap angle. Para ello, se mantienen constantes el canal meridional, la distribuci´on de ´angulos y el rake angle. A partir de ahora, se mencionar´an estas geometr´ıas como muestra la Tabla 2.16: Geometr´ıa Wrap angle GW1 Wrap 80º GW2 Wrap 70º GW3 Wrap 60 º GW4 Wrap 50 º Tabla 2.16: Nomenclatura geometr´ıas con diferentes wrap angles. Finalmente, se simular´an geometr´ıas con diferentes rake angles. Para ello, se dejar´an constantes la distribuci´on de ´angulos, wrap angle y canal meridional. Se mencionar´an como muestra la Tabla 2.17: Geometr´ıa Wrap angle GR1 Rake 30º GR2 Rake 0º GR3 Rake -30 º Tabla 2.17: Nomenclatura geometr´ıas con diferentes rake angles. Todas las geometr´ıas se han dise˜nado a partir de las variables geom´etricas unidimensionales del dise˜no de l´ınea de media de Shiffmann et al. La distribuci´on de espesores que se ha utilizado en todas las geometr´ıas es la siguiente (Figura 2.34): 50 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Figura 2.34: Distribuci´on del espesor del ´alabe Antonio Orti Hortelano 51 3. Resultados y discusi´on 3. Resultados y discusi´on El an´alisis de resultados se divide en tres partes. Se ha automatizado la generaci´on y simulaci´on de geometr´ıas con diferentes canales meridionales, distribuciones de ´angulo, wrap y rake angles para estudiar el efecto que tienen en las diferentes variables termodin´amicas durante la compresi´on, as´ı como de los dos par´ametros fundamentales que son la relaci´on de compresi´on y la eficiencia del compresor. 3.1. Efecto del canal meridional La variaci´on del canal meridional se ha realizado mediante la variaci´on de los puntos de control de las curvas del canal meridional del rodete, definidas por curvas de B´ezier, dando lugar a las curvas que se muestran en la Figura 3.1 : Figura 3.1: Distintos canales meridionales La evoluci´on del canal meridional es la que se muestra en la Figura 3.1, siendo la primera geometr´ıa que se simula la que tiene canal meridional 1, la segunda geometr´ıa tiene el canal meridional 2, la tercera el canal meridional 3 y la cuarta geometr´ıa el canal meridional 4. Se trata, por tanto, de una evoluci´on donde el canal meridional se va curvando cada vez m´as. Todas estas geometr´ıas tienen la misma distribuci´on de ´angulos, que se definen mediante puntos de control de las curvas de B´ezier que definen las distribuciones, as´ı como un wrap angle de 60ºy un rake angle de 0º. Se analizar´an la evoluci´on de la entrop´ıa, entalp´ıa, n´umero de Mach y presi´on en el compresor, y c´omo la forma del canal meridional influye sobre la relaci´on de compresi´on y la eficiencia. 52 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR A continuaci´on se muestra la evoluci´on de la presi´on (Figura 3.2) en las diferentes geometr´ıas que se han simulado. (a) (b) (c) (d) Figura 3.2: Evoluci´on de la presi´on en el canal meridional (a) GCM1, (b) GCM2 ,(c) GCM3 y (d) GCM4. Como se puede observar, conforme se va curvando el canal meridional se produce una disminuci´on de la presi´on en la salida del rodete. La evoluci´on de la presi´on es homog´enea a lo largo del canal meridional en todos los casos, lo que garantiza un buen comportamiento del compresor. Los valores de la relaci´on de compresi´on obtenidos son la relaci´on de presiones integradas y promediadas en los planos de entrada y salida del compresor, y se muestran en la Figura 3.3: Antonio Orti Hortelano 53 3. Resultados y discusi´on Figura 3.3: Valores de la relaci´on de compresi´on total a total. La evoluci´on de entrop´ıa se muestra en la Figura 3.4. (a) (b) (c) (d) Figura 3.4: Evoluci´on de la entrop´ıa en el canal meridional (a) GCM1, (b) GCM2 ,(c) GCM3 y (d) GCM4. 54 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Se puede observar una zona de mayor entrop´ıa en la parte superior debido al tip clearance que va aumentando por la mayor curvatura del canal. Asimismo, se incrementa ligeramente la entrop´ıa a la salida del rodete, ya que el fluido experimenta un cambio m´as brusco en su direcci´on. El valor de la entrop´ıa integrado y promediado en el plano de salida del rodete obtenido para cada geometr´ıa se muestra en la Figura 3.5. Figura 3.5: Entrop´ıa a la salida del rodete Antonio Orti Hortelano 55 3. Resultados y discusi´on En cuanto a la evoluci´on de la entalp´ıa, se muestra en la Figura 3.6 (a) (b) (c) (d) Figura 3.6: Evoluci´on de la entalp´ıa en el canal meridional (a) GCM1, (b) GCM2 ,(c) GCM3 y (d) GCM4. La evoluci´on de la entalp´ıa disminuye conforme se aumenta la curvatura del canal meridional, por lo que la relaci´on de compresi´on que se alcanza es menor. La transferencia de energ´ıa centr´ıfuga a energ´ıa del fluido en forma de entalp´ıa disminuye, lo que afectar´a negativamente a la eficiencia del compresor. La evoluci´on del n´umero de Mach a lo largo del canal meridional de las diferentes geometr´ıas se muestra en la Figura 3.7: 56 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Se puede apreciar la estela de desprendimiento que se va reduciendo al reducir la curvatura del ´alabe. Al tener mayor curvatura, como en la geometr´ıa GDA1, el fluido se desprende en la cara de succi´on, pues le es m´as dif´ıcil seguir el contorno del ´alabe, tal y como se puede ver en la Figura 3.13a y c´omo se va reduciendo a medida que el ´alabe reduce dicha curvatura. Los valores de entrop´ıa integrados y promediados en el plano de salida del rodete son los siguientes (Figura 3.14): Figura 3.14: Evoluci´on de la entrop´ıa a la salida del rodete No hay una tendencia clara respecto al valor promediado de la entrop´ıa en el plano de salida. Mientras que en las figuras se puede observar que los focos de alta entrop´ıa disminuyen con la disminuci´on de la curvatura del ´alabe, tambi´en se puede observar que el ´area verde (≈38,65[J/kgK]) aumenta debido a que al haber menos curvatura, la variaci´on de secci´on que experimenta el fluido es mayor entre la entrada y salida del rodete, resultando en un guiado m´as deficiente y aumento de la entrop´ıa. Asimismo, los ´alabes de mayor curvatura tienen una mayor superficie de rozamiento con las paredes, lo que implica mayores p´erdidas por fricci´on. En la zona de mayor entrop´ıa cerca de la zona de deslizamiento se puede constatar que la modificaci´on de la curvatura del ´alabe modifica los flujos secundarios, lo que implica cambios en el deslizamiento a nivel entr´opico. Antonio Orti Hortelano 63 3. Resultados y discusi´on Respecto a la evoluci´on de la entalp´ıa (Figura 3.15): (a) (b) (c) (d) (e) Figura 3.15: Evoluci´on de la entalp´ıa (a) GDA1, (b) GDA2 ,(c) GDA3 ,(d) GDA4 y (e) GDA5. 64 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR La entalp´ıa aumenta en la salida y se uniformiza, lo que supondr´a un aumento de la relaci´on de compresi´on. En la entrada del rodete, se puede observar que hay una zona de baja entalp´ıa debido a la mayor curvatura del ´alabe que obstaculiza la entrada de fluido. Conforme disminuye la curvatura, esta zona desaparece y se homogeniza, aumentando la entalp´ıa. Asimismo, existe una zona de baja entalp´ıa en el borde de ataque debido a los grandes ´angulos de incidencia que tiene el rodete. Antonio Orti Hortelano 65 3. Resultados y discusi´on La evoluci´on del n´umero de Mach se presenta a continuaci´on (Figura 3.16): (a) (b) (c) (d) (e) Figura 3.16: Evoluci´on del n´umero de Mach (a) GDA1, (b) GDA2 ,(c) GDA3 ,(d) GDA4 y (e) GDA5. 66 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR La zona de aceleraci´on a la salida del rodete, donde existen zonas supers´onicas, se va reduciendo a medida que disminuye la curvatura del ´alabe; adem´as, se puede observar c´omo en la zona de desprendimiento previamente comentada hay una zona de bajo n´umero de Mach donde se frena el flujo. El hecho de disminuir la curvatura reduce esta zona de frenado, homogeneizando el n´umero de Mach a lo largo del canal. Asimismo, se puede observar el efecto del difusor, donde disminuye de forma homog´enea el n´umero de Mach, reduciendo la velocidad del fluido y aumentando la presi´on del mismo. Figura 3.17: Eficiencia total a total En cuanto la eficiencia total a total obtenida (Figura 3.17), se puede ver que aumenta a medida que disminuye la curvatura del ´angulo, pero que en las dos ´ultimas geometr´ıas y por aumento de la zona de alta entrop´ıa, disminuye la eficiencia. Es decir, disminuir la curvatura del ´alabe ha hecho aumentar la relaci´on de compresi´on, pero llegada cierta curvatura, la eficiencia cae por un aumento de la entrop´ıa, que es m´as significativo que el aumento de entalp´ıa. En la Figura 3.18 se muestra el valor promediado de la entalp´ıa en el plano de salida del rodete. Figura 3.18: Entalp´ıa a la salida del rodete Frente al comportamiento de la entrop´ıa, hay un aumento mon´otono de la entalp´ıa que sigue la tendencia de la compresi´on. A medida que se reduce la curvatura del ´alabe, el paso entre la entrada y la salida del rodete es m´as corto y el ´area transversal cambia m´as r´apidamente, por lo que el fluido no ser´a guiado correctamente, dando lugar a flujos secundarios. En cambio, puesto que la l´ınea de corriente es m´as larga en ´alabes m´as curvados, se producen m´as p´erdidas por fricci´on en los perfiles. Antonio Orti Hortelano 67 3. Resultados y discusi´on 3.3. Efecto del wrap angle Se define wrap angle como el ´angulo θformado entre el borde de ataque y el borde de salida del ´alabe en el hub (Figura 2.10). Se trata de una variable crucial en el dise˜no del compresor, pues delimita los valores de las distribuciones de ´angulo. Para estudiar la influencia que tiene el wrap angle sobre las variables termodin´amicas, la relaci´on de compresi´on y la eficiencia, se han simulado cuatro geometr´ıas con diferentes wrap angles (80o, 70o, 60oy 50o), pero con un mismo canal meridional, igual distribuci´on de ´angulos y un rake angle de 0º. La evoluci´on de la presi´on (Figura 3.19): (a) (b) (c) (d) Figura 3.19: Evoluci´on de la presi´on (a) GW1 (80o), (b) GW2 (70o) ,(c) GW3 (60o) y (d) GW4 (50o). 68 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Conforme el canal reduce su longitud, se produce una transferencia m´as eficiente de energ´ıa centr´ıfuga a energ´ıa del fluido en forma de entalp´ıa al haber menos p´erdidas por fricci´on y turbulencia.Destacan las diferencias en la evoluci´on de la presi´on entre las geometr´ıas GW1 y GW4, donde en la primera geometr´ıa la velocidad no comienza a convertirse en presi´on hasta que el fluido no ha atravesado pr´acticamente el 50 % del canal, mientras que en la geometr´ıa GW4, antes de alcanzar la mitad del canal, ya ha transformado esa velocidad en presi´on. Por estos motivos, la presi´on a la salida del rodete aumentar´a, dando lugar a relaciones de compresi´on m´as altas, tal y como se muestra en la Figura 3.20, resultados obtenidos de la simulaci´on para la relaci´on de compresi´on (valores integrados y promediados en los planos de entrada y salida del compresor). Figura 3.20: Relaci´on de compresi´on total a total Antonio Orti Hortelano 69 3. Resultados y discusi´on La evoluci´on de la entrop´ıa para las diferentes geometr´ıas es la que se muestra a continuaci´on (Figura 3.21): (a) (b) (c) (d) Figura 3.21: Evoluci´on de la entrop´ıa (a) GW1 (80o), (b) GW2 (70o) ,(c) GW3 (60o) y (d) GW4 (50o). Se puede observar c´omo la entrop´ıa se va reduciendo conforme disminuye el wrap angle. El canal va reduciendo la longitud meridional de la l´ınea media del flujo conforme disminuye este ´angulo. En canales de mayor longitud, se produce un mayor gradiente de presi´on, dando lugar a mayores p´erdidas irreversibles por recirculaciones. Asimismo, en canales con m´as longitud se da m´as tiempo a que el fluido desarrolle efectos turbulentos y cizalladura, lo que genera un aumento de entrop´ıa debido a la disipaci´on viscosa. No obstante, en canales cortos el fluido es guiado de manera m´as deficiente, dando lugar a un aumento del ´area de alta entrop´ıa, pero siendo esta de menor intensidad que 70 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR en los casos con mayor longitud en el canal. En cuanto a la evoluci´on de la entalp´ıa (Figuras 3.22): (a) (b) (c) (d) Figura 3.22: Evoluci´on de la entalp´ıa (a) GW1 (80o), (b) GW2 (70o) ,(c) GW3 (60o) y (d) GW4 (50o). Antonio Orti Hortelano 71 3. Resultados y discusi´on Los valores de entalp´ıa y entrop´ıa integrados y promediados en el plano de salida del rodete son los mostrados en la Figura 3.23 y 3.24 respectivamente: Figura 3.23: Entalp´ıa a la salida del rodete Figura 3.24: Entrop´ıa a la salida del rodete Geometr´ıas con ´angulos wrap mayores implican ´alabes m´as largos y se producen mayores p´erdidas por fricci´on con las paredes, dando lugar a una transformaci´on menos eficiente entre energ´ıa mec´anica y energ´ıa en forma de entalp´ıa al fluido. Por tanto, a medida que se reduce el canal, se tienen menos p´erdidas por rozamiento y la entalp´ıa a la salida aumenta. 72 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR La evoluci´on de la entalp´ıa obtenida (Figura 3.31): (a) (b) (c) Figura 3.31: Evoluci´on de la entalp´ıa (a) GR1, (b) GR2 y (c) GR3. El comportamiento es similar al de la entrop´ıa, con unos mapas muy parecidos entre las diferentes geometr´ıas. Los valores obtenidos de entalp´ıa promediados en el plano de salida del rodete son los siguientes (Figura 3.32): Figura 3.32: Valor de la entalp´ıa a la salida del rodete. Antonio Orti Hortelano 79 3. Resultados y discusi´on Se obtiene menor entalp´ıa promediada en el plano de salida en la geometr´ıa con rake positivo, donde se producen mayores flujos secundarios que afectan la eficiencia de la transmisi´on de energ´ıa mec´anica al fluido en forma de entalp´ıa; no obstante, las diferencias son muy bajas. Esta disminuci´on de la entalp´ıa afectar´a tanto a la eficiencia del compresor como a la relaci´on de compresi´on. La geometr´ıa con un rake negativo es la que ofrece una mayor relaci´on de compresi´on, siguiendo una mayor entalp´ıa. Asimismo, es la que menor generaci´on de entrop´ıa presenta a lo largo del canal durante el proceso de compresi´on, debido a que con un rake negativo, se controlan mejor los flujos secundarios. Al sufrir este tipo de rodetes mayores tensiones, estas tensiones se equilibran con mayores gradientes de presi´on, obteniendo as´ı mayores relaciones de compresi´on. La eficiencia ser´a mayor en este tipo de rodetes, aunque las tensiones disminuir´an la vida ´util del rodete. Un rake nulo se presenta siempre como una soluci´on intermedia tanto para el valor de las variables termodin´amicas como para la eficiencia. Los valores de eficiencia se muestran en la siguiente Figura 3.33. Figura 3.33: Eficiencia total a total 80 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR 3.5. Influencia sobre los par´ametros fundamentales En esta secci´on se muestra el impacto de cada par´ametro tridimensional en la eficiencia total a total (ηtt) y relaci´on de compresi´on total a total (PRtt) fruto de la modificaci´on de cada uno de los par´ametros. En la Tabla 3.1 se muestran los valores obtenidos de relaci´on de compresi´on y eficiencia para cada geometr´ıa, dentro de la variaci´on que se ha hecho a cada par´ametro tridimensional. Canal Meridional PRtt ηtt GCM1 2.3679 0.749 GCM2 2.3253 0.746 GCM3 2.2578 0.738 GCM4 2.1728 0.722 Distribuci´on de ´angulos PRtt ηtt GDA1 2.1746 0.752 GDA2 2.2 0.754 GDA3 2.24 0.755 GDA4 2.3 0.754 GDA5 2.368 0.748 Wrap angle PRtt ηtt GW1 2.1434 0.715 GW2 2.2242 0.73 GW3 2.287 0.74 GW4 2.332 0.745 Rake angle PRtt ηtt GR1 2.175 0.722 GR2 2.2242 0.73 GR3 2.25 0.742 Tabla 3.1: Valores de PRtt yηtt obtenidos. En las Tablas 3.2 y 3.3 se muestra la diferencia absoluta as´ı como la variaci´on porcentual de la relaci´on de compresi´on total a total (PRtt) y la eficiencia total a total (ηtt) tomando como referencia para calcular estos valores la primera geometr´ıa simulada en cada uno de los casos (GCM1, GDA1, GW1). Los valores porcentuales y diferencias absolutas del rake angle se calculan tomando como referencia la geometr´ıa con rake nulo (GR2). Los valores negativos y positivos indican la disminuci´on o aumento de PRtt yηtt, respectivamente, en la modificaci´on de los par´ametros tridimensionales. Par´ametro tridimensional Diferencia absoluta ηtt Variaci´on porcentual ( %) Canal Meridional 0.0261 -3.48 Distribuci´on de ´angulos 0.004 0.51 Wrap angle 0.03 4 Rake angle positivo 0.0086 -1.1 Rake angle negativo 0.012 1.68 Tabla 3.2: Influencia sobre ηtt. Antonio Orti Hortelano 81 3. Resultados y discusi´on Par´ametro tridimensional Diferencia absoluta PRtt Variaci´on porcentual ( %) Canal Meridional 0.195 -8.2 Distribuci´on de ´angulos 0.193 8.16 Wrap angle 0.19 8 Rake angle positivo 0.05 -2 Rake angle negativo 0.023 1.16 Tabla 3.3: Influencia sobre PRtt. 82 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR 4. Conclusiones El presente proyecto cumple con los objetivos establecidos, desarrollando un c´odigo que automatiza la definici´on y el dise˜no de geometr´ıas tridimensionales de peque˜nos compresores centr´ıfugos, as´ı como la simulaci´on num´erica con fluidos no ideales. El c´odigo implementado en MATLAB permite definir y dise˜nar diferentes geometr´ıas a partir de las variables unidimensionales obtenidas en el dise˜no de la l´ınea media y las variables tridimensionales. Esto se logra mediante la modificaci´on de puntos de B´ezier, asegurando adem´as la automatizaci´on del proceso de simulaci´on para cada una de las geometr´ıas generadas. A continuaci´on, se presentan las conclusiones recogidas sobre el impacto del canal meridional, la curvatura de los ´alabes y los ´angulos wrap yrake. Se ha observado que, en geometr´ıas con canales meridionales m´as curvados, la entalp´ıa tiende a disminuir, lo que conlleva una reducci´on en la relaci´on de compresi´on. Adem´as, el aumento de entrop´ıa asociado a una mayor curvatura del canal meridional afecta negativamente a la eficiencia del compresor. La influencia sobre la relaci´on de compresi´on y la eficiencia es considerable; la modificaci´on del canal meridional ha supuesto una disminuci´on del 3.5 % en la eficiencia del compresor, as´ı como de m´as de un 8 % en la relaci´on de compresi´on, por lo que un buen dise˜no del canal meridional tiene una gran influencia sobre estos par´ametros. En cuanto a la influencia de las distribuciones de ´angulos, se concluye que tienen un impacto crucial en la relaci´on de compresi´on, donde se ha observado una variaci´on porcentual de m´as del 8 % en este par´ametro. A medida que disminuye la curvatura del ´angulo, la relaci´on de compresi´on aumenta debido a la reducci´on del rozamiento con las paredes, lo que favorece una transferencia m´as eficiente de energ´ıa mec´anica al fluido en forma de entalp´ıa. La capa de desprendimiento se reduce con la curvatura, aunque se ha observado un incremento de entrop´ıa a la salida en las ´ultimas geometr´ıas. Este aumento de entrop´ıa tiene un efecto m´as significativo que el de la entalp´ıa, lo que provoca una disminuci´on del rendimiento en las geometr´ıas GDA4 y GDA5; no obstante, la influencia de la curvatura de los ´alabes en la eficiencia del compresor es inferior al 1 %. El efecto del wrap angle es determinante tanto en el rendimiento como en la relaci´on de compresi´on. Un wrap angle mayor implica una mayor longitud del canal, lo que incrementa las p´erdidas por fricci´on con las paredes y reduce la eficiencia. Disminuir el valor de este ´angulo ha supuesto aumentos de eficiencia de entorno a un 4 %. Respecto a la relaci´on de compresi´on, los canales m´as cortos favorecen un cambio de presi´on m´as abrupto y una transici´on m´as eficiente entre energ´ıa centr´ıfuga y entalp´ıa, lo que genera un incremento progresivo de la relaci´on de compresi´on, aumentando entorno a un 8 %. Se ha constatado que un rake angle negativo mejora tanto la eficiencia (1.68 %) como la relaci´on de compresi´on (1.16 %), debido a la reducci´on de los flujos secundarios. Sin embargo, este tipo de configuraciones est´an sometidas a mayores tensiones, lo que reduce la vida ´util del compresor. La geometr´ıa con rake negativo presenta la menor generaci´on de entrop´ıa y el mejor rendimiento, as´ı como una mayor relaci´on de compresi´on. La geometr´ıa con rake nulo se sit´ua en un punto intermedio tanto en los valores de las variables termodin´amicas como en la relaci´on de compresi´on y la eficiencia. Antonio Orti Hortelano 83 5. L´ıneas futuras 5. L´ıneas futuras Este proyecto ha sentado las bases para la automatizaci´on del proceso de simulaci´on de compresores centr´ıfugos mediante la integraci´on de MATLAB y ANSYS. Sin embargo, existen diversas ´areas que podr´ıan explorarse para mejorar y ampliar el alcance del trabajo. En la parte de generaci´on de la malla en TurboGrid, se podr´ıa automatizar la comprobaci´on de las estad´ısticas de la malla desde el c´odigo en Matlab y establecer un bucle de manera que se modifiquen los par´ametros de la malla hasta que cumplan con los criterios de estas estad´ısticas sin necesidad de entrar en la interfaz gr´afica de Ansys y asegurando un mallado de calidad. Asimismo, ser´ıa interesante la aplicaci´on de inteligencia artificial o algoritmos para la variaci´on de los puntos de B´ezier de las diferentes curvas y, seg´un los resultados obtenidos de la simulaci´on en cada geometr´ıa, que se modifiquen hasta alcanzar la geometr´ıa ´optima. Para ello, ser´ıa necesario que el c´odigo en MATLAB pueda tomar los resultados de la simulaci´on. Se podr´ıa hacer mediante tablas exportadas a Excel. Finalmente, se podr´ıa extender la automatizaci´on al an´alisis tensional del compresor en Ansys Mechanical, realizando as´ı un estudio detallado del rodete, teniendo en cuenta las limitaciones mec´anicas. Estas propuestas representan una evoluci´on natural del trabajo realizado y su implementaci´on permitir´a avanzar hacia un dise˜no m´as eficiente y vers´atil de compresores centr´ıfugos automatizados. 84 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR Referencias [1] J. Schiffmann y D. Favrat. ((Experimental Investigation of a Direct Driven Radial Compressor for Domestic Heat Pumps)). En: International Journal of Refrigeration 32 (2009), p´ags. 1918-1928. doi:10.1016/j.ijrefrig.2009.07.006. [2] E. E. Agency. ((Heating and Cooling: Overview of EU Policies and Market Trends)). En: European Environment Agency Reports 2020/12 (2020). url:https://www. eea.europa.eu/publications/heating-and-cooling-overview. [3] F. H. Toledo-Quevedo, B. A. L´opez-Rosero, D. C. Velasco-Ni˜no y G. E. Mu˜nozCorena. ((Rendimiento de bombas de calor con ciclo de compresi´on de simple etapa y ciclo de compresi´on de doble etapa con inyecci´on de vapor para producci´on de agua caliente sanitaria)). En: Informaci´on Tecnol´ogica 33.1 (2022), p´ags. 215-222. doi: 10.4067/S0718-07642022000100215.url:https://www.scielo.cl/scielo. php?script=sci_arttext&pid=S0718-07642022000100215. [4] J. L. Abrisqueta. ((Dise˜no de Turbocompresores Radiales con Fluido no Ideal para Aplicaciones de Bomba de Calor)). Trabajo de fin de M´aster. Universidad Polit´ecnica de Madrid, 2024. [5] M. V. Casey, D. Kr¨ahenbuhl y C. Zwyssig. ((The Design of Ultra-High-Speed Miniature Centrifugal Compressors)). En: Proceedings of the 10th European Conference on Turbomachinery Fluid Dynamics & Thermodynamics (ETC10). Lappeenranta, Finland: Euroturbo, abr. de 2013, p´ags. 1-14. url:https://www.euroturbo.eu. [6] S. A. Wright, T. M. Conboy y A. J. Suo-Anttila. Summary of the Sandia Supercritical CO2Development Program. Technical Report. Albuquerque, NM, United States: Sandia National Laboratories (SNL-NM), 2011. [7] A. Sebasti´an, R. Abbas y M. Vald´es. ((Analytical prediction of Reynolds-number effects on miniaturized centrifugal compressors under off-design conditions)). En: Energy 227 (2021), p´ag. 120477. doi:10.1016/j.energy.2021.120477. [8] A. Sebasti´an, R. Abbas y M. Vald´es. ((Experimental investigation of key aerothermal phenomena in micro-scale radial turbocompressors)). En: Thermal Science and Engineering Progress 39 (2023), p´ag. 101748. issn: 2451-9049. doi:10.1016/j. tsep.2023.101748. [9] A. Meroni, B. Z¨uhlsdorf, B. Elmegaard y F. Haglind. ((Design of centrifugal compressors for heat pump systems)). En: Applied Energy 232 (2018), p´ags. 139-156. doi:10.1016/j.apenergy.2018.09.201. [10] M. Casey y C. Robinson. ((A method to estimate the performance map of a centrifugal compressor stage)). En: Journal of Turbomachinery 135.2 (2013), p´ag. 021034. doi:10.1115/1.4006586. [11] P.-Y. Li, C.-W. Gu e Y. Song. ((A new optimization method for centrifugal compressors based on 1D calculations and analyses)). En: Energies 8.5 (2015), p´ags. 4317-4334. doi:10.3390/en8054317. [12] M. Casey y C. Robinson. Radial Flow Turbocompressors: Design, Analysis, and Applications. Cambridge, UK: Cambridge University Press, 2021. [13] C. Yang, J. Du y H. Zhang. ((Time-Marching Throughflow Analysis of Centrifugal Compressors)). En: Applied Sciences 12.13 (2022), p´ag. 6576. doi:10.3390/ app12136576. Antonio Orti Hortelano 85 Referencias [14] M. Casey y C. Robinson. ((A New Streamline Curvature Throughflow Method for Radial Turbomachinery)). En: Journal of Turbomachinery 12.13 (2010), p´ag. 6576. doi:10.1115/1.3151601. [15] C. Xu y R. S. Amano. ((Meridional Considerations of the Centrifugal Compressor Development)). En: International Journal of Rotating Machinery 2012 (2012), p´ags. 1-11. doi:10.1155/2012/518381. [16] S. Pakle y K. Jiang. ((Design of a transonic centrifugal compressor for high-speed turbomachinery)). En: Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science 0.0 (2020), p´ags. 1-14. doi:10.1177/ 0954406220953171.url:https://journals.sagepub.com/home/pic. [17] M. V. del Fresno. Turbocompresores centr´ıfugos. Material acad´emico. Apuntes del Tema 6, Universidad Polit´ecnica de Madrid. 2024. [18] S. L. Dixon y C. A. Hall. Fluid Mechanics and Thermodynamics of Turbomachinery. 7th. Oxford, UK: Butterworth-Heinemann, 2014. [19] The MathWorks, Inc. MATLAB 2024: The MathWorks, Inc. Available at: https: //www.mathworks.com. Natick, Massachusetts, USA, 2024. [20] M. Mojaddam y S. A. Moussavi Torshizi. ((Design and optimization of meridional profiles for the impeller of centrifugal compressors)). En: Journal of Mechanical Science and Technology 31.10 (2017), p´ags. 4853-4861. doi:10.1007/s12206-0170933-3.url:https://link.springer.com/article/10.1007/s12206-0170933-3. [21] M. V. Casey. ((A Computational Geometry for the Blades and Internal Flow Channels of Centrifugal Compressors)). En: ASME Turbo Expo: Power for Land, Sea, and Air. ASME Paper 82-GT-153. Sulzer Brothers Ltd., Z¨urich, Switzerland. New York, NY, USA: American Society of Mechanical Engineers (ASME), 1982. [22] D. J. Ferguson. ((Interactive Interpolation and Approximation by B´ezier Polynomials)). En: SIAM Journal on Numerical Analysis 21.2 (1984), p´ags. 318-331. doi: 10.1137/0721021. [23] S. Kim, J. Park y J. Baek. ((A numerical study of the effects of blade angle distribution on the performance and loss generation of centrifugal compressor impellers)). En: Proceedings of the Institution of Mechanical Engineers, Part A: Journal of Power and Energy 226.2 (2011), p´ags. 208-216. doi:10.1177/0957650911429000. [24] ANSYS, Inc. ANSYS TurboSystem User’s Guide. Release 2021 R1. Contains proprietary and confidential information of ANSYS, Inc. and its subsidiaries and affiliates. ANSYS, Inc. Southpointe, Canonsburg, PA, USA, 2021. url:https://www. ansys.com/. [25] F. Merli y P. Gaetani. ((Numerical Investigation of the Lean Effects in Centrifugal Compressors)). En: Proceedings of 13th European Conference on Turbomachinery Fluid Dynamics and Thermodynamics (ETC13). Lausanne, Switzerland: European Turbomachinery Society, 2018, ETC2019-143. url:http://www.lfm.polimi.it. [26] R. H. Aungier. ((Centrifugal compressor stage preliminary aerodynamic design and component sizing)). En: Turbo Expo: Power for Land, Sea, and Air. Vol. 78781. American Society of Mechanical Engineers. 1995, V001T01A010. [27] H. Cohen, G. Rogers y H. Saravanamuttoo. Gas Turbine Theory. 4th. Essex, UK: Longman Group Limited, 1996. [28] J. M. T. Pulido. Turbom´aquinas. An´alisis, Dise˜no y Simulaci´on Num´erica. Formato impreso y digital, 596 p´aginas. Espa˜na: Ediciones Ingenier´ıa Energ´etica y Fluidomec´anica, 2019. isbn: 978-84-16898-95-4. 86 Escuela T´ecnica Superior de Ingenieros Industriales (UPM) AUTOMATIZACI´ ON DEL DISE˜ NO CFD TRIDIMENSIONAL DE PEQUE˜ NOS TURBOCOMPRESORES RADIALES PARA BOMBA DE CALOR [29] A. A. Ortiz. ((Caracterizaci´on y Modelaci´on Termodin´amica de un Turbocompresor Mediante CFD)). Tesis para obtener el grado de Maestro en Ingenier´ıa Mec´anica. Tesis de Maestr´ıa. Pachuca de Soto, Hidalgo, M´exico: Instituto Tecnol´ogico de Pachuca, nov. de 2015, p´ag. 127. [30] J. Oro. T´ecnicas num´ericas en ingenier´ıa de fluidos: Introducci´on a la din´amica de fluidos computacional (CFD) por el m´etodo de vol´umenes finitos. Reverte, 2012. isbn: 9788429192773. url:https://books.google.es/books?id=x9zeDwAAQBAJ. [31] I. Simulations. Turbulence Models in CFD. Accessed: 2025-01-02. 2025. url:https: //www.idealsimulations.com/resources/turbulence-models-in-cfd/. [32] F. R. Menter. ((Zonal Two Equation k−ωTurbulence Models for Aerodynamic Flows)). En: AIAA Paper. 93-2906. Orlando, FL: American Institute of Aeronautics y Astronautics (AIAA), 1993. doi:10.2514/6.1993-2906.url:https://doi. org/10.2514/6.1993-2906. [33] A. Sebasti´an y M. Vald´es. Tema 8: Introducci´on a la mec´anica de fluidos computacional (CFD) en turbom´aquinas t´ermicas. Turbom´aquinas T´ermicas - Grado Ing. Tecnolog´ıas Industriales. Email de contacto: [email protected]. 2023. [34] ANSYS, Inc. ANSYS TurboGrid User’s Guide. Release 2021 R1. ANSYS, Inc. and ANSYS Europe, Ltd are ISO 9001:2015 certified. ANSYS, Inc. Southpointe, 2600 ANSYS Drive, Canonsburg, PA 15317, USA, ene. de 2021. url:http: / / www . ansys.com. [35] ANSYS, Inc. Appendix A: Mesh Quality, ANSYS Meshing Application Introduction. Inventory #002645. ANSYS, Inc. April 28, 2009, 2009. url:https://www.ansys. com. [36] I. ANSYS. ANSYS CFX-Solver Theory Guide. 2021 R1. ISO 9001:2015 certified. ANSYS, Inc. Canonsburg, PA, USA, ene. de 2021. url:http://www.ansys.com. [37] J. Schiffmann y D. Favrat. ((Experimental investigation of a direct driven radial compressor for domestic heat pumps)). En: International Journal of Refrigeration 32.7 (2009). Published online 25 July 2009, p´ags. 1918-1928. doi:10.1016 /j . ijrefrig . 2009 . 07 . 008.url:https : / / www . sciencedirect . com / science / article/pii/S0140700709001578. [38] O. F. Atac, J.-E. Yun y T. Noh. ((Aerodynamic Design Optimization of a Micro Radial Compressor of a Turbocharger)). En: Energies 11.7 (2018), p´ag. 1827. doi: 10.3390/en11071827. Antonio Orti Hortelano 87 Referencias ANEXOS A. Planificaci´on temporal El proyecto de fin de grado ha tenido una duraci´on de seis meses, comenzando en julio de 2024 con las reuniones previas y definici´on de objetivos. Durante ese mes se realiz´o un estudio de los trabajos previos relacionados con el presente TFG. A continuaci´on, se presenta la planificaci´on temporal del proyecto. Figura 5.1: Planificaci´on temporal 88 Escuela T´ecnica Superior de Ingenieros Industriales (UPM)