Simulación bidimensional del flujo de convección libre en una cavidad cuadrada utilizando el método de la ecuación de Boltzmann en redes¹

Two-Dimensional Free Convection Flow Simulations on a Square Cavity Through the Use of Lattice Boltzmann Equation Method²

Simulação bidimensional do fluxo de convecção livre em uma cavidade quadrada utilizando o método da equação de Boltzmann em redes³

Elkin Flórez-Serrano⁴
Raquel Laguado Ramírez⁵

¹ Este artículo se deriva de un proyecto de investigación denominado Desarrollo de un software libre para simulación bidimensional de flujos de fluido, de los grupos de Investigación, Gestión, Administración de Producción y Operaciones (INGAPO) y de Ingeniería Mecánica de la Universidad de Pamplona (GIMUP) y fue financiado por la Universidad de Pamplona, Norte de Santander, Colombia.
²Submitted on: January 6, 2011. Accepted on: May 2, 2011. This article results from the research project Development of Free Software for Bidimensional Simulation of Fluid Flows, developed by the research groups INGAPO and GIMUP and financed by the Universidad de Pamplona, Norte de Santander, Colombia.
³Data de recepção: 6 de janeiro de 2011. Data de aceitação: 2 de maio de 2011. Este artigo se deriva de um projeto de pesquisa denominado Desenvolvimento de um software livre para simulação bidimensional de fluxos de fluido, dos grupos Pesquisa, Gestão, Administração de Produção e Operações (INGAPO) e o grupo de Pesquisa em Engenharia Mecânica da Universidade de Pamplona (GIMUP) e foi financiado pela Universidade de Pamplona, Norte de Santander, Colômbia.
⁴Ingeniero mecánico, Universidad Francisco de Paula Santander, Cúcuta, Colombia. Magíster en Ingeniería Mecánica, Universidad de los Andes, Bogotá, Colombia. Magíster en Ingeniería Química y de Procesos, Universidad Rovira i Virgili, España. Doctor en Ingeniería Mecánica, Fluidos y Aeronáutica, Universidad Politécnica de Cataluña, España. Profesor asociado Programa de Ingeniería Mecánica, Universidad de Pamplona, Norte de Santander, Colombia. Correo electrónico: eflorez@unipamplona.edu.co. ]]> ⁵Ingeniera Industrial, Universidad de Santander, Bucaramanga, Colombia. Magíster en Organización Industrial, Universidad Rovira i Virgili, España. Profesora de Ingeniería Industrial, Universidad de Pamplona, Norte de Santander, Colombia. Correo electrónico: raquell@unipamplona.edu.co.

Fecha de recepción: 6 de enero de 2011. Fecha de aceptación: 2 de mayo de 2011.

Resumen

La convección es un fenómeno fundamental encontrado en varias aplicaciones atmosféricas e industriales. En el presente artículo se presenta la simulación numérica del flujo en una cavidad, calentada en una de sus paredes, en régimen laminar, números de Rayleigh (Ra) de 10³ ≤ Ra ≤ 10⁶ y para un número de Prandtl de Pr = 0,7, mediante el método de la ecuación de Boltzmann en redes (LBE). Esta técnica ha demostrado ser muy eficiente y poderosa en la dinámica de fluidos computacional (CFD). La velocidad del flujo se calcula a partir del uso de la tradicional función densidad de distribución en el modelo D2Q9 (de dos dimensiones y nueve velocidades) y la temperatura se obtiene a partir de la inclusión de una función de densidad de distribución de energía interna en un modelo D2T5 (dos dimensiones y cinco temperaturas). El nuevo modelo térmico utilizado demostró ser estable y los resultados obtenidos demuestran una gran exactitud al ser comparados con los resultados obtenidos experimental y numéricamente, mediante otros métodos de CFD.

Palabras clave: Simulación numérica, ecuación de Boltzmann en redes, convección natural.

Abstract

Convection is a key phenomenon found in different atmospheric and industrial applications. This work presents the numerical simulation of flow in a cavity which has been heated on one wall, in laminar regime, with Rayleigh numbers (Ra) of 103 ≤ Ra ≤ 106 and a Prandtl number of Pr = 0.7, through the lattice Boltzmann equation method -LBE. This technique has proven to be very efficient and powerful in computational fluid dynamics -CFD- (Y. W. Know, 2006). The flow velocity is calculated based on the use of density distribution function in the model D2Q9 -, with two dimensions and nine velocities. Temperature readings are obtained through the model D2T5 -, with two dimensions and five temperatures. The new thermal model used showed to be stable, and the results are highly accurate compared to the experimental and numerical results obtained through other CFD methods.

Key words: Numerical simulation, Lattice Boltzmann equation, natural convection.

Resumo

A convecção é um fenômeno fundamental encontrado em varias aplicações atmosféricas e industriais. En o presente artigo se apresenta a simulação numérica do fluxo em uma cavidade, aquecida em uma de suas paredes, em regime laminar, números de Rayleigh (Ra) de 10³ ≤ Ra ≤10⁶ e para um número de Prandtl de Pr = 0,7, mediante o método da equação de Boltzmann em redes (LBE). Esta técnica demonstrou ser muito eficiente e poderosa na dinâmica de fluidos computacional (CFD). A velocidade do fluxo é calculada a partir do uso da tradicional função densidade de distribuição no modelo D2Q9 (de duas dimensões e nove velocidades) e a temperatura é obtida a partir da inclusão de uma função de densidade de distribuição de energia interna em um modelo D2T5 (duas dimensões e cinco temperaturas). O novo modelo térmico utilizado demonstrou ser estável e os resultados obtidos demonstram uma grande exatidão quando comparados com os resultados obidos experimental e numericamente, mediante outros métodos de CFD.

Introducción

Desde la última década del siglo pasado, una gran cantidad de problemas que involucran diferentes aplicaciones de flujos de fluidos se han simulado utilizando el método de la ecuación de Boltzmann en redes (LBEM), por mencionar algunos (Quian et ál., 1992; Filippova y Hannel, 2000; Dazhi et ál., 2003; Flórez et ál., 2008). Para los flujos isotérmicos, el método ha demostrado ser de gran exactitud, estable y "económico" computacionalmente comparado con los métodos clásicos, utilizados en la mecánica de fluidos computacional (CFD).

La LBEM representa la mínima forma de la ecuación cinética de Boltzmann (Higuera y Jiménez, 1989), y el resultado es una elegante y simple ecuación, para una función de distribución de densidad discreta . Esta función de distribución de densidad representa la probabilidad de encontrar una partícula en un sitio de la red , en un instante t y a una velocidad c_i.

Originalmente, solo fueron consideradas la conservación de masa y momento en los modelos LBEM. No obstante, en muchas aplicaciones los efectos térmicos en flujos de fluidos son importantes y algunas veces crítico. Por lo tanto, se han desarrollado algunos modelos que permiten incluir dichos efectos. En general, la inclusión de los efectos térmicos en los LBEM encajan en tres modelos: el de multivelocidad, el de un escalar pasivo y el de doble función de distribución.

El modelo de multivelocidad consiste en aumentar o expandir la función de distribución, con el fin de obtener la temperatura macroscópica (Chen et ál., 1994). Sin embargo, este requiere un gran esfuerzo computacional, debido a la adición de velocidades discretas en cada nodo de la red y sufre de inestabilidades (Kuznik et ál., 2007). El modelo de un escalar pasivo consiste en resolver el campo de velocidad mediante el LBEM y la ecuación de temperatura macroscópica independientemente.

La ecuación para la temperatura macroscópica es similar a una ecuación de evolución de un escalar pasivo, sin tener en cuenta la disipación de calor viscoso y el trabajo de compresión realizado por la presión. El acople de la ecuación de la temperatura con el LBEM se realiza por medio de un potencial que se adhiere a la ecuación de la función de distribución. Este modelo presenta una mayor estabilidad numérica, comparado con el modelo multivelocidad (D'Orazio et ál., 2004); pero tiene la desventaja de no permitir incluir en las soluciones la disipación de calor viscosa y el trabajo de compresión hecho por la presión; además, al incluir una ecuación para la temperatura, se pierde la simplicidad del LBEM.

El tercer modelo, el de doble función de distribución, utilizado en el presente trabajo, se basa en el trabajo de He et ál. (1998). En este modelo se introduce una función de distribución de energía interna para simular el campo de temperatura, la cual es análoga a la función de distribución de densidad para simular el campo de velocidad. La estabilidad y la precisión del modelo de doble función de distribución han sido verificadas en sendos estudios (He et ál., 1998; Guo et ál., 2002; Kuznik et ál., 2004).

La simulación del problema de la convección natural en una cavidad cuadrada se trabajó en régimen laminar, números de Ra ≤ 10⁶, y para Pr = 0,71. Para la incorporación de la fuerza de cuerpo, existente en los problemas de convección natural, se utiliza la aproximación de Boussinesq.

1. Modelo 2D para simular el flujo térmico incompresible

En la Figura 1, T_h y T_c describen la temperatura de la pared caliente y fría, respectivamente; las paredes horizontales no tienen pérdida de calor (paredes adiabáticas). u y v describen la velocidad horizontal y vertical del flujo en las fronteras, respectivamente. Entre tanto, g describe la aceleración de la gravedad. Las principales hipótesis utilizadas en el modelo son:

• El término de colisión en la ecuación de Boltzmann se expresa en función de un simple tiempo de relajación para el equilibrio local (Bhatnagar et ál., 1954).
• El número de Knudsen se asume pequeño (teoría del continuo).
• El flujo es incompresible.

De esta forma, la función de distribución, tanto para la función de distribución de densidad f como para la de energía g, obedecen a su respectiva ecuación de transporte de Boltzmann, definidas por:

(1)

(2)

Para el campo de temperatura.

Donde, Π_v y Π_t son los tiempos de relajación de la función de distribución de la densidad y de la energía interna, respectivamente; F_i es la cantidad de movimiento debida a la fuerza de cuerpo o boyante. Para obtener las ecuaciones (1) y (2) de la ecuación de Boltzmann, la velocidad debe ser discretizada: durante un dt, la función de distribución se mueve a lo largo de una dirección discreta de la red dx_i = e_idt. En las simulaciones realizadas se utilizó un modelo de dos dimensiones y nueve direcciones discretas (D2Q9) para la función de distribución de densidad, y otro de dos dimensiones y cinco direcciones discretas (D2T5) para la energía interna (Figura 2).

Las funciones de distribución de equilibrio para la densidad y para la energía interna están dadas por las expresiones (Guo et ál., 2002):

(3)

(4)

Respectivamente. Las variables macroscópicas, densidad ρ, velocidad u y temperatura θ son definidas por: ρ = Σ f_i= u Σ f_i e_iα/ρ, θ = Σ g_i, y donde i y α representan las coordenadas cartesianas. Utilizando la expansión de Chapman-Enskog, la ecuación de continuidad y las ecuaciones de Navier-Stokes pueden obtenerse con una aproximación de segundo orden de la ecuación (1), como se muestra en (Quian y Orszag, 1993). Similarmente, la ecuación de transporte para la temperatura se puede obtener de la ecuación (2), como se muestra en (Guo et ál., 2002).

La fuerza de cuerpo se obtiene a partir de la aproximación de Boussinesq, donde todas las propiedades del fluido son consideradas constantes, excepto el término de la fuerza de cuerpo. De esta forma, la densidad del fluido queda definida por ρ = ρ_m [1 - ß (T-T_m)}, donde ρ_m es la densidad promedio del fluido, T_m es la temperatura promedio del fluido y ß es el coeficiente de expansión térmico. Así, la fuerza de cuerpo de la ecuación (2) queda definida por:

(5)

Donde ω es el factor de peso, que debe determinars para lograr la isotropía del tensor de cuarto orden de las velocidades y la invariancia Galileana (Quian et ál., 1992), g es la aceleración de la gravedad y c_s = l/√3 es la velocidad del sonido. La viscosidad para simular la convección natural es necesario definir, a partir de los parámetros del problema, una velocidad característica apropiada (V_c = √ß gΔθ H ), donde H es la longitud característica de la cavidad, para mantener el flujo dentro del régimen incompresible. Además, se definen los números de Prandtl y Rayleigh, que permiten contar con dos ecuaciones más en el modelo: una para la viscosidad cinemática (v = √V ²H ² Pr/Ra) y otra para la difusividad térmica (α = v/Pr). En la simulación se utilizaron condiciones de Dirichlet para paredes con temperatura fija y condiciones de Newman para las paredes adiabáticas. Los criterios de convergencia para todos los casos estudiados son:

(6)

2. Resultados

Para la convección natural en la cavidad cuadrada se obtuvieron resultados para números de Rayleigh de 10³, 10⁴, 10⁵ y 10⁶, y un número de Prandtl Pr = 0,71. La temperatura promedio utilizada en el estado de equilibrio es T_p = (T_h + T_c)/2. Los principales valores numéricos que caracterizan el flujo están dados por el número de Nusselt promedio en la pared caliente, el máximo valor de velocidad vertical y horizontal y las respectivas posiciones donde estas se producen. El número de Nusselt se calcula mediante:

Donde el gradiente de temperatura se obtuvo utilizando un esquema de segundo orden. Las comparaciones con los resultados de la literatura existente se llevan a cabo para 10³ < Ra < 10⁶ y Pr = 0,71. Se analizaron los valores de las cantidades que caracterizan el flujo: el número de Nusselt promedio en la pared caliente y las máximas velocidades horizontal y vertical del flujo, con las respectivas ubicaciones donde se producen. En la Tabla 1 se muestran dichos valores. Los resultados obtenidos con el LBEM, descrito en este trabajo, son muy similares a los existentes en la literatura (un error relativo promedio alrededor del 1,5%).

En la Figura 3 se muestran las isotermas y las líneas de corriente del problema de convección natural de la cavidad cuadrada para distintos números de Ra y para un Pr = 0,71. A medida que aumenta el número de Rayleigh, el movimiento del fluido tiene lugar principalmente cerca de las paredes diferencialmente climatizada y el flujo en el núcleo de la cavidad se convierte en casi inmóvil: las características de este flujo son bien capturadas por el método numérico propuesto.

3. Conclusiones

Se revisó el comportamiento de un flujo laminar por convección libre, simulado por medio del modelo LBM-BGK con simple tiempo de relajación. Los resultados obtenidos permiten observar el buen comportamiento del método propuesto. Respecto a otros métodos, el propuesto ha demostrado tener una buena estabilidad numérica y permitir la inclusión de efectos de calentamiento viscoso. Se han implementado condiciones de frontera que permiten establecer directamente la temperatura en la pared (condiciones del tipo Dirichlet) o condiciones que permiten imponer un flujo de calor en la pared (condiciones de Noumman) transformadas utilizando una aproximación de segundo orden.

Además, en el método utilizado en el presente trabajo ninguna ecuación de Poisson se debe resolver y las condiciones de contorno irregular son fáciles de manejar. Esto hace los métodos de LBEM sean herramientas competitivas, en relación con los enfoques teóricos habituales y con las técnicas numéricas estándar, para la simulación de fenómenos hidrodinámicos complejos.

Referencias

BHATNAGAR, E L.; GROSS, E. P y KROOK. M. A model for collision processes in gases. I Small amplitude processes in charged and neutral one-component system. Physical Review. 1954, vol. 94, pp. 511-525.        [ Links ]

CHEN, Y; OHASHI, H. y AKIYAMA, M. A. Thermal lattice Bhatnagar-Gross-Krook model without nonlinear deviations in macrodynamic equations. Physical Review E. 1994; vol. 50, núm. 4, pp. 2776-2783.        [ Links ]

D'ORAZIO, A.; MASSIMO, C. y GIAN PIERO, C. Application to natural convection enclosed flow of a lattice Boltzmann BGK model coupled with a general purpose thermal boundary conditions. International Journal of Thermal Sciences. 2004, vol. 43, pp. 575-586.         [ Links ]

DAZHI, Y; RENWEI, M.; LI-SHI, L. y WEI, S. Viscous flow computations with the method of lattice Boltzmann equations. Progress in Aerospace Sciences. 2003, vol. 39, núm. 5, pp. 329-367.        [ Links ]

FILIPPOVA, O. y PANEL, D. A novel BGK approach for low Mach number combustión. Journal of Computations Physics. 2000, vol. 158, núm. 2, pp. 139-160.        [ Links ]

FLOREZ, S. E.; CUESTA, I. y SALUEÑA, C. Flujo de Poiseuille y la cavidad con pared móvil calculado usando el método de la ecuación de lattice Boltzmann. Ingeniería & Desarrollo. 2008, vol. 1, núm. 24, pp. 117-132.        [ Links ]

GUO, Z.; SHI, B. y ZHENG, C. A coupled lattice BGK model for the Boussinesq equations. International Journal for Numerical Methods in Fluids. 2002, vol. 39, pp. 325-342.        [ Links ]

HE, X.; CHEN, S. y DOOLEN, G. D. A novel thermal model for the lattice Boltzman method in incompresible limit. Journal of Computational Physics. 1998, vol. 146, núm. 1, pp. 282-300.         [ Links ]

HIGUERA, F. J. y JIMÉNEZ, J. Boltzmann approach to lattice gas simulation. Europhysics Letters. 1989, vol. 9, núm. 7, pp. 663-668.        [ Links ]

KNOW, Y W Development of coupling technique for LBM and FEM for FSI application. InternationalJournal for Computer-Aided Engineering and Software. 2006, vol. 23, núm. 8, pp. 860-875.        [ Links ]

KUZNICK, F.; VAREILLES, J.; RUSAOUEN, G. y KRAUSS, G. A double-population lattice Boltzmann method whit non-uniform mesh for the simulation of natural convection in square cavity. International Journal of Heat and Fluid Flow. 2007, vol. 28, pp. 862-870.         [ Links ]

KUZNIK, F. ; VAREILLES, J. ; RUSAOUEN, G. y KRAUSS, G. Application to natural convection enclosed flows of a lattice Boltzmann BGK model coupled with a general purpose thermal boundary conditions. InternationalJournal of Thermal Sciences. 2004, vol. 43, pp. 575-586.         [ Links ]

QUIAN, Y; D'HUMIERES, D. y LALLEMAND P Lattice BGK models for Navier-Stokes Equation. Europhysics Letters. 1992, vol. 17, pp. 479-84.        [ Links ]

QUIAN, Y H. y ORSZAG, S. A. Lattice BGK models for Navier-Stokes equation: non-linear deviation in compressible regimes. Europhysics Letters. 1993, vol. 21, pp. 255-259.        [ Links ] ]]> 1954 94 511-525 1994 50 4 4 2776-2783 2004 43 575-586 2003 39 5 5 329-367 2000 158 2 2 139-160 2008 1 24 24 117-132 2002 39 325-342 1998 146 1 1 282-300 1989 9 7 7 663-668 2006 23 8 8 860-875 2007 28 862-870 2004 43 575-586 1992 17 479-84 1993 21 255-259