Melhorando a Convergência no ANSYS CFX – Parte 3: Avaliando um Caso Prático

Autor: João Victor Barbosa Alves

Email: joaovictorbr@hotmail.com

Este é o primeiro artigo da InFinite Simulation idealizado e executado por um leitor! João Victor Barbosa Alves é Engenheiro Químico e de Segurança, Especialista em Projetos Avançados de CFD. Também é Professor de Pós-Graduação no Instituto ESSS e Consultor de CFD.

Introdução

Este artigo é uma continuação da série sobre convergência no ANSYS CFX. Nesta terceira postagem, serão vistos na prática os conceitos apresentados na postagem anterior e os seus impactos na convergência da simulação.

O exemplo escolhido para avaliarmos a convergência de uma simulação foi um caso simples de um tubo circular em regime laminar. Os resultados desta análise serão comparados com a solução analítica quanto à perda de carga, ao comprimento equivalente de entrada, ao perfil de velocidades e ao valor da velocidade máxima teórica.

Considerações

As equações presentes neste artigo seguem as seguintes considerações:

  • Escoamento laminar;
  • Escoamento permanente;
  • Escoamento incompressível;
  • Escoamento completamente desenvolvido (ex. perda de carga);
  • Tubo horizontal.

Número de Reynolds

O Número de Reynolds nos permite determinar em qual regime de escoamento um fluido se encontra. Para escoamentos internos são apresentadas as seguintes relações:

  • Escoamento laminar: Re ≤ 2.300
  • Escoamento em transição: 2.300 ≤ Re ≤ 4.000
  • Escoamento turbulento: Re ≥ 4.000

O número de Reynolds para escoamentos internos em dutos circulares é definido pela relação:

Re =\frac{Dv_{med}\rho}{\mu}           Eq. (1)

Onde D = diâmetro; vméd = velocidade média; ρ = densidade; μ = viscosidade e Re = número de Reynolds.

Comprimento de entrada e Perfil de velocidade

O perfil de velocidade em um escoamento laminar varia desde a entrada até um certo comprimento, chamado de comprimento de entrada, Equação 2. No momento em que o perfil de velocidade deixa de variar significa que entramos na região onde o escoamento se encontra completamente desenvolvido. Podemos observar que a camada limite vai se expandindo até encontrar a linha de centro da tubulação (FOX e MCDONALDS, 2001), conforme apresentado na  Figura 1.

figura 1

Figura 1: Escoamento na região de entrada de um tubo com escoamento laminar.
Fonte: Adaptado de CENGEL e CIMBALA (2018).

O comprimento de entrada pode ser obtido com a relação:

L_e = 0.06ReD           Eq. (2)

Onde Le é comprimento de entrada.

O perfil de velocidade completamente desenvolvido, para escoamento laminar, apresenta velocidade máxima na linha de centro do tubo e o seu valor equivale a duas vezes a velocidade média, conforme pode ser verificado na Equação 3.

v(r) = 2v_{med}\left[ 1 - \left( \frac{r}{R} \right)^2 \right]           Eq. (3)

Onde r = posição radial; R = raio interno da tubulação; v = perfil de velocidade radial.

Queda de pressão e Perda de carga

A queda de pressão em escoamento laminar pode ser calculada analiticamente para um tubo linear, com escoamento completamente desenvolvido, com área da seção reta constante através das Equações 4 e 5 (FOX e MCDONALDS, 2001; CENGEL e CIMBALA, 2018).

O fator de atrito de Darcy para um tubo circular em regime laminar é representado pela Equação 5.

\Delta P = f \frac{L}{D}\frac{\rho v^2_{med} }{2}           Eq. (4)

f = \frac{64}{Re}           Eq. (5)

Onde f= fator de atrito de Darcy.

A queda de pressão também pode ser expressa em termos de altura de coluna de fluido equivalente (normalmente água) e com isso passa a ser chamada de perda de carga (CENGEL e CIMBALA, 2018), conforme apresentado na Equação 6.

h_t = \frac{\Delta P}{\rho g} = f \frac{L}{D}\frac{v^2_{med} }{2g}           Eq. (6)

Geometria

A geometria escolhida para o problema consiste de um tubo de 0,05 m de diâmetro e com 6 metros de comprimento, Figura 2.

figura 2

Figura 2: Geometria do tubo

Malha

Quatro diferentes malhas foram preparadas para a avaliação da convergência deste problema. Com isso, poderemos analisar como alguns dos fatores apresentados no artigo anterior interferem na convergência e na qualidade do resultado em caso prático. São eles:

  • Alinhamento dos elementos com o escoamento;
  • Elementos mais regulares possíveis;
  • Elementos com razão de aspecto (relação entre dimensões, aspect ratio em inglês) próxima a 1;
  • Elementos prismáticos próximos às paredes.

A primeira malha é composta por elementos tetraédricos sem a presença da camada de prismas (MTSP), Figura 3.

figura 3

Figura 3: Malha com elementos tetraédricos sem a presença de camada de prismas (MTSP). Visão isométrica (esquerda) e visão frontal (direita).

A segunda malha é composta por elementos tetraédricos com a presença de camada de prismas (MTCP), Figura 4.

figura 4

Figura 4: Malha com elementos tetraédricos com a presença de camada de prismas (MTCP). Visão isométrica (esquerda) e visão frontal (direita).

A terceira malha é composta por elementos hexaédricos e foi construída com o programa meshing utilizando o método sweep (MH-Meshing), Figura 5.

figura 5

Figura 5: Malha com elementos hexaédricos feita no meshing com o método sweep (MH-Meshing). Visão isométrica (esquerda) e visão frontal (direita).

A quarta malha é composta por elementos hexaédricos construídos com o programa ICEM utilizando o método de blocos (MH-ICEM), Figura 6.

figura 6

Figura 6: Malha com elementos hexaédricos feita no ICEM com o método de blocos (MH-ICEM). Visão isométrica (esquerda) e visão frontal (direita).

Os elementos das malhas hexas são menores que os elementos das malhas tetra. Isso ocorreu pelo fato das malhas tetras gastarem uma quantidade muito maior de elementos, o que superou a quantidade disponível na licença acadêmica. A malha hexa gasta menos elementos pelo fato de podermos alonga-los mais na direção do escoamento e manter a qualidade do resultado.

Modelos e Condições de Contorno

Na Tabela 1 e Tabela 2 são apresentados os principais modelos, propriedades e condições de contorno.

Tabela 1: Condições de contorno e propriedades utilizadas nas simulações.

Parâmetro Informação
Velocidade na entrada 0,03 m/s
Pressão na saída 0 Pa
Paredes Condição de aderência
Fluido água
Densidade 997 Kg/m3
Viscosidade dinâmica 8,899×10-4 Kg/(m.s)
Temperatura 25 °C

Tabela 2:  Principais modelos empregados nas simulações numéricas

Parâmetro Opção
Modelo de transferência de calor Isotérmico
Modelo de turbulência Laminar
Esquema de advecção High Resolution
Passo de tempo 2,27 s
Regime temporal Estacionário

(pseudo-transiente)

Convergência

Como qualquer método numérico iterativo, existem dois critérios de convergência básico para as simulações de CFD em regime permanente. A primeira é definida pelo número de iterações e a segunda pelo valor do resíduo, que pode ser o resíduo médio (RMS) ou resíduo máxima (MAX) no CFX.

Geralmente os usuários costumam avaliar a convergência de um problema pelos resíduos. Entretanto, essa pode não ser a melhor escolha, já que o nível do resíduo pode variar em relação as físicas e equações utilizadas.

A solução deve ser considerada convergida quando os balanços estiverem fechados (ex. massa, energia, espécies, etc.). Além disso, devemos monitorar algumas variáveis pertinentes ao problema avaliado, como queda de pressão, velocidade, temperatura, turbulência ou concentração em um certo ponto. A partir do momento em que observarmos que estas propriedades não mais variam podemos acreditar que a simulação atingiu a convergência de uma forma mais segura.

Mesmo obtendo a convergência de um problema isso não significa dizer que o problema está resolvido adequadamente ou possui solução precisa.

A precisão da solução depende, entre outros, dos seguintes fatores:

  • Malha adequada (testada a independência)
  • Condições de contorno corretas (ex. perfil desenvolvido nas entradas)
  • Limitação dos modelos utilizados (ex. turbulência, combustão, …)
  • Simplificações (ex. modelos, geometria, propriedades, …)
  • Precisão numérica (precisão simples ou precisão dupla)
  • Ordem de discretização (ex. segunda ordem)

Por essas razões a validação do problema é de grande importância para as simulações de CFD. Isso nos permite saber se estamos utilizando os modelos, regimes de escoamento, condições de contorno corretas, entre outras coisas. A validação pode utilizar problemas similares, soluções analíticas ou alguma medição de campo em um modelo físico.

Quando não conseguimos convergir uma solução podemos utilizar alguma das estratégias abaixo:

  • Iniciar com uma condição inicial mais próxima da solução
  • Aumentar a complexidade do problema aos poucos (ex. rodar laminar para simulação turbulenta; resolver o escoamento antes da troca térmica)
  • Melhorar a qualidade da malha utilizada
  • As altas ordens de discretização são mais precisas, mas podem apresentar dificuldade de convergência (lentidão ou divergência). Por isso, pode inicializar com um caso com Upwind, por exemplo.
  • Utilizar um caso estacionário como condição inicial para uma simulação transiente.

Resultados

Malha com elementos tetraédricos sem a presença de camada de prismas (MTSP)

Podemos observar na Figura 7 que a queda de pressão está convergida na iteração 50 e neste ponto o pior resíduo se encontra em 5,4×10-4. Já o gráfico da velocidade máxima na saída atinge a convergência na iteração 80 onde o pior resíduo apresenta o valor de
1,0×10-5.

figura 7

Figura 7: Gráficos de convergência da malha com tetraédricos se a presença de prismas (MTSP). Resíduo das equações (a), Queda de pressão (b) e Velocidade máxima na saída do tubo (c).

Malha com elementos tetraédricos com a presença de camada de prismas (MTCP)

Podemos observar na Figura 10 que a queda de pressão está convergida na iteração 50 e neste ponto o pior resíduo se encontra em 6,3×10-5. Já o gráfico da velocidade máxima na saída atinge a convergência na iteração 70 onde o pior resíduo apresenta o valor de
1,0×10-5.

figura 8

Figura 8: Gráficos de convergência da malha com tetraédricos com a presença de prismas (MTCP). Resíduo das equações (a), Queda de pressão (b) e Velocidade máxima na saída do tubo (c).

Malha com elementos hexaédricos feita no meshing com o método sweep (MH-Meshing)

Podemos observar na Figura 9 que a queda de pressão está convergida na iteração 40 e neste ponto o pior resíduo se encontra em 4,9×10-5. Já o gráfico da velocidade máxima na saída atinge a convergência na iteração 80 onde o pior resíduo apresenta o valor de
6,5×10-7.

figura 9

Figura 9: Gráficos de convergência da malha com elementos hexaédricos feita no meshing com o método sweep (MH-Meshing). Resíduo das equações (a), Queda de pressão (b) e Velocidade máxima na saída do tubo (c).

Malha com elementos hexaédricos feita no ICEM com o método de blocos (MH-ICEM)

Podemos observar na Figura 10 que a queda de pressão está convergida na iteração 40 e neste ponto o pior resíduo se encontra em 1,0×10-4. Já o gráfico da velocidade máxima na saída atinge a convergência na iteração 80 onde o pior resíduo apresenta o valor de
8,6×10-7.

figura 10

Figura 10: Gráficos de convergência da malha com elementos hexaédricos feita no ICEM. Resíduo das equações (a), Queda de pressão (b) e Velocidade máxima na saída do tubo (c).

Comparação entre os resultados

As malhas hexaédricas apresentam resultados independentes da malha. Entretanto, foram utilizados mais elementos nas faces para garantir uma boa resolução do perfil de velocidades.

A malha feita no ICEM apresenta menor número de elementos devido à facilidade de distribuir os elementos na face de forma mais adequada. Já as malhas com elementos tetraédricos ainda não se encontram com resultados independentes, em função da quantidade de elementos disponível na licença acadêmica. Entretanto, a experiência em outras casos nos mostra que para atingirmos o mesmo resultado demandaríamos de muitos mais elementos deste tipo.

Na Tabela 3 são apresentados os dados sobre a malha e o tempo de simulação.

Tabela 3: Número de elementos, número de nós e tempo gasto em cada simulação utilizando 4 processadores

Malha Número de elementos Número de nós Tempo de simulação

(min:s)

Hexa-ICEM 91.476 95.300 3:37
Hexa-Meshing 160.020 167.567 5:40
Tetra com prisma 1.031.166 347.118 34:5
Tetra sem prisma 833.894 163.052 23:35

A velocidade máxima em um escoamento laminar atinge o seu valor máximo na linha de centro da tubulação, sendo o seu valor igual ao dobro da velocidade média. Contudo, esse valor é assintótico, ou seja, ele tende a este valor e se aproxima aos poucos.

Na Figura 11 podemos observar que as malhas hexaédricas conseguiram se aproximar do valor de velocidade máxima teórico. Entretanto, as malhas tetraédricas não conseguem atingir este valor. A malha tetra sem prisma já se torna diferente deste valor logo no início da simulação, demonstrando a importância da camada de prisma mesmo para um problema laminar. Já a malha tetra com prisma apresenta bons resultados para os trechos iniciais, que tende a uma variação linear, mas com o aumento da curvatura do gráfico os resultados distam do teórico, apresentando uma difusão numérica.

figura 11

Figura 11: Velocidade máxima na linha de centro da tubulação

O perfil de velocidade analítico, utilizando a Equação 3, foi comparado com os resultados da simulação numérica na saída da tubulação. O resultado dessa comparação é apresentado na Figura 12.

Assim como o que já ocorreu nos resultados anteriores, apenas as simulações feitas com elementos hexaédricos foram capazes de capturar o perfil analítico.

figura 12

Figura 12: Comparação do perfil de velocidade radial na saída da tubulação das simulações numéricas e o perfil analítico.

O comprimento de entrada é uma variável bem difícil de ser obtido, devido ao fato de a velocidade máxima ser assintótica. O valor analítico foi obtido com a Equação 1 e o seu valor foi de 5,04 m, ou seja, a partir deste valor não é esperada variação no perfil de velocidade.

Uma outra forma de avaliar o comprimento de entrada é traçar o gráfico do perfil de velocidade para vários trechos da tubulação e verificar em que momento o perfil deixa de variar, conforme apresentado na Figura 12. Essa avaliação foi feita apenas para a simulação utilizando malha hexaédrica MH-ICEM, em virtude da quantidade de dados gerados e para dar maior clareza a avaliação.

É possível verificar na Figura 12 que na entrada o perfil de velocidade é representado pela média da velocidade. A medida que o escoamento avança podemos verificar que o perfil vai se tornando cada vez mais parabólico. Além disso, verifica-se que apartir de 5 m de distância da entrada o perfil é praticamente o mesmo da saída. Este valor se aproxima bastante do valor previsto pela Equação 1.

figura 13

Figura 13: Perfil de velocidades para diferentes trechos da tubulação.

O último resultado a ser avaliado foi a queda de pressão. Ao utilizar as Equações 4 e 5 no trecho total da geometria, representado pelos pontos 1 e 2 na Figura 2 obtém-se o valor de 2,05 Pa. Mas, a simulação numérica para esse trecho encontrou valores de 2,61 Pa com as malhas Hexaédricas, o que corresponde a um erro de aproximadamente 27%. Entretanto, a condição necessária para utilização das Equações 4 e 5 considera que o escoamento deva ser completamente desenvolvido.

Como o comprimento da entrada (Le) é igual a 5,04 m, a queda de pressão teórica passa a ser igual a 0,33 Pa e a obtida com a melhor simulação é igual à 0,35 Pa. Um erro relativo de 5,50%, conforme apresentado na Tabela 4.

Tabela 4: Valores analíticos e numéricos para a velocidade máxima, comprimento de entrada e queda de pressão.

Teórico MTSP MTCP MH-MESHING MH-ICEM
vmáx (m/s) 0,06 0,0429 0,0457 0,0595 0,0598
Le (m) 5,04 5,6 4,73
ΔP* (Pa) 2,05 2,76 2,99 2,61 2,61
ΔP** (Pa) 0,327 0,430 0,430 0,344 0,346

* Queda de pressão referente ao comprimento total de 6 m; ** Queda de pressão referente ao trecho plenamente desenvolvido de 0,96 m.

Conclusões

Os resultados desta análise foram comparados com a solução analítica quanto à perda de carga, ao comprimento equivalente de entrada, ao perfil de velocidades e ao valor da velocidade máxima teórica.  E chegou-se as seguintes conclusões:

Os elementos hexaédricos apresentaram excelentes resultados, utilizando menos elementos, já que estes estão alinhados com o escoamento;

Dependendo do escoamento, os gradientes podem ser maiores na direção radial e menores na direção longitudinal. Com isso a relação entre dimensões (aspect ratio) podem ser bem maiores que um para casos similares a estes;

Os elementos prismáticos são de grande importância nas proximidades da parede, mesmo para casos laminares.;

Alguns casos podem apresentar resíduos baixos e não estarem adequadamente convergidos. Por essa razão, é de grande importância monitorarmos o fechamento dos balanços (ex. massa, energia, espécies, etc.) e algumas variáveis pertinentes ao problema avaliado, como queda de pressão, velocidade, temperatura, turbulência ou concentração em um certo ponto;

A partir do momento em que observarmos que estas propriedades não mais variam podemos acreditar que a simulação atingiu a convergência de uma forma mais segura;

Verificamos que mesmo obtendo a convergência do problema isso não significa dizer que o problema está resolvido adequadamente ou possui solução precisa;

A precisão da solução depende, entre outros, dos seguintes fatores:

  • Malha adequada (testada a independência)
  • Condições de contorno corretas (ex. perfil desenvolvido nas entradas)
  • Limitação dos modelos utilizados (ex. turbulência, combustão, …)
  • Simplificações (ex. modelos, geometria, propriedades, …)
  • Precisão numérica (precisão simples ou precisão dupla)
  • Ordem de discretização (ex. segunda ordem)

Por essas razões a validação do problema é de grande importância para as simulações de CFD;

Referências

Fox, R.W; McDonald, A.T.; Introdução à Mecânica dos Fluidos. 5aed. Editora LTC. 2001.

Cengel, Y.; Cimbala, J.; Fluid Mechanics: Fundamentals and Applications. 4thed. Editora  McGRAW-HILL. 2018.

 

Autor: João Victor Barbosa Alves

Email: joaovictorbr@hotmail.com

Este é o primeiro artigo da InFinite Simulation idealizado e executado por um leitor! João Victor Barbosa Alves é Engenheiro Químico e de Segurança, Especialista em Projetos Avançados de CFD. Também é Professor de Pós-Graduação no Instituto ESSS e Consultor de CFD.

Deixe um comentário

Preencha os seus dados abaixo ou clique em um ícone para log in:

Logotipo do WordPress.com

Você está comentando utilizando sua conta WordPress.com. Sair /  Alterar )

Foto do Google+

Você está comentando utilizando sua conta Google+. Sair /  Alterar )

Imagem do Twitter

Você está comentando utilizando sua conta Twitter. Sair /  Alterar )

Foto do Facebook

Você está comentando utilizando sua conta Facebook. Sair /  Alterar )

Conectando a %s