Uma formulação de volumes finitos para a solulção de problemas convectivos-difusivos, para quaisquer números de Peclet

Revista Internacional de Métodos Numéricos para Cálculo y Diseño en Ingeniería. Vol. 11,4, 669-681(1995) UMA FORMULACÁO DE VOLUMES FINITOS PARA A SOLUCÁO DE PROBLEMAS CONVECTIVOS-DIFUSIVOS, PARA QUAISQUER NÚMEROS DE PECLET CARLOS M. FERNANDES e MAURI FORTES Departamento de Engenharia Mecanica Universidade Federal de Minas Gerais Belo Horizonte - Minas Gerais - Brazil Apresenta-se uma metodologia de volumes finitos para solucáo de problemas convectivosdifusivos, que pode ser aplicada a escoamentos envolvendo quaisquer valores de números locais de Peclet. O método tem por base a integracáo do termo de adveccáo ao longo da linha de corrente e uma integras50 especial dos termos difusivos, pela utilizaciio de médias de fluxos nas fronteiras dos volumes de controle. O método proposto mostrou-se bastante eficaz e acurado para resolver problemas complexos que normalmente levam a dispersoes e difusóes numéricas. SUMMARY A finite volume method for the solution of convective-diffusive problems is presented, which can be applied to flows with any Peclet numbers. The method is based on the integration of the advective term along a stream line and a special integration of the diffusive terms, using average flow at the boundary at the control volumes. The proposed method has shown to be efficient and acurate in complex applications normally leading to numeric dispersion and diffusion. O método de volumes finitos, MVF, (Spalding17), desenvolvido visando a solucáo numérica das equacóes de Navier-Stokes e de problemas convectivos-difusivos, mostrouse bastante competitivo em relacáo aos métodos de diferencas finitas (MDF) e de elementos finitos (MEF). As solucóes numéricas obtidas por meio de técnicas do MVF (Patankar7, Raithbyl1~l2, Sharif e Busnaina15; Patel et MEF (Rice e Schnipke13, Brooks e Hughesl) e MDF (Hoffman4), apresentam problemas (erros) numéricos, chamados de falsa difusáo e dispersáo numérica. Sharif e Busnaina15 fazem distincáo Recibido: Enero 1995 OUniversitat Politecnica de Catalunya (España) ISSN 0213-1315 670 C.M. FERNANDES E M. FORTES entre difusao e dispersa0 numéricas: pode-se definir difusa0 numérica como qualquer efeito que tende a suavizar ou amortecer gradientes ou descontinuidades presentes na soluciio exata de um problema, enquanto que dispersa0 numérica é o efeito que resulta em oscilacoes na soluqao aproximada. Um grande número de métodos ou esquemas visando a diminuicao de difusiio e dispersiio numéricas encontram-se disponíveis, na literatura, em revisoes abrangentes aplicadas ao MDF (Hoffman4), MVF (Patankar7, Patel et al.', Marchi5) e MEF (Rice e Schnipke13, Brooks e Hughesl; Zienkiewicz e Taylor18). As técnicas de MDF, MVF e MEF mais citadas na literatura referente a soluqGes de equacoes de convecciiodifusao, levaram aos esquemas de Diferencas a Montante (Upwind Difference Scheme ou UDS) (Brooks e Hughesl, Patankar7; Patel et al.') , de Diferencas Centrais (Central Difference Scheme ou CDS) (Brooks e Hughesl, Patankar7, Patel et al.'; Nieckele6), Exponencial (Locally-Exact Difference Scheme ou LEDS) (Spalding17, Patankar7) ; Híbrido (Hybrid Difference Scheme ou HDS) (Spalding17); Lei de Potencia (Power Difference Scheme ou PDS) (Patankar7, Patel et al.'); de Diferencas a Montante Assimétrico (Ske Upstream Difference Scheme ou SUDS) (Raithbyl1,l2); de Diferencas Ponderadas a Montante (Weighted-Upstream Difference Scheme ou WUDS) (Raithby e Torrancel0); de Interpolaqk Quadrática a Montante para Cinemática Convectiva (Quadratic Upstream Interpolation for Convective Kinematics ou QUICK); Diferencial Adaptável (Adaptable Difference Scheme ou ADS) (Marchi5); MEF com formulacao de Galerkin (Brooks e Hughesl); MEF do tipo montante (Upwind) na direeso das linhas de corrente, que utilizam a fomulacao de Petrov-Galerkin, mínimos quadrados, esquemas de interpolacao de alta ordem e outras aproximacoes (Brooks e Hughesl, Rice e Schnipke13) Como tem sido evidenciado na literatura (Brooks e Hughesl, Smith e Hutton16, Patel et Marchi5; Zienkiewicz e Taylor18), praticamente todos os esquemas numéricos de MDF, MEF e MVF introduzem difusa0 e dispersiio numéricas nas solucoes. Tendo por base a revisa0 sucinta apresentada acima e os seguinte argumentos: nao foi encontrado, na literatura revista, um esquema que pudesse ser considerado ótimo, em termos de minimizacao de falsa difusao e dispersa0 numérica para uma vasta gama de números de Peclet; m que existe a possibilidade de explora$io de técnicas envolvendo todos ou dois dos MVF, MEF e MDF; m que existem problemas que t6m sido aceito, pela literatura de diferencas, volumes e elementos finitos como padroes (benchmark), que, embora nao permitam inferir que um método é "ótimo", permitem, por outro lado, comparar os resultados a serem obtidos por um determinado esquema com vários outros propostos e em uso corrente, os objetivos específicos deste trabalho sao: - apresentar uma metodologia de volumes finitos para a soluciio de problemas convectivos-difusivos, tendo por base um esquema a montante, na direeso das linhas de corrente; - apresentar um conjunto de solucoes de problemas padroes (benchmark) e comparar as solucoes obtidas pelo método proposto com outras soluq6es numéricas disponíveis na literatura. O objetivo deste trabalho é a soluciio numérica da equaciio de difusiio-convecciio ou equaciio de transporte em que u e v s5o as componentes da velocidade nas diresoes das coordenadas x e y, p é a massa especifica, 4 é a variável dependente, S é o termo fonte e r é o coeficiente difusivo. Obtencao de variáveis a montante Neste trabalho, desenvolveu-se uma técnica de discretizaciio dos termos de advecciio, que possui algumas semelhancas com a técnica desenvolvida por Rice e Schnipke13, num esquema do MEF, visando a sua aplicaciio por meio do MVF. Os termos convectivos da equaciio (1) podem ser reescritos em coordenadas de linha de corrente (s,n), Figura 1 , na forma em que S é uma coordenada de linha de corrente e u,, o vetor velocidade tangente a linha de corrente, é dado por: u, = dw Figura 1. Coordenadas de linhas de corrente Substituindo a equac5o (2) na equac5o (l), tem-se 672 C.M. FERNANDES E M. FORTES A malha computacional mostrada na Figura 2 define quatro volumes de controle, que possuem, em comum, o ponto nodal P. Admite-se, a priori, que o ponto P seja um ponto nodal "downwind" , isto é, um ponto a jusante de uma linha de corrente do campo de escoamento. Supóe-se ser conhecido o vetor velocidade no ponto P. O critério de escolha do volume de controle referente ao ponto nodal P pode ser exemplificado como segue. Se o vetor velocidade no ponto P estiver dirigido para o quadrante 1 da malha computacional, Figura 2, as coordenadas (x',y') estar20, necessariamente, nos lados S ou w do quadrante 111 (ou coincidir20 com S ou W). Neste caso, o quadrante 111 torna-se o volume de controle para a integras50 da equa~iio (2). Raciocínio análogo pode ser aplicado no caso das outras direcóes plausíveis de V. LID01 Figura 2. Esquema para visualizac50 do ponto a montante (upwind) de um ponto nodal P A localizaqáo do ponto (x',y') a montante da linha de corrente pode ser avaliada a partir de fatores de interpolac20, como feito por Schnipke e Rice13. Os valores de Fp e Fn s2o definidos pelas expressóes abaixo, em func2o das vazoes mássicas nas fronteiras do volume de controle (Figura 3): F4 F2 Lado 1: Fp = max(min(abs(-, l), 0) e Fn = max(min(abs(-, l), 0) (4) F3 Fl As vazóes mássicas podem ser avaliadas por expressoes da forma (Figura 3) Y S YS F~ = J pu dy - J pv dx YSW YSW Neste trabalho supóe-se que toda propriedade (variável) dependente (u,v,+ por exemplo) pode ser interpolada por um polinomio linear em todo o volume de controle e em sua superfície. Para o caso de fluidos incompressíveis, a integraciio da express50 acima leva a FORMULACAO DE V. F. PARA PROBL. CONV.-DIF., PARA QUAISQUER N. DE PECLET 673 A Figura 3 mostra como Fp e Fn sáo orientados no volume de controle em relaqiio ao nó "downwind" P. Os fatores de interpolaqáo variam de zero a um. Se o ponto a montante estiver no lado 1, o valor Fn será um número entre zero e um e Fp será igual a um. Se o ponto a montante estiver no lado quatro Fn será igual a um e Fp terá valor entre zero e um. Se o ponto a montante coincidir com o ponto nodal SW, entáo, Fp e Fn seráo iguais a um. rt" ', 4 " SW. x'."', a' Figura 3. Vazoes mássicas e fatores de interpolacao relativos a um volume de controle no terceiro quadrante Após calcular os fatores de interpolaqáo, podem-se determinar as coordenadas (x,y) do ponto "upstream" e o valor da propriedade qi neste ponto. Se o ponto nodal P for um ponto "downwind" para a linha de corrente que passa pelo terceiro quadrante, entiio, sup6e-se, neste trabalho, que 4' = (1 - Fp)qiw + (1 - Fn)qis + FpFnqisw (9) Expressoes similares podem ser obtidas para outros nós situados em quadrantes diferentes. Discretizaciio da equacáo de transporte A fim de obter um algoritmo de discretizaqáo da equaqáo (4), deve-se definir o quadrante no qual será feita a integraqáo. Para isso analisa-se a direqiio e o sentido do vetor velocidade. Considere a malha computacional apresentada na Figura 2. Integrando a equaqáo (4), no terceiro quadrante dessa figura, tem-se e n e n e n dxdy = a aqi aqi J J -(r-) dxdy + J J aay(r-) dxdy + JJ S dxd~ (10) dx dx W S W S ws ay WS 674 C.M. FERNANDES E M. FORTES Nesta equaqiio, o termo advectivo pode ser aproximado por sendo Ax e Ay siio as dimensoes do volume de controle e Integrando a equaqiio (lo), tem-se onde e n S' = //S; dxdy W S Admitindo-se que as propriedades nas fronteiras em e, w, n e s sejam obtidas pela média aritmética das propriedades nos pontos nodais adjacentes, a equaqiio (13) pode ser reescrita na forma Fazendo Dn = Ax(~P + rw) Ds = ( wrs + rsw) 4(S~)n ~(SY), onde (sx), = XE - xp e (sx), = xp - xw (18) Substituindo as equacóes (16) e (17) em (15) tem-se Reagrupando os termos nesta equaciio, obtém-se ap4p = UE~E + aw4w + UN~N + as+s + b em que b = (AFpFn + DW + D.s)~sw + D~~Nw + De4SE- - (Dn + Ds)4~ - (De + DW)4s + S' A fim de validar o modelo proposto, solucionaram-se dois problemas que tem sido considerados como padróes (benchmark) para se avaliar a capacidade de um método qualquer (MDF, MEF ou MVF) minimizar falsa difusiio e dispersa0 numéricas. 676 C.M. FERNANDES E M. FORTES Escoamento laminar puramente convectivo Este problema teve o objetivo de analisar a capacidade de o método proposto neste trabalho predizer o transporte puramente convectivo de uma propriedade ao longo de uma linha de corrente. O domínio do problema é uma regia0 quadrada, como mostrado na Figura 4. Presume-se que o campo de velocidade seja uniforme no d~míniol~~~~~; ii entrada da regiáo, em y < y,, a propriedade 4 do fluido vale 1 e em y 2 y,, 4 = 0. Figura 4. Esquema do problema de escoamento puramente convectivo Sendo nulo o coeficiente difusivo r , a equacáo de transporte para o problema, equacao (4) torna-se 84 PUS- = O em que u = cos 8 e v = seno as (23) em que u e v si50 os componentes do vetor velocidade e 8 é o angulo de inclinacáo do vetor velocidade em relacáo eixo de coordenadas x. As condicoes de contorno sáo As condicoes de contorno nos lados de saída do escoamento nao devem, em princípio, ser especificadas, pois a equacao (23) é uma equacáo diferencial hiperbólica de primeira ordem. O problema foi simulado para valores de 8 iguais a 22.5, 45 e 67.5'. A regia0 do escoamento foi discretizada numa malha com 21x21 pontos nodais. Os resultados estáo plotados nos gráficos das Figuras 5 e 6; os dados para 67.5' nao foram mostrados pois coincidiram (simetricamente) com os dados para 22.5'. Vários autores (Brooks e Hughesl, Rice e Schnipkel4, Fernandes e Lemos2, Zienkiewicz e Taylor18; Patankar7) publicaram solucoes numéricas para este problema, obtidas pela aplicacao de vários esquemas dos MEF e MVF para o tratamento do termo convectivo da equaqao de transporte. A maior parte dos esquemas, neles incluidos o de Galerkin e Petrov-Galerkin (Zienkiewicz e Taylor18) e SCDS, apresentaram FORMULA~AO DE V. F. PARA PROBL. CONV.-DIF., PARA QUAISQUER N. DE PECLET 677 oscilacoes espaciais náo físicas e os modelos CDS, UDS, LEDS e seus assemelhados (Patankar7, Brooks e Hughesl) apresentaram problemas de difusáo numérica. O modelo, apresentado neste trabalho, náo mostrou problemas de oscilacoes espaciais nao físicas. Para o Angulo O = 45O, a solucáo numérica coincidiu com a exata, náo ocorrendo qualquer difusáo numérica. Para O = 22.5O e 67.5O ocorreu uma pequena falsa difusáo, que, no entanto, foi igual ou menor que o apresentado por todos os esquemas citados. Figura 5. SolucCio do problema puramente convectivo para = 22.5' Figura 6. Soluciio do problema puramente convectivo para = 45' O problema de Smith e Hutton16 Este caso-teste foi proposto por Smith e Hutton16 com a finalidade de avaliar a performance de vários métodos numéricos ao serem aplicados para obter aproximacoes para a soluqáo de problemas convectivos-difusivos com linhas de corrente curvas (recirculacáo), condicoes complexas de contorno e sob toda a faixa de números de Peclet (Figura 7) O problema de Smith e Hutton, portanto, permite avaliar o desempenho de esquemas computacionais quando da soluciio de problemas altamente difusivos e problemas puramente convectivos. A equacáo (3) foi avaliada no domínio mostrado na Figura 7, onde o campo de velocidade é dado pelas relaqoes abaixo: 2 2 u = 2y(lx ) e v = -2x(19 ) para - 15 x 5 1 e O 5 y 5 1 (25)