Una contribución al análisis de problemas de flujo compresible sin utilización de viscosidad artificial

Revista Internacional de Métodos Numéricos para Cálculo y Diseño en Ingeniería. Vol. 10,3, 187207 (1994) UNA CONTRIBUCIÓN AL ANÁLISIS DE PROBLEMAS DE FLUJO COMPRESIBLE SIN UTILIZACI~N DE VISCOCIDAD ARTIFICIAL FERNANDO QUINTANA División Mecánica Computacional Centro Atómico Bariloche. CNEA. Argentina. RESUMEN En este trabajo se presentan unas técnicas desarrolladas con el objeto de obtener un algoritmo para la solución de problemas de flujo compresible sin la aplicación de términos de viscosidad artificial. Esto es motivado por el hecho de que si bien la calidad de la solución obtenida mediante procedimientos conocidos es muy buena -en particular los choques y las capas límite son capturados de manera muy precisa y sin oscilaciones-, en muchas ocasiones es preciso recurrir a la pericia del usuario para ajustar unos coeficientes de difusión artificial para preservar la monotonicidad de la solución, en particular cuando en la misma aparecen choques como es el caso de problemas transónicos y supersónicos. Las ideas propuestas se basan en aumentar la viscosidad numérica propia de un algoritmo conocido (Taylor Galerkin), teniendo en cuenta las velocidades de transporte de cada ecuación. Este tipo de técnica resulta adecuada para problemas subsónicos, pero en problemas supersónicos e hipersónicos, las discontinuidades (choques) que aparecen en la solución son muy severas y hacen que la viscosidad numérica sea insuficiente para eliminar las oscilaciones. Se resuelven algunos ejemplos de flujo subsónico para un perfil NACA0012 mediante las ecuaciones de Euler. Se estudia el mismo problema a número de Mach creciente hasta llegar a un problema transónico con aparición de un choque. SUMMARY The present study is devoted to the development of an algorithm for the solution of compressible flow problems without artificial viscosity. The study is motivated by the fact that the solution obtained using known algorithms is accurate and free of ascillations close to shocks and boundary layers, but some tuning is necessary for the artificial viscosity coefficients if monotonicity is to be preserved when the solution presents discontinuities. The algorithm used is one-step Taylor-Galerkin. The proposed ideas are based on increasing the numerical viscosity according with the transport velocities. Although the techniques here presented work for subsonic problems, in supersonic and hypersonic cases, the strong discontinuities that the solution may present induce oscillations that are not supressed with the numerical viscosity that the algorithm introduces. Some subsonic examples for a NACA0012 are studied, and to analyze the behaviour of the algorithm, the Mach number is increased for the same geometry, until a shock is obtained in a transonic problem. Recibido: Septiembre 1993 OUniversitat Politbcnica de Catalunya (España) ISSN 0213-1315 188 F. QUINTANA A. ECUACIONES DIFERENCIALES DEL PROBLEMA Se trata de modelar problemas de flujo compresible sin tener en cuenta los términos viscosos. Para ello se utiliza el sistema de ecuaciones de conservación, lineales, de primer orden y acopladas, conocido como ecuaciones de Euler, que escritas en forma conservativa se escriben: donde Sij es la delta de Kronecker, E es la energía total por unidad de volumen (cinética más interna), las u, son las componentes cartesianas del vector velocidad V y p es la presión, que para el caso del aire, se calcula en base a la ecuación de estado para un gas perfecto: donde R es una constante del gas, que está relacionada con los calores específicos por medio de: y T es la temperatura T = w, con lo cual, la ecuación de estado (2) puede escribirse: B. IMPORTANCIA DEL PASO TEMPORAL La extensión a elementos finitos del esquema de Lax-Wendroff en diferencias finitas, es conocida como algoritmo de Taylor-GalerkinlO. En su versión de un paso, el método permite escribir el sistema discreto de ecuaciones como: Puede observarse, tal como lo hiciera notar el autor de este algoritmo2, que el término disipativo que introduce el mismo está vinculado a.1 incremento de tiempo At. La ecuación (5) puede escribirse nuevamente como: Se observa aquí que con el incremento de tiempo At varía linealmente la magnitud del término disipativo. La idea consiste en incrementar este término y una forma de hacerlo es a través del incremento de tiempo: PROBLEMAS DE FLUJO COMPRESIBLE SIN VISCOCIDAD ARTIFICIAL Se han colocado los supraíndices n y b para destacar que no necesariamente se trata del mismo At. En efecto: si se incrementa Atb para aumentar el término disipativo, dado el carácter explícito del algoritmo, será también necesario disminuir Ata para preservar la estabilidad del. esquema. El siguiente análisis permite determinar el valor de Ata que garantiza estabilidad para el caso de incrementos de tiempo diferentes. Según se demostrara en1', para el caso unidimensional lineal el algoritmo de LaxWendroff (LW) en diferencias finitas, resulta totalmente equivalente a Taylor Garlerkin (TG). Se analiza aquí la estabilidad del algoritmo de Lax-Wendroff siguiendo el rnét,odo de Von Neumann o de Fourier en el caso en que se escoge un incremento de tiempo mayor para el término que aplica la disipación numérica. Siendo 4 una variable escalar, para el caso unidimensional, en una malla uniforme de tamaño elemental Ax, la ecuación tipo que resulta para un nodo cualquiera de la malla luego de efectuar las integraciones correspondientes es: Donde C = e es el número de Courant, el subíndice indica número de nod6.y el supraíndice paso temporal. Si se utiliza un incremento de tiempo a veces mayor para el segundo término: Se inserta un modo de Fourier: Para prevenir el crecimiento no acotado de la solución, debe verificarse que1 Al 5 1, que es equivalente a: Con lo que la condición de estabilidad resulta en una limitación del número de Courant, y por tanto en el paso de tiempo: 190 F. QUINTANA En el caso de trabajar con matriz de masas consistente, la condición resultante es: Se propone como caso de prueba, resolver el flujo no viscoso alrededor de un perfil NACA0012 con un número de Mach igual a 0.5 en el flujo no perturbado. El mismo se resuelve usando la técnica de TG de dos pasos con viscosidad artificial, la cual será usada como referencia. A continuación, se resuelve el mismo problema introduciendo las modificaciones propuestas en el Apartado B, y se comparan los resultados. La malla utilizada en ambos casos tiene 2550 elementos y 1226 nodos y fue generada con un generador automático que utiliza la técnica de avance frontal. En las Figuras 1 a 3 se pueden apreciar los resultados obtenidos utilizando la técnica de TG de dos pasos con viscosidad artificial. $as Figuras 1 y 2 corresponden a las isolíneas de presión alrededor del perfil NACA0012 y un detalle del aspecto de las mismas en proximidades del punto de estancamiento. Para la interpretación de los valores numéricos debe tenerse en cuenta que se ha efectuado una normalización (adimensionalización) de las ecuaciones utilizando como variables adimensionales las siguientes: Donde el subíndice r indica valores de referencia. La barra destaca que se trata de magnitudes adimensionales, aunque en adelante y para simplificar la notación se omitirá el uso de la misma. En este sentido, la presión adimensional calculada en el flujo no perturbado es p, = 2.857. Este valor surge de considerar que la velocidad y la densidad en el flujo no perturbado tienen valor unitario. Como la velocidad del sonido es e, = E, de aquí puede despejarse p, en función de magnitudes conocidas como el número de Mach y el cociente de calores específicos y. En la Figura 3 se observa la variación de la presión a lo largo de la línea de - corriente correspondiente al punto de estancamiento, en los nodos próximos al mismo. l PROBLEMAS DE FLUJO COMPRESIBLE SIN VISCOCIDAD ARTIFICIAL ,191 Figura l. Isolíneas de presión para Mm = 0.5 Figura 2. Detalle de las isolíneas de presión en la zona próxima al punto de estancamiento. 192 F. QUINTANA Figura 3. Variación de la presión a lo largo de la línea de corriente de estancamiento Figura 4. Detalle de las isolíneas de densidad en la zona cercana al punto de estancamiento. PROBLEMAS DE FLUJO COMPRESIBLE SIN VISCOCIDAD ARTIFICIAL 193 La Figura 4 muestra las isolíneas de densidad en la zona próxima a dicho punto. Debe tenerse en cuenta que la densidad para el flujo no perturbado es p, = 1. La Figura 5 representa el comportamiento de la densidad a lo largo de la línea de corriente del punto de estancamiento, en los nodos situados cerca del extremo anterior 'del perfil. En las siguientes Figuras (6 a 10) se presentan los mismos resultados obtenidos utilizando la técnica presentada anteriormente. Asimismo se ha indicado el valor "teórico". Estos valores fueron obtenidos en base a las relaciones válidas para un gas ideal que sufre una compresión isoentrópicat y considerando que en el punto de parada la velocidad es cero. Analizando la Figura 6 puede verse que el campo de presiones es coincidente de forma general con el obtenido con el método usado como referencia por lo que la solución puede considerarse como buena. Figura 5. Variación de la densidad a lo largo de la línea de corriente del punto de estancamiento. En la Figura 7 el detalle de las isolíneas de presión demuestra que no existen comportamientos oscilatorios ni picos espúreos y lo mismo puede decirse de la Figura 8 que muestra una variación suave de la presión a lo largo de la línea de corriente que pasa por el punto de parada. Sin embargo, puede verse en las Figuras 9 y 10 que la densidad sufre un pequeño salto (overshoot) localizado en dos o tres elementos vecinos. Los autores han podido comprobar que este comportamiento persiste en caso que se haga un remallado y se utilice en esta región un tamaño de elemento más pequeño y en todos los casos estará limitado a los dos o tres elementos que tienen en común el nodo situado exactamente en el punto de parada. Si se modifica la malla de manera que no haya un nodo coincidente con dicho punto, el comportamiento de las isolíneas de densidad'resulta más suave. En ningún caso estas oscilaciones ponen en Peligro la estabilidad del método propuesto. t Ver Apartado 2.10 de la referencia 3. 194 F. QUINTANA En este Apartado se ha mostrado que en base a las consideraciones hechas anteriormente es posible obtener soluciones sin la adición de viscosidad artificial para problemas que no presenten choques. Se ha mostrado que la calidad de las soluciones es muy buena y solo persiste un problema de oscilaciones en una región muy limitada. Mediante el siguiente estudio, se muestra que las velocidades con que se transportan cada una de las ecuaciones son diferentes, se deducen unos incrementos de tiempo óptimos, se establece un criterio para elegir el tamaño elemental y se presentan los correspondientes ejemplos numéricos. Figiira 6. Isolíneas de presión usando TG de 1111 paso sin viscosidad art.ific.iii1 Figura 7. Detalle de las isolíneas de presión usando la técnica de TG de un paso sin viscosidad artificial. PROBLEMAS DE FLUJO COMPRESIBLE SIN VISCOCIDAD ARTIFICIAL Figura 8. Variación de la presión a lo largo de la línea de corriente del punto de estancamiento. Figura 9. Detalle de las isolíneas de densidad en la zona próxima al punto de estancamiento. F. QUINTANA Figura 16. Comparación de las isolíneas de presión y densidad para el caso M, = 0.8 versus la solución por Taylor Galerkin usando viscosidad artificial. Figura 17. Variación de los residuos (norma Lp) para el caso M, = 0.8. PROBLEMAS DE FLUJO COMPRESIBLE SIN VISCOCIDAD ARTIFICIAL Figura 18. Variación de la presión alrededor del perfil para el caso M, = 0.8. lX 1 DENSIDAD Figura 19. Variación de la densidad alrededor del perfil para el caso M, = 0.8. 208 F. QUINTANA AGRADECIMIENTOS El autor agradece a los profesores Eugenio Oñate y Juan Miquel Canet de la Universidad politécnica de Cataluña por su asesoriariiiento durante la realización de la tesis doctoraliO~ll de la cual se han extraído los resultados de este artículo. Se agradece también la ayuda prestada por el Centro Internacional de Métodos Numéricos en Ingeniería ,para la realizar este trabajo. REFERENCIAS 1. F.P. Brueckner & J.C. Heinrich, "Petrov Galerkin finite element model for compressible flows", Department of Aerospace and Mechanical Engineering, University of Arizona, Tucson, Arizona 85721, (1990). 2. J. Donea , "A Taylor-Galerkin Method for the convective transport problems," Int. Journal .Num. Methods in Engineering, Vol. 20, pp. 101-119, (1984). 3. H: Liepmiinn & A. Roshko, "Elements of Gasdyñamics,)' Galcit Aeronautical Series, John Wiley & Sons Ed, (1957). 4. J. Miquel, E. Oñate, & F. Quintana, "A Formulation for the Solution of Euler Equations for Compressible Flow using Finite Elements. ", Proceedings of the 5th Int. Symp. on Num. Methods in Engineering, Springer Verlag Ed.(1989), Lausanne, Switzerland, (1989). 6. H. Liepmann & A. Roshko, "Elements of Gasdynamics," Galcit Aeronautical Series, John Wiley & Sons Ed. (1957). 7. D. Potter, Computational Physics, John Wiley & Sons, (1980), (First Edition 1972). 8. W. J. Usab & E. Murman, "Embedded Mesh Solutions of the Euler Equations Using a Multiple Grid Method", en Advances in Computational Transonics, (W. G. Habashi, Ed). Pineridge Press, pp. 447-472, (UK), (1985). 9. O.C. Zienkiewicz & R.L. Taylor, "The Finite Element Method", Fourth edition, Mc. GrawHill, (1989). 10. J.F. Quintana, Análisis de problemas de flujo compresible de advección dominante por el método de elementos finitos, Tesis Doctoral, E.T.S.E.C.C.P.B., Universidad Politécnica de Cataluña, España, Marzo, (1993). 11. F. Quintana, "Análisis de problemas de flujo compresible de advección dominante por el método de elementos finitos", Tesis doctoral Universidad Politécnica de Cataluña, Marzo (1993). 12. E. Oñate, F. Quintana, R. Codina y J. Miquel, "Finite element formulations for compressible and incompressible fluids", Publicación de Investigación N" 14, CIMNE, Barcelona, (1991).