VALIDAÇÃO DE UM CÓDIGO NUMÉRICO PARA O ESTUDO DO

ESCOAMENTO BIDIMENSIONAL AO REDOR DE CILINDROS CIRCULARES

Bruno Correa Ferreira

Projeto de Graduação apresentado ao Curso de Engenharia

Naval e Oceânica da Escola Politécnica, Universidade

Federal do Rio de Janeiro, como parte dos requisitos

necessários à obtenção do título de Engenheiro Naval e

Oceânico.

Orientador: Prof. Juan Bautista Villa Wanderley, Ph.D.

Rio de Janeiro

Março de 2015

VALIDAÇÃO DE UM CÓDIGO NUMÉRICO PARA O ESTUDO DO

ESCOAMENTO BIDIMENSIONAL AO REDOR DE CILINDROS CIRCULARES

Bruno Correa Ferreira

PROJETO DE GRADUAÇÃO SUBMETIDO AO CORPO DOCENTE DO CURSO

DE ENGENHARIA NAVAL E OCEÂNICA DA ESCOLA POLITÉCNICA DA

UNIVERSIDADE FEDERAL DO RIO DE JANEIRO COMO PARTE DOS

REQUISITOS NECESSÁRIOS PARA A OBTENÇÃO DO GRAU DE

ENGENHEIRO NAVAL E OCEÂNICO.

Examinada por:

Prof. Juan Bautista Villa Wanderley, Ph.D.

Prof. Sérgio Hamilton Sphaier, Ph.D.

Prof. Richard David Schachter, Ph.D.

RIO DE JANEIRO, RJ - BRASIL

MARÇO de 2015

iii

Ferreira, Bruno Correa

Validação de um Código Numérico para o Estudo do

Escoamento Bidimensional ao Redor de Cilindros

Circulares/Bruno Correa Ferreira – Rio de Janeiro:

UFRJ/Escola Politécnica, 2015.

XI, 68 p.: il.; 29,7 cm.

Orientador: Juan Bautista Villa Wanderley

Projeto de Graduação – UFRJ/POLI/Engenharia Naval e

Oceânica, 2015.

Referência Bibliográficas: p. 66-68.

1. Método dos Volumes Finitos 2. Vibração Induzidas

por Vórtices 3. Influência dos Limites Computacionais 4.

Euler Explícito 5. Runge-Kutta de 2ª Ordem Explícito

I. Wanderley, Juan Bautista Villa. II. Universidade Federal

do Rio de Janeiro, UFRJ, Engenharia Naval e Oceânica.

III. Validação de um Código Numérico para o Estudo do

Escoamento Bidimensional ao Redor de Cilindros

Circulares.

iv

Dedicatória

Dedico este trabalho a Deus primeiramente. A minha esposa e filha, Fernanda e

Isabel. Aos meus pais, Rosimere Correa e Valtecir Ferreira Maciel. Aos meus avós,

Eliete e Ilton. E a todos os meus familiares e amigos que me apoiaram durante essa

empreitada.

v

Agradecimentos

Agradeço primeiramente a Deus, pois tenho a certeza de que Ele me

acompanhou por todo este tempo.

Agradeço a minha esposa e filha, não apenas por serem o meu maior motivo de a

cada dia ser melhor do que ontem, mais também por compreenderem as minhas

inúmeras ausências ao longo de todo esse período da graduação.

Agradeço aos meus pais, pelos infindáveis conselhos e apoio incondicional

durante essa jornada.

Agradeço ao meu orientador e amigo, Professor Ph.D. Juan Wanderley, pela sua

amizade, confiança e muita paciência ao longo da minha graduação e elaboração do

presente trabalho.

vi

Resumo do Projeto de Graduação apresentado à Escola Politécnica/ UFRJ como parte

dos requisitos necessários para a obtenção do grau de Engenheiro Naval e Oceânico.

VALIDAÇÃO DE UM CÓDIGO NUMÉRICO PARA O ESTUDO DO

ESCOAMENTO BIDIMENSIONAL AO REDOR DE CILINDROS CIRCULARES

Bruno Correa Ferreira

Março/2015

Orientador: Prof.. Juan Bautista Villa Wanderley, Ph.D

Curso: Engenharia Naval e Oceânica

No presente trabalho as equações de Navier-Stokes levemente compressíveis são

resolvidas numericamente pelo método dos volumes finitos para o escoamento

bidimensional ao redor de um cilindro. Para tanto, um esquema de segunda ordem no

espaço é utilizado para discretizar tanto os vetores de fluxos invíscidos quanto os

viscosos do escoamento, e quatro métodos explícitos de marcha no tempo são utilizados

para resolver as equações governantes escritas na forma conservativa e em coordenadas

cartesianas, para comparação.

A validação do código computacional para investigações de VIV foi realizada pela

comparação dos resultados numéricos com os dados encontrados na literatura. Eles

mediram o coeficiente de arrasto e de sustentação tanto para os números de Reynolds

40, 100 e 200 quanto a frequência de emissão de vórtices para os números de Reynolds

100 e 200. A comparação com estes resultados mostrou que o código computacional foi

capaz de reproduzir o coeficiente de arrasto, de sustentação e a frequência de emissão

de vórtices.

Palavras-chave: Método dos Volumes Finitos, Vibração Induzidas por Vórtices,

Influência dos Limites Computacionais, Método de Euler Explícito, Métodos de Runge-

Kutta Explícitos.

vii

Abstract of Undergraduate Project presented to POLI/UFRJ as a partial fulfillment of

the requirements for the degree of Engineer.

VALIDATION OF A NUMERICAL CODE FOR TWO-DIMENSIONAL FLOW

STUDY THE CIRCULAR CYLINDER AROUND

Bruno Correa Ferreira

March/2015

Advisor: Juan Bautista Villa Wanderley

Course: Naval and Ocean Engineering

In this study, the slightly compressible Navier-Stokes equations are solved numerically

by the finite volume method for the two-dimensional flow around a circular cylinder.

Second-order space approximations are used to discretize both inviscid and viscous flux

vectors. The discretized governing equations written in the conservative form and in

Cartesian coordinates are integrated in time using four different explicit methods, for

comparison.

The validation of the computational code for VIV investigation was conducted by

comparing the numerical results with the data found in the literature. They measured the

drag and lift coefficients for Reynolds numbers 40, 100 and 200, and the vortex

shedding frequency for Reynolds numbers 100 and 200. The comparison showed that

such computer code was able to reproduce the drag and lift coefficients and the

frequency of vortex shedding.

Keywords: Finite Volumes Methods, Vortex-Induced Vibrations, Influence of

Computational Limits, Explicit Euler Scheme, Explicit Runge-Kutta Scheme.

viii

Sumário

1 INTRODUÇÃO .......................................................................................... 1

2 FORMULAÇÃO MATEMÁTICA ................................................................... 2

2.1 Equações de Navier-Stokes ................................................................................................ 2

2.2 Coeficientes de Força .......................................................................................................... 8

2.3 Coeficientes de Pressão....................................................................................................... 8

3 FORMULAÇÃO NUMÉRICA ....................................................................... 8

3.1 Geração do Domínio e da Malha Computacional ............................................................ 9

3.2 Volume Finito .................................................................................................................... 11

3.3 Propriedades Geométricas dos Elementos da Malha Computacional ......................... 12

3.3.1 Cálculo do Volume ................................................................................................... 13

3.3.2 Cálculo do Vetor Área............................................................................................... 14

3.4 Métodos de Marcha no Tempo ou de Integração........................................................... 14

3.4.1 Método de Euler Explícito ........................................................................................ 15

3.4.2 Métodos de Runge-Kutta de 2ª Ordem Explícitos..................................................... 15

3.5 Discretização do Fluxo Invíscido ..................................................................................... 17

3.6 Discretização do Fluxo Viscoso ........................................................................................ 18

3.7 Dissipação Artificial Não Linear ..................................................................................... 20

3.8 Condições Iniciais e de Contorno .................................................................................... 21

3.8.1 Condições Iniciais ..................................................................................................... 21

3.8.2 Condições de Contorno ............................................................................................. 22

4 ESTUDO DO EFEITO DOS LIMITES DO DOMÍNIO COMPUTACIONAL NO

CAMPO DE ESCOAMENTO ........................................................................................ 25

4.1 Definição de um Domínio Computacional de Referência .............................................. 26

4.2 Definição das Características da Malha Computacional a ser Gerada........................ 34

ix

4.3 Estudo da Influência da Localização das Fronteiras Laterais (Condição de Contorno

de “Freestream”) no Campo de Escoamento .......................................................................................... 35

4.4 Estudo da Influência da Localização da Fronteira de Entrada (Condição de Contorno

de “Inflow”) no Campo de Escoamento .................................................................................................. 39

4.5 Estudo da Influência da Localização da Fronteira de Saída (Condição de Contorno

de “Outflow”) no Campo de Escoamento ............................................................................................... 44

4.6 Estudo do Efeito do Refinamento da Malha no Campo de Escoamento...................... 49

5 ESTUDO DO EFEITO DOS MÉTODOS DE MARCHA NO TEMPO NO CAMPO DE

ESCOAMENTO ......................................................................................................... 53

5.1 Resultados para o Número de Reynolds 40 .................................................................... 53

5.2 Resultados para o Número de Reynolds 100 .................................................................. 57

5.3 Resultados para o Número de Reynolds 200 .................................................................. 61

6 CONCLUSÃO ......................................................................................... 65

7 TRABALHOS FUTUROS ........................................................................... 66

8 REFERÊNCIAS BIBLIOGRÁFICAS ............................................................... 66

x

Lista de Símbolos

𝑎 Velocidade do som

𝐶𝐿 Coeficiente de sustentação

𝐶𝐷 Coeficiente de arrasto

𝐶𝑃 Coeficiente de pressão

𝐷 Diâmetro do cilindro

𝐸𝑒 Vetor de fluxo invíscido na direção x

𝐸𝑣 Vetor de fluxo viscoso na direção x

𝐹𝑒 Vetor de fluxo invíscido na direção y

𝐹𝑣 Vetor de fluxo viscoso na direção y

𝐹𝑥 Força na direção x

𝐹𝑦 Força na direção y

𝑀 Número de Mach

�⃗� Vetor normal à superfície

𝑝 Pressão dinâmica do fluido

𝑃 Vetor definido pelos vetores de fluxo

𝑃𝑒 Vetor definido pelos vetores de fluxo invíscidos

𝑃𝑣 Vetor definido pelos vetores de fluxo viscosos

𝑄 Vetor de variáveis conservadas

𝑅𝑒 Número de Reynolds

𝑆 Superfície da malha computacional

𝑆𝑡 Número de Strouhal

𝑇 Temperatura

xi

𝑡 Tempo

𝑢 Velocidade do fluido na direção x

𝑣 Velocidade do fluido na direção y

𝑉 Volume do elemento da malha computacional

𝑥 Coordenada cartesiana horizontal

𝑦 Coordenada cartesiana vertical

∇ Operador divergente

𝜇 Viscosidade dinâmica

𝜌 Massa específica do fluido

𝜏 Coeficiente de compressibilidade isotérmica

∞ Condição de Free-stream (Escoamento não-perturbado)

1

1 Introdução

“Escoamentos ao redor de cilindros são encontrados em muitas aplicações na

engenharia, como por exemplo: risers, estruturas offshore, pilares de pontes, etc.

Dependendo do número de Reynolds, pode ocorrer a emissão de vórtices no bordo de

fuga do corpo. A emissão de vórtices causa uma distribuição de pressão não-simétrica

sobre o cilindro que força o movimento do corpo. Quando a frequência de emissão de

vórtices coincide com a frequência natural da estrutura, oscilações ressonantes podem

ocorrer com altas amplitudes e possivelmente com a falha da estrutura. Portanto, é

muito importante predizer, entender e possivelmente suprimir a vibração induzida por

vórtices.”, Wanderley et al. (2008b).

No presente trabalho, um estudo numérico do escoamento bidimensional ao

redor do cilindro fixo é investigado por meio da solução das equações de Navier-Stokes,

considerando-se a formulação para escoamentos levemente compressíveis propostas por

Wanderley (2001) e usando-se um algoritmo baseado no método dos volumes finitos.

Uma malha não-estruturada composta de elementos triangulares é gerada ao redor da

seção circular do cilindro. Este tipo de malha é extremamente eficiente para aproximar

geometrias complexas. Para otimizar a solução, foi modelado um domínio

computacional que se adapta a geometria do corpo e assim possui uma densidade de

elementos maior apenas nas regiões de interesse próximas da seção circular do cilindro

com o intuito de capturar tanto a camada limite da superfície do corpo quanto ao

fenômeno de desprendimento de vórtices. A forma integral das equações governantes

são discretizadas usando-se aproximações espaciais de segunda ordem e os métodos

explícitos de Euler de primeira ordem e Runge-Kutta de segunda ordem são usados para

a integração no tempo. É adicionado um termo de dissipação artificial não-linear

proposto por Mavriplis (1987) para garantir a estabilidade do método numérico

implementado.

Para verificar a precisão do método numérico, os resultados para números de

Reynolds 40, 100 e 200 foram comparados com diversos dados numéricos e

experimentais encontrados na literatura técnica, como por exemplo: os dados numéricos

obtidos por Wanderley et al. (2008a), Rengel (1999) e Herfjord (1995); e os dados

experimentais obtidos por Tritton (1959), Constanceau (1977) e Wieselsberger (1921).

2

A concordância dos resultados se mostrou satisfatória, validando assim o código para o

estudo da vibração induzida por vórtices a baixo número de Reynolds.

2 Formulação Matemática

A presente seção tem como objetivo apresentar a formulação matemática para

escoamentos levemente compressíveis proposta por Wanderley (2001). Ao final deste

capítulo, são apresentadas as equações para a obtenção dos coeficientes de arrasto e

sustentação do cilindro e o coeficiente de pressão do escoamento.

2.1 Equações de Navier-Stokes

Na mecânica dos fluidos, um escoamento incompressível é definido como um

escoamento onde a densidade 𝜌 é considerada constante. Tendo esta definição em

mente, as equações de Navier-Stokes são utilizadas na engenharia para governar

escoamentos onde a compressibilidade do mesmo é considerada muito baixa (𝑀𝑎𝑐ℎ <

0.30). Mesmo sendo uma boa aproximação para escoamentos com baixa

compressibilidade, estas equações governantes são muito difíceis de serem resolvidas

numericamente, requerendo assim um esforça tanto matemático quanto computacional

muito grande.

A formulação matemática para escoamentos levemente compressíveis,

desenvolvida por Wanderley (2001) e apresentada aqui, propõe uma equação para o

campo de pressão apropriada para estes escoamentos considerando sua baixa

compressibilidade. A equação do campo de pressão junto com a equação da quantidade

de movimento linear para escoamentos incompressíveis, geram uma formulação

apropriada que não somente representa corretamente a física destes escoamentos, mas é

também aplicável para soluções numéricas, vide Wanderley et al. (2008a).

A equação do campo de pressão é obtida quando se combina a equação da

continuidade para escoamentos compressíveis, Eq. (3.1), com a definição de

compressibilidade isotérmica, Eq. (3.2), segundo Anderson (1990).

𝜕𝜌

𝜕𝑡+ ∇. (𝜌�⃗� ) = 0 (3.1)

3

𝜏 =1

𝜌

𝜕𝜌

𝜕𝑝|𝑇

(3.2)

onde

𝜏 = {10−5𝑚2 𝑁⁄ , 𝑝𝑎𝑟𝑎 𝑜 𝑎𝑟 (𝑝 = 1𝑎𝑡𝑚);

5−10𝑚2 𝑁⁄ , 𝑝𝑎𝑟𝑎 𝑎 á𝑔𝑢𝑎 (𝑝 = 1𝑎𝑡𝑚). (3.3)

Restringindo o presente estudo para escoamentos isotérmicos com coeficiente de

compressibilidade isotérmica constante, 𝜏 = 𝑐𝑜𝑛𝑠𝑡., podemos integrar a Eq. (3.2) para

obter a Eq. (3.5).

∫ 𝜏 𝑑𝑝𝑝∞

𝑝

= 𝜏(𝑝∞ − 𝑝) = ∫1

𝜌𝑑𝜌

𝜌∞

𝜌

= ln(𝜌∞) − ln(𝜌) = 𝑙𝑛 (𝜌∞

𝜌) (3.4)

𝜌 = 𝜌∞𝑒𝜏(𝑝−𝑝∞) (3.5)

Fazendo uma expansão em série de Taylor da Eq. (3.5), obtemos a Eq. (3.6).

𝜌 = 𝜌∞ + 𝜌∞𝜏(𝑝 − 𝑝∞) +1

2𝜌∞𝜏

2(𝑝 − 𝑝∞)2+. .. (3.6)

Considerando que a compressibilidade isotérmica é muito pequena, vide Eq.

(3.3), serão considerados apenas os dois primeiros termos da série de Taylor da Eq.

(3.6). Com isso, substituindo estes dois termos na equação da continuidade para

escoamentos compressíveis, Eq. (3.1), obtém-se como resultado as Eqs. (3.8), (3.9) e

(3.10).

𝜕𝜌

𝜕𝑡=𝜕𝜌∞

𝜕𝑡+ 𝜌∞𝜏

𝜕𝑝

𝜕𝑡+ 𝜌∞𝜏

2(𝑝 − 𝑝∞)𝜕𝑝

𝜕𝑡+. .. (3.7)

𝜕𝜌

𝜕𝑡=𝜕𝜌∞

𝜕𝑡+ 𝜌∞𝜏

𝜕𝑝

𝜕𝑡 (3.8)

4

𝜕(𝜌𝑢)

𝜕𝑥= 𝜌∞𝜏 (𝑢

𝜕𝑝

𝜕𝑥+ 𝑝

𝜕𝑢

𝜕𝑥) + 𝜌∞(1 − 𝜏𝑝∞)

𝜕𝑢

𝜕𝑥 (3.9)

𝜕(𝜌𝑣)

𝜕𝑦= 𝜌∞𝜏 (𝑣

𝜕𝑝

𝜕𝑦+ 𝑝

𝜕𝑣

𝜕𝑦) + 𝜌∞(1 − 𝜏𝑝∞)

𝜕𝑣

𝜕𝑦 (3.10)

Somando e rearranjando os termos das Eqs. (3.9) e (3.10), obtém-se as Eqs.

(3.11) e (3.12).

𝜌∞𝜏 (𝑢𝜕𝑝

𝜕𝑥+ 𝑝

𝜕𝑢

𝜕𝑥) + 𝜌∞𝜏 (𝑣

𝜕𝑝

𝜕𝑦+ 𝑝

𝜕𝑣

𝜕𝑦) = 𝜌∞𝜏∇. (𝑝�⃗� ) (3.11)

𝜌∞(1 − 𝜏𝑝∞)𝜕𝑢

𝜕𝑥+ 𝜌∞(1 − 𝜏𝑝∞)

𝜕𝑣

𝜕𝑦= 𝜌∞(1 − 𝜏𝑝∞)∇. �⃗� (3.12)

Substituindo as Eqs. (3.8), (3.11) e (3.12) na Eq. (3.1) ou (3.13), resulta a Eq.

(3.14).

𝜕𝜌

𝜕𝑡+ ∇. (𝜌�⃗� ) =

𝜕𝜌

𝜕𝑡+∂(𝜌𝑢)

∂x+∂(𝜌𝑣)

∂y= 0 (3.13)

𝜌∞𝜏 (𝜕𝜌

𝜕𝑡+ ∇. (𝑝�⃗� )) + 𝜌∞(1 − 𝜏𝑝∞)∇. �⃗� = 0 (3.14)

Como na formulação incompressível assume-se que a compressibilidade

isotérmica é nula, a Eq. (3.14) se reduz à Eq. (3.15), que é a equação da continuidade no

caso de fluidos incompressíveis. Tornando assim, muito onerosa a solução numérica das

equações governantes, conforme dito anteriormente.

∇. �⃗� = 0 (3.15)

O que no caso da formulação matemática proposta por Wanderley (2001), uma

equação mais conveniente pode ser obtida se considerarmos a compressibilidade

isotérmica diferente de zero. Uma vez que o objetivo aqui é obter os coeficientes de

sustentação, arrasto e pressão, o valor adotado para a pressão do escoamento livre pode

5

ser arbitrado. Um valor conveniente para a pressão do escoamento livre é 𝑝∞ =1

𝜏, pois

simplifica a Eq. (3.14), conforme Eq. (3.16).

𝜕𝜌

𝜕𝑡+ ∇. (𝑝�⃗� ) = 0 (3.16)

Nas equações da quantidade de movimento linear, não se leva em consideração a

compressibilidade isotérmica do fluido.

Como mencionado acima, a Eq. (3.16) junto com as equações da quantidade de

movimento linear para escoamentos incompressíveis formam um sistema de equações

que representam convenientemente escoamentos levemente compressíveis. A derivada

no tempo da pressão e as componentes da velocidade facilitam substancialmente a

implementação de qualquer método de marcha no tempo, conforme pode ser visto nas

Eqs. (3.17), já adimensionalizadas.

𝜕𝑝∗

𝜕𝑡∗+∂(𝑝∗𝑢∗)

∂x∗+∂(𝑝∗𝑣∗)

∂y∗= 0

(3.17) 𝜕𝑢∗

𝜕𝑡∗+

𝜕

𝜕𝑥∗(𝑢∗2 + 𝑝∗ −

𝑀∞

𝑅𝑒

𝜕𝑢∗

𝜕𝑥∗) +

𝜕

𝜕𝑦∗(𝑢∗𝑣∗ −

𝑀∞

𝑅𝑒

𝜕𝑢∗

𝜕𝑦∗) = 0

𝜕𝑣∗

𝜕𝑡∗+𝜕

𝜕𝑥(𝑢∗𝑣∗ −

𝑀∞

𝑅𝑒

𝜕𝑣∗

𝜕𝑥∗) +

𝜕

𝜕𝑦∗(𝑣∗2 + 𝑝∗ −

𝑀∞

𝑅𝑒

𝜕𝑣∗

𝜕𝑦∗) = 0

A forma adimensional das equações é importante numericamente porque garante

que as variáveis terão valores da ordem da unidade e posteriormente facilitarão as

análises dos resultados. A Eq. (3.17) foi adimensionalizada de acordo com as seguintes

definições, onde o símbolo de * significa que a variável é adimensional:

𝑢∗ =𝑢

𝑎∞, 𝑣∗ =

𝑣

𝑎∞, 𝑥∗ =

𝑥

𝐷, 𝑦∗ =

𝑦

𝐷

(3.18)

𝑡∗ =𝑡

𝐷𝑎∞⁄

, 𝑝∗ = 𝑝

𝜌∞𝑎∞2

6

𝑀∞ =𝑈∞

𝑎∞, 𝑅𝑒 =

𝜌∞𝑈∞𝐿

𝜇∞

Com a intenção de simplificar a escrita das equações, será eliminado o símbolo

de asterisco utilizado para representar as variáveis adimensionais.

O mesmo sistema de equações mostrados na Eq. (3.17) é apresentado a seguir

em coordenadas cartesianas, na forma conservativa e vetorial.

𝑄𝑡 + (𝐸𝑒 − 𝐸𝑣)𝑥 + (𝐹𝑒 − 𝐹𝑣)𝑦 = 0 (3.19)

onde 𝑄 é o vetor de variáveis conservadas , 𝐸𝑒 e 𝐹𝑒 são os vetores de fluxo invíscidos e

𝐸𝑣 e 𝐹𝑣 são os vetores de fluxo viscosos, conforme Eq. (3.20).

𝑄𝑡 = {𝑝𝑢𝑣} , 𝐸𝑒 = {

𝑝𝑢

𝑢2 + 𝑝𝑢𝑣

} , 𝐹𝑒 = {

𝑝𝑣𝑢𝑣

𝑣2 + 𝑝},

(3.20)

𝐸𝑣 =𝑀∞

𝑅𝑒

{

0𝜕𝑢

𝜕𝑥𝜕𝑣

𝜕𝑥}

, 𝐹𝑣 =𝑀∞

𝑅𝑒

{

0𝜕𝑢

𝜕𝑦𝜕𝑣

𝜕𝑦}

Para aplicarmos o método dos volumes finitos, a Eq. (3.19) é transformada para

a sua forma integral. Utilizando a definição do vetor �⃗� dada na Eq. (3.21), a Eq. (3.19)

pode ser escrita na forma apresentada na Eq. (3.22).

�⃗� = (𝐸𝑒 − 𝐸𝑣)𝑖̂ + (𝐹𝑒 − 𝐹𝑣)𝑗 ̂ (3.21)

𝑄𝑡 =𝜕𝑄

𝜕𝑡= −∇. �⃗� = (𝐸𝑒 − 𝐸𝑣)𝑥 + (𝐹𝑒 − 𝐹𝑣)𝑦 (3.22)

Integrando a Eq. (3.22) num volume de controle 𝑉, conforme a Eq. (3.23), e

aplicando o teorema de Gauss ao lado direito da Eq. (3.23), resulta a Eq. (3.24).

7

∫𝜕𝑄

𝜕𝑡𝑑𝑉

𝑉

= −∫ (∇. �⃗� )𝑑𝑉𝑉

(3.23)

∫ (∇. �⃗� )𝑑𝑉𝑉

= ∫ (�⃗� . �⃗� )𝑑𝑆𝑉

(3.24)

Substituindo a Eq. (3.24) na Eq. (3.23), a forma integral da equação de Navier-

Stokes para escoamentos levemente compressíveis é obtida, conforme mostra a Eq.

(3.25).

∫𝜕𝑄

𝜕𝑡𝑑𝑉

𝑉

= −∫ (�⃗� . �⃗� )𝑑𝑆𝑆

(3.25)

Para o caso de volumes finitos rígidos e não deformáveis, que é o presente caso,

podemos passar o símbolo de derivada temporal para fora do sinal de integração. O

resultado é mostrado na Eq. (3.26).

𝜕

𝜕𝑡∫ 𝑄𝑑𝑉𝑉

= −∫ (�⃗� . �⃗� )𝑑𝑆𝑉

(3.26)

O vetor médio das variáveis conservadas em cada volume finito é definido

conforme a Eq. (3.27).

�̅�𝑉 = ∫ 𝑄𝑑𝑉𝑉

(3.27)

Substituindo a Eq.(3.27) na Eq.(3.26), obtemos a formulação integral em termos

do vetor médio das variáveis conservadas, como mostra a Eq. (3.28).

𝜕�̅�

𝜕𝑡= −

1

𝑉∫ (�⃗� . �⃗� )𝑑𝑆𝑉

(3.28)

A Eq. (3.28) é a forma mais conveniente da equação de Navier-Stokes para a

utilização do método dos volumes finitos. Sendo assim, esta equação será avaliada em

cada volume finito presente no domínio computacional discretizado.

8

2.2 Coeficientes de Força

Os coeficientes de força são calculados pela soma da integral de superfície tanto

da pressão quanto da tensão cisalhante. O campo de pressão é conhecido, porém faz-se

necessário o cálculo da tensão cisalhante a partir das variáveis conhecidas 𝑢 e 𝑣. Estes

cálculos serão explicados nas seções 3.5 e 3.6. Com estes valores determinados, os

coeficientes são avaliados nas direções 𝑥 e 𝑦 do escoamento, fornecendo as forças 𝐹𝑥 e

𝐹𝑦 , adimensionalizadas da seguinte forma:

𝐶𝐷 =𝐹𝑥

12⁄ 𝜌∞𝑈∞

2𝐷 (3.29)

𝐶𝐿 =𝐹𝑦

12⁄ 𝜌∞𝑈∞

2𝐷 (3.30)

2.3 Coeficientes de Pressão

O coeficiente de pressão é calculado através da equação abaixo:

𝐶𝑃 =𝑝 − 1

12⁄ 𝜌∞𝑈∞

2𝐷 (3.31)

3 Formulação Numérica

Para capturar o campo de escoamento ao redor da seção circular do cilindro e

reproduzir o fenômeno de desprendimento de vórtices, faz-se necessário um domínio e

uma malha computacional adequados. A geração do domínio e da malha computacional

serão descritos no presente capítulo de uma forma simplificada, deixando os maiores

detalhes para as seções 4.1 e 4.2, respectivamente. Aqui também serão apresentados:

como cada volume finito pode ser obtido através da discretização do domínio

computacional em elementos triangulares (malha não-estruturada); e uma descrição

detalhada da implementação numérica do método de volumes finitos usando um

esquema espacial centrado de segundo ordem para discretizar tanto o fluxo invíscido

quanto o viscoso das propriedades do escoamento.

9

3.1 Geração do Domínio e da Malha Computacional

A ferramenta escolhida tanto para a modelagem da geometria quanto para a

geração da malha, foi o software Gmsh. Este é um gerador de malhas com um sistema

CAD embutido. A escolha de tal ferramenta se deu por este ser um software livre,

portável tanto para as plataformas Linux e Windows, intuitivo e de fácil aprendizado. E

como a documentação de tal software é vasta na internet, no presente trabalho não se

fará menção de como os domínios ou as malhas foram geradas, apenas os pontos mais

importantes serão mencionados, como por exemplo: as características geométricas do

domínio e das malhas computacionais geradas.

Tendo como referência Behr et al. (1991) e Behr et al. (1995), foi definido

inicialmente um domínio computacional de partida, conforme Figura 3.1. Domínio este,

que será utilizado como referência para a construção de todos os domínios

computacionais utilizados nos estudos subsequentes. Este domínio computacional pode

ser visualizado na Figura 3.2, logo abaixo.

Figura 3.1 – Domínio computacional de referência.

Figura 3.2 – Domínio computacional de referência gerado no software Gmsh.

10

A malha computacional foi gerada através do método de avanço de fronteira,

conforme mostrado na Figura 3.3. Malha esta composta por elementos triangulares.

Como se pode observar, foram criadas regiões no entorno da seção circular do

cilindro com o intuito de capturar tanto a camada limite nesta superfície quanto o

fenômeno de desprendimento de vórtices ou a esteira de vórtices de von Kárman no

bordo de fuga do mesmo.

Figura 3.3 – Malha gerada no domínio computacional de referência através do software Gmsh.

Esta malha contém elementos ao redor do domínio físico e no interior do corpo

chamados elementos imagens ou “fantasmas”. O objetivo deles é simplificar a

implementação das condições de contorno. Estes elementos podem ser visualizados nas

Figuras abaixo.

Figura 3.4 – Elementos imagem ou “fantasmas” ao redor do domínio computacional.

11

Figura 3.5 – Detalhe dos elementos imagem ou “fantasmas” ao redor da seção circular do cilindro.

Vale ressaltar, que tais elementos imagem ou “fantasmas” não foram gerados

através do software Gmsh, ou seja, não é uma funcionalidade do mesmo. Para tal, foi

desenvolvido um programa também em linguagem Fortran que gera estes elementos.

3.2 Volume Finito

Uma vez discretizado o domínio computacional com uma malha não-

estruturada, tratando-se de elementos triangulares no presente trabalho, o próximo passo

é obter os volumes finitos através destes elementos. Existem diversas formas que podem

ser utilizadas para criar os volumes finitos a partir de uma malha não-estruturada. Uma

delas é a criação destes através do método das medianas, que simplesmente consiste em

unir os centroides dos elementos com os pontos médios de seus lados, conforme Figura

3.6. Esta é uma forma bem genérica, podendo ser aplicado a uma discretização mista,

envolvendo triângulos e quadriláteros.

(a)

(b)

Figura 3.6 – Volumes de controle obtidos através do método das medianas. (a) Elemento triangular 1; e (b) Volumes

de controle associados aos nós 1, 0 e 3.

Uma outra forma seria quando uma discretização apenas com triângulos fosse

utilizada, e formando uma triangulação de Delanay, os volumes de controle já estariam

12

diretamente determinados e se constituem no dual da triangulação, os chamados

diagramas de Voronoi, conforme Figura 3.7.

Figura 3.7 – Volumes de controle obtidos através do dual da triangulação de Delanay.

Mas para o presente trabalho, os volumes finitos serão definidos como os

próprios elementos triangulares resultantes da geração da malha pelo método de avanço

de fronteira. É neste volume de controle ou elemento que a Eq. (3.28) será avaliada.

Todas as propriedades do escoamento serão armazenadas no centroide destes elementos.

Aqui cada volume de controle possui uma profundidade unitária e seu volume é dado

pela área do triângulo ou elemento.

Figura 3.8 – Discretização do espaço bidimensional em elementos triangulares (volumes de controle).

3.3 Propriedades Geométricas dos Elementos da Malha Computacional

A presente seção tem como objetivo mostrar as propriedades geométricas dos

elementos da malha computacional ou de cada volume finito que serão utilizadas

quando a Eq. (3.28) for avaliada em cada um destes.

13

3.3.1 Cálculo do Volume

Para calcular o volume 𝑉𝑖 de cada elemento pertencente à malha computacional,

definem-se os vetores 𝐴 𝑖,2 e 𝐴 𝑖,3 e os nós 1, 2 e 3 conforme regra da mão direita e

Figura 3.9. O volume de cada elemento triangular é dado pela Eq. (3.32).

Figura 3.9 – Definição dos vetores 𝐴 𝑖,2 e 𝐴 𝑖,3 para o cálculo do volume de cada elemento ou volume de controle

presente na malha computacional.

𝑉𝑖 =1

2|𝐴 𝑖,3 × 𝐴 𝑖,2| (3.32)

𝐴 𝑖,3 = (𝑥𝑖,2 − 𝑥𝑖,1)𝑖̂ + (𝑦𝑖,2 − 𝑦𝑖,1)𝑗 ̂ (3.33)

𝐴 𝑖,2 = (𝑥𝑖,3 − 𝑥𝑖,1)𝑖̂ + (𝑦𝑖,3 − 𝑦𝑖,1)𝑗 ̂ (3.34)

Substituindo as Eqs. (3.33) e (3.34) na Eq. (3.32), obtém-se a Eq. (3.36) que

representa o volume 𝑉𝑖 de cada elemento 𝑖.

𝑉𝑖 =1

2|𝐴 𝑖,3 × 𝐴 𝑖,2| = |

𝑖̂ 𝑗̂ �̂�

(𝑥𝑖,2 − 𝑥𝑖,1) (𝑦𝑖,2 − 𝑦𝑖,1) (0)

(𝑥𝑖,3 − 𝑥𝑖,1) (𝑦𝑖,3 − 𝑦𝑖,1) (0)

| (3.35)

14

𝑉𝑖 = |(𝑥𝑖,2 − 𝑥𝑖,1). (𝑦𝑖,3 − 𝑦𝑖,1)�̂� − (𝑦𝑖,2 − 𝑦𝑖,1). (𝑥𝑖,3 − 𝑥𝑖,1)�̂�| (3.36)

3.3.2 Cálculo do Vetor Área

Para aplicar o método dos volumes finitos faz-se necessário calcular o vetor área

perpendicular a cada face dos volumes de controle. As Eqs. de (3.37) a (3.39), mostram

o vetor área e suas componentes nas direções 𝑥 e 𝑦 do volume de controle 𝑖 em todas as

suas faces, conforme Figura 3.10.

𝑆 𝑖,1 = 𝑆𝑥𝑖,1𝑖̂ + 𝑆𝑦𝑖,1𝑗̂ = (𝑦𝑖,3 − 𝑦𝑖,2)𝑖̂ + (𝑥𝑖,2 − 𝑥𝑖,3)𝑗̂ (3.37)

𝑆 𝑖,2 = 𝑆𝑥𝑖,2𝑖̂ + 𝑆𝑦𝑖,2𝑗̂ = (𝑦𝑖,1 − 𝑦𝑖,3)𝑖̂ + (𝑥𝑖,3 − 𝑥𝑖,1)𝑗̂ (3.38)

𝑆 𝑖,3 = 𝑆𝑥𝑖,3𝑖̂ + 𝑆𝑦𝑖,3𝑗̂ = (𝑦𝑖,2 − 𝑦𝑖,1)𝑖̂ + (𝑥𝑖,1 − 𝑥𝑖,2)𝑗̂ (3.39)

Figura 3.10 – Vetor área perpendicular de cada face do elemento ou volume de controle 𝑖.

3.4 Métodos de Marcha no Tempo ou de Integração

A presente seção tem como objetivo apresentar os métodos de marcha no tempo

que serão alvos de estudo no presente trabalho para a análise do efeito destes no campo

de escoamento. Todos estes métodos serão aplicados à Eq. (3.28).

15

3.4.1 Método de Euler Explícito

Dada uma condição inicial e condições de contorno, que serão apresentadas na

seção 3.8, a Eq. (3.28) será resolvida pelo método explícito de Euler de primeira ordem,

conforme mostrado nas equações abaixo.

�̅�𝑖𝑛+1 = �̅�𝑖

𝑛 +1

2∆𝑡 (

𝜕�̅�𝑖𝑛

𝜕𝑡) (3.40)

Figura 3.11 – Exemplo de elemento ou volume de controle pertencente a malha computacional onde a Eq. (3.28) será

avaliada.

Aplicando a Eq. (3.40) ao volume finito ou elemento mostrado na Figura 3.11,

obtemos as Eqs. (3.41) e (3.42).

�̅�𝑖𝑛+1 = �̅�𝑖

𝑛 −1

2

∆𝑡

𝑉𝑖∑(�⃗� 𝑖,𝑗

𝑛 . 𝑆 𝑖,𝑗)

3

𝑗=1

(3.41)

�̅�𝑖𝑛+1 = �̅�𝑖

𝑛 −1

2

∆𝑡

𝑉𝑖(�⃗� 𝑖,1

𝑛 . 𝑆 𝑖,1 + �⃗� 𝑖,2𝑛 . 𝑆 𝑖,2 + �⃗� 𝑖,3

𝑛 . 𝑆 𝑖,3) (3.42)

A Eq. (3.42) será avaliada em cada volume finito 𝑖 pertencente ao domínio

computacional.

3.4.2 Métodos de Runge-Kutta de 2ª Ordem Explícitos

Existem diversos métodos de Runge-Kutta de 2ª ordem, logo abaixo são

apresentados os três métodos também utilizados no presente trabalho.

16

“Heun Method with a Single Corrector”:

𝑄𝑖𝑛+1

= 𝑄𝑖𝑛+ ∆𝑡

𝜕𝑄𝑖𝑛

𝜕𝑡 (3.43)

𝑄𝑖𝑛+1

= 𝑄𝑖𝑛+ (

1

2

𝜕𝑄𝑖𝑛

𝜕𝑡+1

2

𝜕𝑄𝑖𝑛+1

𝜕𝑡) (3.44)

Este é o método de “Heun Method with a Single Corrector”.

“The Midpoint Method”:

𝑄𝑖𝑛+1 2⁄ = 𝑄𝑖

𝑛+1

2∆𝑡𝜕𝑄𝑖

𝑛

𝜕𝑡 (3.45)

𝑄𝑖𝑛+1

= 𝑄𝑖𝑛+ ∆𝑡

𝜕𝑄𝑖𝑛+1 2⁄

𝜕𝑡 (3.46)

Este é o método “The Midpoint Method”.

“Ralston’s Method”:

𝑄𝑖𝑛+3 4⁄ = 𝑄𝑖

𝑛+3

4∆𝑡𝜕𝑄𝑖

𝑛

𝜕𝑡 (3.47)

𝑄𝑖𝑛+1

= 𝑄𝑖𝑛+ (

1

3

𝜕𝑄𝑖𝑛

𝜕𝑡+2

3

𝜕𝑄𝑖𝑛+3 4⁄

𝜕𝑡) (3.48)

Este é o método “Ralston’s Method”.

Então, estes métodos de marcha no tempo também serão aplicados à Eq. (3.28)

de forma semelhante a demonstrada para o método de Euler explícito.

O passo de tempo que aparece na Eq. (3.28) em sua forma dimensional, ∆𝑡∗, é

definido pela Eq. (3.49).

17

∆𝑡∗ = 𝐶𝐹𝐿.𝑙∗

𝑈∗, (3..49)

onde o 𝐶𝐹𝐿 é o número de Courant-Friedrichs-Lewy, 𝑈∗ e 𝑙∗ são respectivamente a

velocidade e o comprimento característicos do escoamento. O passo de tempo na forma

adimensional é dado pela Eq. (3.50).

∆𝑡 = 𝐶𝐹𝐿. 𝑙, (3.50)

onde

𝑙 = √𝑉𝑚𝑖𝑛 (3.51)

𝑉𝑚𝑖𝑛 representa o volume do menor elemento presente na malha, calculado pela

Eq. (3.52).

𝑉𝑚𝑖𝑛 = 𝐴𝑚𝑖𝑛, (3.52)

onde 𝐴𝑚𝑖𝑛 representa a menor área encontrada na malha computacional.

3.5 Discretização do Fluxo Invíscido

Os vetores de fluxo invíscidos nas faces dos volumes de controle são obtidos a

partir dos valores das propriedades médias definidas em seus centroides. Uma forma

utilizada para obter esses valores é por interpolação das propriedades de dois volumes

de controle adjacentes. Conhecendo-se as propriedades nas interfaces, pode-se

determinar os termos invíscidos da Eq. (3.42), conforme a Eq. (3.53).

�⃗� 𝑖,𝑗𝑛 . 𝑆 𝑖,𝑗 = {

(𝑝𝑢)𝑖,𝑗𝑛 (𝑆𝑥)𝑖,𝑗

(𝑢2 + 𝑝)𝑖,𝑗𝑛 (𝑆𝑥)𝑖,𝑗

(𝑢𝑣)𝑖,𝑗𝑛 (𝑆𝑥)𝑖,𝑗

} +

{

(𝑝𝑣)𝑖,𝑗

𝑛 (𝑆𝑦)𝑖,𝑗

(𝑢𝑣)𝑖,𝑗𝑛 (𝑆𝑦)𝑖,𝑗

(𝑣2 + 𝑝)𝑖,𝑗𝑛 (𝑆𝑦)𝑖,𝑗}

(3.53)

18

�⃗� 𝑖,𝑗𝑛 . 𝑆 𝑖,𝑗 = {

(𝑝)𝑖,𝑗𝑛 (𝑞𝑠)𝑖,𝑗

𝑛

(𝑢𝑞𝑠)𝑖,𝑗𝑛 + (𝑝)𝑖,𝑗

𝑛 (𝑆𝑥)𝑖,𝑗

(𝑣𝑞𝑠)𝑖,𝑗𝑛 + (𝑝)𝑖,𝑗

𝑛 (𝑆𝑦)𝑖,𝑗

} (3.54)

onde

𝑞𝑠𝑖,𝑗𝑛 = (𝑢)𝑖,𝑗

𝑛 (𝑆𝑥)𝑖,𝑗 + (𝑣)𝑖,𝑗𝑛 (𝑆𝑦)𝑖,𝑗

(3.55)

Uma forma bastante utilizada, e aqui praticada, para a obtenção dos vetores de

fluxo não-viscosos 𝐸𝑒 e 𝐹𝑒 nas faces dos volumes de controle é pela média aritmética

dos vetores de fluxo definidos nos centroides dos volumes adjacentes. As Eq. (3.56) e

(3.57) mostram o cálculo dos vetores de fluxo na face 1 do elemento 4497 pela média

aritmética dos vetores de fluxo definidos nos centroides dos elementos 4497 e 4142,

mostrados na Figura 3.12.

(𝐸𝑒)4497,1𝑛 =

1

2[(𝐸𝑒)4497

𝑛 + (𝐸𝑒)4142𝑛 ] (3.56)

(𝐹𝑒)4497,1𝑛 =

1

2[(𝐹𝑒)4497

𝑛 + (𝐹𝑒)4142𝑛 ] (3.57)

Este procedimento de interpolação equivale a uma aproximação centrada de

segunda ordem de precisão em diferenças finitas quando a malha é suave.

3.6 Discretização do Fluxo Viscoso

Os termos viscosos nas faces dos volumes de controle são obtidos calculando-se

a média aritmética dos vetores de fluxo viscosos dos volumes de controle adjacentes

àquela face, conforme a Eq. (3.58).

(�⃗� 𝑣𝑛. 𝑆 )

4497,1=1

2(𝐸𝑣4497

𝑛 + 𝐸𝑣4142𝑛 )𝑆𝑥4497,1

+1

2(𝐹𝑣4497

𝑛 + 𝐹𝑣4142𝑛 )𝑆𝑦4497,1

(3.58)

19

Mas, para este fim, faz-se necessário a obtenção das velocidades médias nos

centros das faces dos volumes de controle. Que também são obtidas através da média

das velocidades dos volumes de controle adjacentes a esta face. Com isso, pode-se

calcular as derivadas das componentes de velocidade 𝑢 e 𝑣 em relação à 𝑥 e 𝑦 para um

volume de controle qualquer. Considerando o volume de controle 4497 da Figura 3.12.

Pelo teorema de Gauss, a derivada de 𝑢 em relação à 𝑥 é mostrada na Eq. (3.59).

Figura 3.12 – Malha computacional não estruturada.

𝜕𝑢

𝜕𝑥)4497

𝑛

=1

𝑉4497∑𝑢4497,𝑗

𝑛 (𝑆𝑥4497,𝑗)

3

𝑗=1

(3.59)

onde

𝑢4497,1𝑛 =

1

2[𝑢4497𝑛 + 𝑢4142

𝑛 ] (3.60)

𝑢4497,2𝑛 =

1

2[𝑢4497𝑛 + 𝑢5064

𝑛 ] (3.61)

𝑢4497,3𝑛 =

1

2[𝑢4497𝑛 + 𝑢4416

𝑛 ] (3.62)

De forma análoga, podemos obter também a derivada de 𝑢 em relação à 𝑦,

conforme mostra a Eq. (3.63).

20

𝜕𝑢

𝜕𝑦)4497

𝑛

=1

𝑉4497∑𝑢4497,𝑗

𝑛 (𝑆𝑦4497,𝑗)

3

𝑗=1

(3.63)

De forma completamente análoga, podemos calcular as derivadas da

componente de velocidade 𝑢 em relação à 𝑥 e 𝑦 no centróide do volume de controle

5064, conforme mostram as Eqs. (3.64) e (3.65).

𝜕𝑢

𝜕𝑥)5064

𝑛

=1

𝑉5064∑𝑢5064,𝑗

𝑛 (𝑆𝑥5064,𝑗)

3

𝑗=1

(3.64)

𝜕𝑢

𝜕𝑦)5064

𝑛

=1

𝑉5064∑𝑢5064,𝑗

𝑛 (𝑆𝑦5064,𝑗)

3

𝑗=1

(3.65)

Finalmente, as derivadas na face entre os volumes de controle 4497 e 5064 são

dadas pela média aritmética dos valores obtidos para os volumes de controle 4497 e

5064, conforme mostra as Eqs. (3.66) e (3.67).

𝜕𝑢

𝜕𝑥)4497,2

𝑛

=1

2[𝜕𝑢

𝜕𝑥)4497

𝑛

+𝜕𝑢

𝜕𝑥)5064

𝑛

] (3.66)

𝜕𝑢

𝜕𝑦)4497,2

𝑛

=1

2[𝜕𝑢

𝜕𝑦)4497

𝑛

+𝜕𝑢

𝜕𝑦)5064

𝑛

] (3.67)

3.7 Dissipação Artificial Não Linear

O processo de interpolação utilizado aqui equivale em diferenças finitas a um

esquema centrado. E esquemas centrados não possuem termos de dissipação artificial

implícitos. Portanto, faz-se necessário adicionar explicitamente termos de dissipação

artificial ao presente esquema para filtrar ondas de instabilidade numéricas de alta

frequência para garantir a estabilidade do método. A Eq. (3.68) mostra a formulação

numérica discretizada espacialmente para elementos triangulares e utilizando o método

de Euler explícito para a integração no tempo com a adição explícita de um termo de

dissipação artificial não-linear típico. Esta equação foi aplicada ao volume de controle

4497, mostrado na Figura 3.12.

21

𝑄𝑖𝑛+1

= 𝑄𝑖𝑛−1

2

∆𝑡

𝑉𝑖(�⃗� 𝑖,1

𝑛 . 𝑆 𝑖,1 + �⃗� 𝑖,2𝑛 . 𝑆 𝑖,2 + �⃗� 𝑖,3

𝑛 . 𝑆 𝑖,3 + 𝐷 (𝑄𝑖𝑛)) (3.68)

onde

𝐷 (𝑄𝑖𝑛) = ∑(𝑑(2))

𝑖,𝑗

𝑛3

𝑗=1

−∑(𝑑(4))𝑖,𝑗

𝑛3

𝑗=1

(3.69)

(𝑑(2))𝑖,𝑗

𝑛=𝜀𝑖,𝑗(2)

2(𝑉𝑖 + 𝑉𝑘

∆𝑡) (𝑄𝑘

𝑛− 𝑄𝑖

𝑛) (3.70)

(𝑑(4))𝑖,𝑗

𝑛=𝜀𝑖,𝑗(4)

2(𝑉𝑖 + 𝑉𝑘

∆𝑡) (∇2𝑄𝑘

𝑛− ∇2𝑄𝑖

𝑛) (3.71)

∇2𝑄𝑖𝑛= ∑𝑄𝑘

𝑛− 3𝑄𝑖

𝑛3

𝑘=1

(3.72)

∇2𝑄𝑘𝑛= ∑𝑄𝑖

𝑛− 3𝑄𝑘

𝑛3

𝑖=1

(3.73)

𝜀𝑖,𝑘(2)= 𝐾(2)𝑚𝑎𝑥(𝜐𝑖, 𝜐𝑘) (3.74)

𝜀𝑖,𝑘(4)= 𝑚𝑎𝑥 [0, (𝐾(4) − 𝜀𝑖,𝑘

(2))] (3.75)

𝜐𝑖𝑛 =

|∑ 𝑝𝑘𝑛3

𝑘=1 − 3𝑝𝑖𝑛|

∑ |𝑝𝑘𝑛|3

𝑘=1 + 3|𝑝𝑖𝑛|

(3.76)

𝐾2 =1

2 𝑒 𝐾4 =

1

256 (3.77)

3.8 Condições Iniciais e de Contorno

3.8.1 Condições Iniciais

Para a simulação do escoamento ao redor da seção circular do cilindro, as

condições iniciais para o escoamento são mostradas abaixo.

22

𝑢 = 𝑈∞ = 0.20, 𝑣 = 0.00 𝑒 𝑝 = 𝑝∞ = 1.00 (3.78)

3.8.2 Condições de Contorno

No domínio computacional, as condições de contorno impostas são as seguintes:

Condição de contorno imposta na fronteira de inflow

A condição de contorno imposta na fronteira de inflow permite a entrada de

fluido com velocidade somente em x (𝑈∞ = 0.20 e 𝑉∞ = 0.00) e a pressão de referência

é igual à 𝑃∞ = 1.00. Estas condições de contorno podem ser escritas matematicamente

da seguinte forma:

𝜕𝑢

𝜕𝑥= 0,

𝜕𝑣

𝜕𝑥= 0 𝑒

𝜕𝑝

𝜕𝑥= 0 (3.79)

Figura 3.13 – Condição de contorno imposta na fronteira de entrada do domínio computacional.

Condição de contorno imposta na fronteira de freestream

As fronteiras superior e inferior (freestream), impõem a condição de escoamento

no infinito (escoamento não-perturbado “freestream” = “𝑈∞ = 0.20, 𝑉∞ = 0.00 e 𝑃∞ =

23

1.00“). Estas condições de contorno podem ser escritas matematicamente da seguinte

forma:

𝜕𝑢

𝜕𝑦= 0, 𝑣 = 0 𝑒

𝜕𝑝

𝜕𝑦= 0 (3.80)

Figura 3.14 – Condição de contorno imposta nas fronteiras de freestream.

Condição de contorno imposta na fronteira de outflow

A fronteira de outflow permite a saída de massa (derivada nula). Esta condição

de contorno pode ser escrita matematicamente da seguinte forma:

𝜕𝑢

𝜕𝑥= 0,

𝜕𝑣

𝜕𝑥= 0 𝑒

𝜕𝑝

𝜕𝑥= 0 (3.81)

24

Figura 3.15 – Condição de contorno imposta na fronteira de saída do domínio computacional.

Condição de contorno imposta na superfície do corpo

Na parede do cilindro, aplicam-se as condições de contorno de não-

escorregamento e impenetrabilidade com velocidade nula. Esta condição de contorno

pode ser escrita matematicamente da seguinte forma:

𝑢 = 0, 𝑣 = 0 𝑒 𝜕𝑝

𝜕𝑛= 0 (3.82)

Figura 3.16 – Condição de contorno imposta no corpo de interesse (seção circular do cilindro).

25

4 Estudo do Efeito dos Limites do Domínio Computacional no

Campo de Escoamento

“Uma das fases mais importantes na simulação numérica, é a definição do

domínio computacional, que nada mais é do que uma aproximação do domínio real

onde o problema físico estudado é observado. Muitos tipos de condições de contorno

utilizadas em aplicações práticas são apenas aplicáveis se for garantido que estas estão

suficientemente removidas da região em que a precisão da solução é importante. O

desejo de limitar o custo computacional muitas das vezes proporcionam uma motivação

de reduzir o domínio computacional. E o que muitas das vezes não é observado, é que

os melhores resultados só podem ser alcançados se estas tendências totalmente

contraditórias forem perfeitamente equilibradas.”, tradução livre de Behr et al. (1995).

É neste contexto, que o presente estudo tem como objetivo encontrar um

domínio computacional dito ótimo, onde as influências das fronteiras destes, ou das

condições de contorno impostas nestas, estão suficientemente removidas da região em

que a precisão da solução é importante. E posteriormente utilizar este domínio

computacional para estudar os efeitos de diferentes métodos de marcha no tempo no

campo de escoamento. Para isso, será utilizado como benchmark os resultados da

literatura obtidos para o escoamento entorno de uma seção circular de um cilindro fixo

com número de Reynolds igual à 40. Para o presente estudo, será utilizado o método

explícito de Runge-Kutta de 2ª ordem “The Midpoint Method”, como método de marcha

no tempo.

Para o estudo em questão, foram seguidas as etapas:

1. Definição de um domínio computacional de referência:

1.1. Aqui será definido um domínio computacional que servirá de referência ou

ponto de partida para o estudo da influência dos limites deste sob o campo de

escoamento.

2. Definição das características da malha computacional a ser gerada:

2.1. Aqui serão definidas as características das malhas a serem geradas.

3. Estudo da influência da localização das fronteiras laterais (condição de contorno de

“freestream”) no campo de escoamento:

26

3.1. O estudo em questão resume-se na análise da influência da localização das

fronteiras laterais (condição de contorno de “freestream”) no campo de

escoamento.

4. Estudo da influência da localização da fronteira de entrada (condição de contorno de

“inflow”) no campo de escoamento:

4.1. O estudo em questão resume-se na análise da influência da localização da

fronteira de entrada (condição de contorno de “inflow”) no campo de

escoamento.

5. Estudo da influência da localização da fronteira de saída (condição de contorno de

“outflow”) no campo de escoamento:

5.1. O estudo em questão resume-se na análise da influência da localização da

fronteira de saída (condição de contorno de “outflow”) no campo de

escoamento.

6. Estudo do efeito do refinamento da malha no campo de escoamento:

6.1. O estudo em questão resume-se na análise da influência do refinamento da

malha, na região próxima à superfície da seção circular do cilindro, no campo

de escoamento.

4.1 Definição de um Domínio Computacional de Referência

Tendo como referência Behr et al. (1995) e Behr et al. (1991), foi definido aqui

um domínio computacional de partida, conforme Figura 4.1. Domínio este, que será

utilizado como referência para a construção de todos os domínios computacionais

utilizados nos estudos subsequentes. Este domínio computacional pode ser visualizado

na Figura 4.2, logo abaixo.

Figura 4.1 – Domínio computacional de referência.

27

As características geométricas deste domínio computacional podem ser vistas na

Tabela 1 a seguir.

Tabela 1 – Características geométricas do domínio computacional de referência.

Características Geométricas do Domínio Computacional de Referência

A [ ] B [ ] C [ ] D [ ] E [ ] F [ ] G [ ]

𝟒 ∗ 𝑫 𝟑 ∗ 𝑫 𝟏𝟎 ∗𝝅 ∗ 𝑫

𝟏𝟎𝟎 1.00 𝟗 ∗ 𝑫 𝟖. 𝟓 ∗ 𝑫 𝟐𝟑. 𝟎 ∗ 𝑫

onde:

A (A+B):

o É o comprimento da Região 2;

B:

o É a largura da Região 2;

C:

o É a largura da Região 1;

D:

o É o diâmetro da seção circular do cilindro. Que para o presente

estudo sempre será 1.0;

E:

o É a largura do domínio computacional ou duas vezes a distância

entre o centro da seção circular do cilindro e as fronteiras laterais

(condição de contorno de “freestream”) do domínio

computacional;

F:

o É a distância entre o centro da seção circular do cilindro e a

fronteira de entrada (condição de contorno de “inflow”) do

domínio computacional;

G:

o É a distância entre o centro da seção circular do cilindro e a

fronteira de saída (condição de contorno de “outflow”) do

domínio computacional;

Região 1:

28

o É a região adjacente à parede da seção circular do cilindro, gerada

com o intuito de capturar a camada limite no entorno desta

superfície.

Região 2:

o É a região no entorno da seção circular do cilindro e da região 1,

gerada com o intuito de capturar a esteira de vórtices de von

Kárman gerada por esta seção.

Região 3:

o É a região restante do domínio computacional.

29

(a)

(b)

(c)

(d)

Figura 4.2 – Domínio computacional de referência gerado no software Gmsh.

Ao se construir os diversos domínios computacionais, o mesmo domínio de

referência é utilizado como ponto de partida, adicionando apenas regiões a este,

conforme o caso ilustrado na Figura abaixo.

30

(a)

(b)

(c)

(d)

(e)

(f)

Figura 4.3 – Construção dos diversos domínios computacionais utilizados no estudo do efeito da localização das

fronteiras laterais no campo de escoamento. (a) Domínio de referência (Domínio 1); (b) Domínio 2; (c) Domínio 3;

(d) Domínio 4; (e) Domínio 5 e (f) Domínio 6.

Todos os domínios computacionais foram gerados com esta sistemática com o

intuito de remover qualquer influência de refinamento da malha devido à presença da

seção circular do cilindro, da região 1 ou 2. Ou seja, cada domínio computacional

imediatamente maior, é construído a partir do domínio computacional imediatamente

31

menor, no limite em que cada domínio computacional menor poderia ser visto como um

truncamento do domínio computacional imediatamente maior. Assim, tenta-se garantir

que a malha na região de interesse, malha computacional 1 (Malha 1), não sofra grandes

modificações de uma malha para outra. Com isso, tenta-se eliminar neste estágio a

influência do refinamento da malha nos resultados, deixando este estudo para a última

análise a ser realizada antes do estudo dos métodos de marcha no tempo.

Um exemplo pode ser visto na Figura 4.4, onde temos as imagens dos domínios

computacionais e malhas 1 e 6. Como se pode observar, a malha 6 não foi gerada

através da sistemática explicada anteriormente, e com isso salta aos olhos a diferença

entre elas. Podendo se observar a “olho nu”, que nas regiões de interesse, regiões 1 e 2,

não se consegui diferir entre uma malha e outra. No entanto, na malha no entorno da

região 2, fica evidente a diferença entre estas. O que não se consegue afirmar, quando a

sistemática anteriormente explicada é aplicada. O que é o caso demonstrado pela Figura

4.5, onde poderia se afirmar, que quando se reduz a malha 6 a malha 1, a “olho nu” não

se conseguiria diferir entre uma malha e outra. Justificando assim, a sistemática

utilizada e apresentada anteriormente para a geração dos domínios, e consequentemente

das malhas computacionais utilizadas nos estudos subsequentes.

32

(a)

(d)

(B)

(e)

(c)

(f)

Figura 4.4 – Comparação entre as malhas computacionais 1 e 6. Onde fica evidente a diferença entra as malhas

computacionais quando não se aplica a sistemática de geração dos domínios computacionais e consequentemente das

malhas anteriormente explicada.

33

(a)

(e)

(c)

(f)

(d)

(g)

Figura 4.5 – Comparação entre as malhas computacionais 1 e 6. Agora, ficando difícil de a “olho nu” diferir entre

uma malha e outra quando se aplica a sistemática de geração dos domínios computacionais e consequentemente das

malhas anteriormente explicada.

34

4.2 Definição das Características da Malha Computacional a ser Gerada

Todas as malhas computacionais foram geradas através do software Gmsh. No

Gmsh, as características dimensionais da malha foram aplicadas nas fronteiras das

regiões modeladas, conforme Figura 4.6.

Figura 4.6 – Características da malha computacional.

Desta forma as características da malha computacional são as seguintes, segundo

as regiões modeladas (Figura 4.2):

Tabela 2 – Características da malha computacional segundo as regiões modeladas (1, 2 ou 3).

Regiões Modeladas Tamanho Característico da Malha

Computacional Presente na Região

Tamanho Característico da Malha

Computacional

Região 1 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 0.031

Região 2 𝟖 ∗ 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 0.251

Região 3 𝟑 ∗ 𝑫𝟓⁄ 0.600

A superfície da seção circular do cilindro é composta por 100 elementos cujos

comprimentos de suas faces são iguais (𝜋 ∗ 𝐷 100⁄ ). A região 1 é caracterizada por

elementos que possuem as mesmas características dos elementos presentes na superfície

da seção circular do cilindro, ou seja, possuem um tamanho característico de suas faces

de 0.031, em média. Enquanto que a região 2, é composta por elementos cujo

comprimento de suas faces são 8 (oito) vezes maior que o comprimento das faces dos

elementos presentes na região 1. Já a região 3, é composta por elementos cujo

comprimento de suas faces é 0.600, em média. Ressaltando, que estas características são

replicadas para todas as malhas aqui utilizadas. E deixando claro que em todas as

35

malhas geradas existem apenas 3 (três) regiões, conforme explicado no subcapítulo

anterior, cujas as características da malha estão presentes na Tabela acima.

4.3 Estudo da Influência da Localização das Fronteiras Laterais (Condição de

Contorno de “Freestream”) no Campo de Escoamento

Para o estudo do efeito da localização das fronteiras laterais no campo de

escoamento, foram utilizadas seis malhas computacionais não estruturadas, conforme

dados fornecidos na Tabela abaixo. As distâncias das regiões de entrada e saída,

também chamadas de inflow e outflow, respectivamente, da seção circular do cilindro,

foram mantidas fixas. O centro desta está localizado a uma distância E/2 das fronteiras

laterais do domínio computacional, e será esta variável o objeto de estudo da presente

seção. Com este propósito, foi explorado um intervalo de E que vai de 9.00 ∗ 𝐷 à

31.50 ∗ 𝐷. Estas malhas serão referenciadas pelos símbolos dados na Tabela 3 logo

abaixo.

Tabela 3 – Características das malhas utilizadas no estudo do efeito da localização das fronteiras laterais no campo de

escoamento.

Malha Computacional

Dimensões do Domínio Computacional

Número de Elementos Número de Nós

E F G

Malha 1 9.00 8.50 23.00 8254 4247

Malha 2 13.50 8.50 23.00 9162 4709

Malha 3 18.00 8.50 23.00 10068 5170

Malha 4 22.50 8.50 23.00 11000 5644

Malha 5 27.00 8.50 23.00 11858 6081

Malha 6 31.50 8.50 23.00 12786 6553

Os domínios e as malhas computacionais utilizadas no presente estudo podem

ser visualizadas nas Figuras abaixo.

36

(a)

(b)

(c)

(d)

(e)

(f)

Figura 4.7 – Domínios computacionais utilizados no estudo do efeito da localização das fronteiras laterais no campo

de escoamento. (a) Domínio de referência (Domínio 1); (b) Domínio 2; (c) Domínio 3; (d) Domínio 4; (e) Domínio 5

e (f) Domínio 6.

37

(a)

(b)

(c)

(d)

(e)

(f)

Figura 4.8 – Malhas computacionais utilizadas no estudo do efeito da localização das fronteiras laterais no campo de

escoamento. (a) Malha 1; (b) Malha 2; (c) Malha 3; (d) Malha 4; (e) Malha 5 e (f) Malha 6.

Os resultados do coeficiente de arrasto médio podem ser vistos na Tabela 4 e na

Figura 4.9.

38

Tabela 4 – Variação dos valores do coeficiente de arrasto médio (CD) com a distância E.

Malha 𝑬 𝑪𝑫 [ ]

Malha 1 9.00 1.760

Malha 2 13.50 1.707

Malha 3 18.00 1.679

Malha 4 22.50 1.666

Malha 5 27.00 1.658

Malha 6 31.50 1.656

Figura 4.9 – Variação dos valores do coeficiente de arrasto médio com a distância E.

Pode ser claramente observado na Figura 4.9, que a localização dos limites

laterais do domínio computacional afetam os valores do coeficiente de arrasto médio, e

por consequência também as propriedades do campo de escoamento. Quando a distância

E é tão pequena quanto 9.00 ∗ 𝐷, podemos observar um aumento de 6.28% no valor do

coeficiente de arrasto médio. É evidente também a partir da mesma Figura, que estes

1,65

1,66

1,67

1,68

1,69

1,70

1,71

1,72

1,73

1,74

1,75

1,76

1,77

9,00 13,50 18,00 22,50 27,00 31,50 36,00

Co

efi

cie

nte

de

Arr

asto

(C D

)

Distância E

39

valores estabilizam dentro dos limites de 0.60% apenas quando a distância E excede

22.50 ∗ 𝐷.

Com isso, chega-se à conclusão que as fronteiras laterais devem estar a uma

distância mínima da seção circular do cilindro de 11.25 ∗ 𝐷. Caso isso não seja

respeitado, não apenas o valor do coeficiente de arrasto médio será artificialmente

aumentado, mais também serão modificadas as características do escoamento no

entorno da seção circular do cilindro. Então, o domínio computacional escolhido como

referência para dar prosseguimento aos estudos de influência da localização das

fronteiras no campo de escoamento, foi o 6 (Domínio 6).

4.4 Estudo da Influência da Localização da Fronteira de Entrada (Condição de

Contorno de “Inflow”) no Campo de Escoamento

Seguindo a mesma sistemática apresentada no capítulo anterior, para o estudo do

efeito da localização da fronteira de entrada no campo de escoamento, foram utilizadas

seis malhas computacionais não estruturadas, conforme dados fornecidos na Tabela

abaixo. Sendo que o domínio computacional 6, é o domínio de referência, ou seja, todos

os outros cinco domínios computacionais serão construídos tendo como referência o

domínio 6 (Domínio 6). As distâncias das fronteiras laterais e da fronteira de saída da

seção circular do cilindro, foram mantidas fixas. O centro desta seção circular está

localizado a uma distância F da fronteira de entrada do domínio computacional, e agora

será esta variável o objeto de estudo do presente subcapítulo. Com este propósito, foi

explorado um intervalo de F que vai de 8.00 ∗ 𝐷 à 23.00 ∗ 𝐷. As malhas serão

referenciadas pelos símbolos dados na Tabela 5 logo abaixo.

40

Tabela 5 – Características das malhas utilizadas no estudo do efeito da localização da fronteira de entrada no campo

de escoamento.

Malha Computacional

Dimensões do Domínio Computacional

Número de Elementos Número de Nós

E F G

Malha 6 31.50 8.50 23.00 12786 6553

Malha 6_1 31.50 11.40 23.00 13454 6891

Malha 6_2 31.50 14.30 23.00 14060 7199

Malha 6_3 31.50 17.20 23.00 14716 7532

Malha 6_4 31.50 20.10 23.00 15378 7868

Malha 6_5 31.50 23.00 23.00 16006 8187

As condições de contorno impostas nas fronteiras do domínio computacional e o

tempo adimensional, são os mesmos que os utilizados no capítulo anterior.

Os domínios e as malhas computacionais utilizadas no presente estudo, podem

ser visualizadas nas Figuras abaixo.

41

(a)

(d)

(b)

(e)

(c)

(f)

Figura 4.10 – Domínios computacionais utilizados no estudo do efeito da localização da fronteira de entrada no

campo de escoamento. (a) Domínio de referência (Domínio 6); (b) Domínio 6_1; (c) Domínio 6_2; (d) Domínio 6_3;

(e) Domínio 6_4 e (f) Domínio 6_5.

42

(a)

(d)

(b)

(e)

(c)

(f)

Figura 4.11 – Malhas computacionais utilizadas no estudo do efeito da localização da fronteira de entrada no campo

de escoamento. (a) Malha 6; (b) Malha 6_1; (c) Malha 6_2; (d) Malha 6_3; (e) Malha 6_4 e (f) Malha 6_5.

Os resultados do coeficiente de arrasto médio podem ser vistos na Tabela 6 e na

Figura 4.12.

43

Tabela 6 – Variação dos valores do coeficiente de arrasto médio (CD) com a distância F.

Malha 𝑭 𝑪𝑫 [ ]

Malha 6 8.50 1.656

Malha 6_1 11.40 1.633

Malha 6_2 14.30 1.621

Malha 6_3 17.20 1.614

Malha 6_4 20.10 1.612

Malha 6_5 23.00 1.609

Figura 4.12 – Variação do valor médio do coeficiente de arrasto com a distância F.

Analisando-se os resultados do coeficiente de arrasto médio apresentados na

Figura 4.12, conclui-se que a localização da fronteira de entrada do domínio

computacional também afeta os valores do coeficiente de arrasto médio, e por

consequência também as propriedades do campo de escoamento. Quando a distância F é

tão pequena quanto 8.50 ∗ 𝐷, observar-se um aumento de 2.92% no valor do coeficiente

1,60

1,61

1,62

1,63

1,64

1,65

1,66

1,67

1,68

1,69

1,70

1,71

1,72

8,50 13,00 17,50 22,00 26,50 31,00 35,50

Co

efi

cie

nte

de

Arr

asto

(C D

)

Distância F

44

de arrasto médio. Sendo evidente também a partir da mesma Figura, que estes valores

estabilizam dentro dos limites de 0.31%, apenas quando a distância F excede 17.20 ∗ 𝐷.

Com isso, chega-se à conclusão que mais uma vez, a região de entrada deve

estar a uma distância mínima da seção circular do cilindro de 17.20 ∗ 𝐷. Caso isso não

seja respeitado, não apenas o valor do coeficiente de arrasto médio será artificialmente

aumentado, mais também as características do escoamento no entorno da seção circular

do cilindro serão modificadas. Então, o domínio computacional escolhido como

referência para dar prosseguimento aos estudos de influência da localização das

fronteiras no campo de escoamento, foi o 6_5 (Domínio 6_5).

4.5 Estudo da Influência da Localização da Fronteira de Saída (Condição de

Contorno de “Outflow”) no Campo de Escoamento

Seguindo a mesma sistemática apresentada nas seções anteriores, para o estudo

do efeito da localização da fronteira de saída no campo de escoamento, foram também

utilizadas seis malhas computacionais não estruturadas, conforme dados fornecidos na

Tabela abaixo. Sendo que o domínio computacional 6_5, é o domínio de referência, ou

seja, todos os outros cinco domínios computacionais serão construídos tendo este como

referência. As distâncias das fronteiras laterais e de entrada da seção circular do

cilindro, foram mantidas fixas. O centro desta seção está localizado a uma distância G

da fronteira de saída do domínio computacional, e agora será esta variável o objeto de

estudo da presente seção. Com este propósito, nós exploramos um intervalo de G que

vai de 23.00 ∗ 𝐷 à 37.50 ∗ 𝐷. Estas malhas serão referenciadas pelos símbolos dados na

Tabela 7 logo abaixo.

45

Tabela 7 – Características das malhas utilizadas no estudo do efeito da localização da fronteira de saída no campo de

escoamento.

Malha Computacional

Dimensões do Domínio Computacional

Número de Elementos Número de Nós

E F G

Malha 6_5 31.50 23.00 23.00 16006 8187

Malha 6_5_1 31.50 23.00 25.90 16650 8513

Malha 6_5_2 31.50 23.00 28.80 17292 8839

Malha 6_5_3 31.50 23.00 31.70 17936 9166

Malha 6_5_4 31.50 23.00 34.60 18588 9497

Malha 6_5_5 31.50 23.00 37.50 19210 9813

As condições de contorno impostas nas fronteiras do domínio computacional e o

tempo adimensional, também são os mesmos que os utilizados nos capítulos anteriores.

Os domínios e as malhas computacionais utilizadas no presente estudo, podem

ser visualizadas nas Figuras abaixo.

46

(a)

(d)

(b)

(e)

(c)

(f)

Figura 4.13 – Domínios computacionais utilizados no estudo do efeito da localização da fronteira de saída no campo

de escoamento. (a) Domínio de referência (Domínio 6_5); (b) Domínio 6_5_1; (c) Domínio 6_5_2; (d) Domínio

6_5_3; (e) Domínio 6_5_4 e (f) Domínio 6_5_5.

47

(a)

(d)

(b)

(e)

(c)

(f)

Figura 4.14 – Malhas computacionais utilizadas no estudo do efeito da localização da fronteira de saída no campo de

escoamento. (a) Malha 6_5; (b) Malha 6_5_1; (c) Malha 6_5_2; (d) Malha 6_5_3; (e) Malha 6_5_4 e (f) Malha

6_5_5.

Os resultados do coeficiente de arrasto médio podem ser vistos na Tabela 8 e na

Figura 4.15.

48

Tabela 8 – Variação dos valores do coeficiente de arrasto médio (CD) com a distância G.

Malha 𝑮 𝑪𝑫 [ ]

Malha 6_5 23.00 1.609

Malha 6_5_1 25.90 1.607

Malha 6_5_2 28.80 1.606

Malha 6_5_3 31.70 1.605

Malha 6_5_4 34.60 1.603

Malha 6_5_5 37.50 1.604

Figura 4.15 – Variação do valor do coeficiente de arrasto médio com a distância G.

Pode ser claramente observado na Figura 4.15, que a localização da fronteira de

saída do domínio computacional afeta muito pouco os valores do coeficiente de arrasto

médio, e por consequência também as propriedades do campo de escoamento. Quando a

distância G é tão maior quanto 23.00 ∗ 𝐷, podemos observar uma queda apenas de

0.31% no valor do coeficiente de arrasto médio.

1,60

1,61

1,62

1,63

1,64

1,65

1,66

1,67

1,68

1,69

1,70

1,71

1,72

23,00 27,50 32,00 36,50 41,00 45,50 50,00

Co

efi

cie

nte

de

Arr

asto

(C D

)

Distância G

49

Com isso, chega-se à conclusão que sob as condições impostas aqui, a região de

saída pode permanecer na posição inicialmente estimada, ou seja, a uma distância de

23.00 ∗ 𝐷 da seção circular do cilindro. A escolha de uma localização mais afastada não

traria maiores benefícios, pois além da diferença de um resultado para o outro não ser

tão significativa assim, vide resultados acima, o tempo computacional consumido seria

maior. Então, o domínio computacional escolhido como referência para o estudo da

influência do refinamento da malha no campo de escoamento foi o 6_5 (Domínio 6_5).

4.6 Estudo do Efeito do Refinamento da Malha no Campo de Escoamento

Uma vez realizado os estudos de influência da localização dos limites do

domínio computacional no campo de escoamento e escolhido o domínio dito ótimo, na

presente seção será analisada o efeito do refinamento da malha no campo de

escoamento. Todavia, este estudo será localizado, ou seja, será realizado um

refinamento da malha apenas na região 2. E este refinamento é ditado pelo

comprimento característico da malha definido na fronteira que delimita esta região,

conforme explicado no subcapítulo 4.2. Sendo assim, agora a variável objeto de estudo

do presente subcapítulo será este tamanho característico. Com este propósito, nós

exploramos um intervalo deste tamanho característico que vai de 8.0 ∗ 𝐷 𝜋⁄ à

1.0 ∗ 𝐷 𝜋⁄ . As malhas obtidas através deste refinamento serão referenciadas pelos

símbolos dados na Tabela 9 logo abaixo.

Tabela 9 – Características das malhas utilizadas no estudo do efeito do refinamento destas no campo de escoamento.

Malha

Computacional

Dimensões do Domínio

Computacional

Tamanho

Característico da

Malha

Computacional na

Região 2

Número de

Elementos

Número de

Nós

E F G

Malha 6_5 31.50 23.00 23.00 𝟖. 𝟎 ∗ 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 16006 8187

Malha 6_5_A 31.50 23.00 23.00 𝟔. 𝟔 ∗ 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 16796 8582

Malha 6_5_B 31.50 23.00 23.00 𝟓. 𝟐 ∗ 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 17980 9174

Malha 6_5_C 31.50 23.00 23.00 𝟑. 𝟖 ∗ 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 20492 10430

Malha 6_5_D 31.50 23.00 23.00 𝟐. 𝟒 ∗ 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 26986 13677

Malha 6_5_E 31.50 23.00 23.00 𝟏. 𝟎 ∗ 𝝅 ∗ 𝑫𝟏𝟎𝟎⁄ 60960 30664

50

As condições de contorno impostas nas fronteiras do domínio computacional e o

tempo adimensional são os mesmos que os utilizados nos estudos anteriores.

As malhas computacionais utilizadas no presente estudo podem ser visualizadas

nas Figuras abaixo.

Figura 4.16 – Malha computacional 6_5.

Figura 4.17 – Malha computacional 6_5_A.

Figura 4.18 – Malha computacional 6_5_B.

51

Figura 4.19 – Malha computacional 6_5_C.

Figura 4.20 – Malha computacional 6_5_D.

Figura 4.21 – Malha computacional 6_5_E.

52

Os resultados do coeficiente de arrasto médio podem ser vistos na Tabela 10 e

na Figura 4.22 abaixo.

Tabela 10 – Variação dos valores do coeficiente de arrasto médio (CD) com o refinamento da malha.

Malha 𝑪𝑫 [ ]

Malha 6_5 1.609

Malha 6_5_A 1.597

Malha 6_5_B 1.587

Malha 6_5_C 1.580

Malha 6_5_D 1.576

Malha 6_5_E 1.574

Figura 4.22 – Variação do coeficiente de arrasto médio com o refinamento da malha computacional na região 2.

Pode ser claramente observado na Figura 4.22, que o refinamento da malha afeta

os valores do coeficiente de arrasto médio, e por consequência também as propriedades

do campo de escoamento. Quando o tamanho do elemento característico na região 2 é

tão maior quanto 0.251, observar-se um aumento de 2.22% no valor do coeficiente de

1,510

1,520

1,530

1,540

1,550

1,560

1,570

1,580

1,590

1,600

1,610

1,620

15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000

Co

efi

cie

nte

de

Arr

asto

(C

D)

Número de Elementos

53

arrasto médio. Sendo evidente também a partir da mesma Figura, que estes valores

estabilizam dentro dos limites de 0.13%, apenas quando o tamanho característico do

elemento é tão menor quanto 0.075.

Com isso, chega-se à conclusão que o resultado obtido pela Malha

Computacional 6_5_D é o ótimo após o presente estudo. Pois, como pode-se observar,

da Malha 6_5_D para a 6_5_E a quantidade de elementos mais que duplicou, porém não

houve um ganho significativo no coeficiente de arrasto médio que justificasse a

utilização de uma malha com 60960 elementos. Assim, chegamos tanto ao Domínio

quanto a Malha computacional ditos ótimos (6_5_D) para iniciarmos os estudos do

efeito dos métodos de marcha no tempo no campo de escoamento.

5 Estudo do Efeito dos Métodos de Marcha no Tempo no Campo

de Escoamento

Uma vez obtido tanto o domínio quanto a malha computacional ditos ótimos, o

próximo passo será a análise do efeito dos métodos de marcha no tempo no campo de

escoamento. Então, aqui serão apresentados, comparados e analisados os diversos

resultados obtidos através dos métodos de marcha no tempo para os números de

Reynolds 40, 100 e 200.

5.1 Resultados para o Número de Reynolds 40

Inicialmente são apresentados e comparados os valores do coeficiente de arrasto

(𝐶𝐷) e sustentação (𝐶𝐿) obtidos numericamente com os valores numéricos e

experimentais encontrados na literatura técnica para estudar tanto a influência quanto a

precisão da solução de cada método de marcha no tempo para este regime de

escoamento, Reynolds 40. Estes resultados podem ser visualizados na Tabela 11, logo

abaixo.

54

Tabela 11 – Comparação dos coeficientes de arrasto (𝐶𝐷), sustentação (𝐶𝐿) e das características geométricas do par de

vórtices obtidos através dos diversos métodos de marcha no tempo para o número de Reynolds 40.

Referência CD 𝑳𝑫⁄ 𝒂

𝑫⁄ 𝒃𝑫⁄ 𝜽𝒔 Δt Tipo de Estudo

Tritton (1959) 1.57 - - - - - Experimental

Constanceau (1977) - 2.13 0.76 0.59 53.50 - Experimental

Rengel (1999) 1.61 2.23 0.72 0.58 54.06 - Numérico (FVM)

Wanderley (2008a) 1.56 2.29 0.73 0.60 53.8 - Numérico (FDM)

Método de Euler

(Presente Estudo)

1.554 2.27 0.73 0.60 52.34 0.002 Numérico (FVM)

“Heun Method with a Single Corrector”

(Presente Estudo)

1.732 2.25 0.72 0.59 51.51 0.005 Numérico (FVM)

“The Midpoint Method”

(Presente Estudo)

1.576 2.28 0.73 0.60 52.23 0.005 Numérico (FVM)

“Ralston’s Method”

(Presente Estudo)

1.805 2.23 0.71 0.59 51.21 0.005 Numérico (FVM)

Figura 5.1 – Dimensões características do par de vórtices presentes na esteira de um cilindro para um escoamento

incidente com número de Reynolds 40.

Pode-se observar uma diferença significativa entre os métodos de marcha no

tempo de Ralston e de Heun quando comparados com os valores do coeficiente de

arrasto (𝐶𝐷) presentes na Tabela 11. Apesar destes métodos representaram a física deste

escoamento com certa fidelidade, conforme Figura 5.2, estes se distanciaram e muito

dos resultados numéricos e experimentais encontrados na literatura técnica. Com isso,

55

estes métodos mostram-se incapazes de serem utilizados como métodos de marcha no

tempo para o estudo numérico do escoamento ao redor de um cilindro para número de

Reynolds 40. Todavia, os métodos de Euler e “The Midpoint”, mostraram-se muito

eficientes quando comparados com os demais resultados encontrados na literatura.

Na Tabela 11, também são apresentadas algumas dimensões características dos

dois vórtices gerados na esteira do cilindro para o número de Reynolds 40. A Figura 5.1

apresenta a definição destas dimensões características contidas na Tabela 11. A

concordância, mais uma vez, entre o coeficiente de arrasto obtido através dos métodos

de Euler e “The Midpoint” com os resultados numéricos de Rengel (1999) e Wanderley

(2008a), e com o resultado experimental obtido por Tritton (1957) são satisfatórios. As

dimensões características dos dois vórtices obtidos por estes métodos também

concordam bem com os obtidos experimentalmente por Constanceau (1977), e

numericamente com Rengel (1999) e Wanderley (2008a), conforme Tabela 11 e Figura

5.2.

As Figuras Figura 5.3 e Figura 5.4 mostram a variação dos coeficientes de

arrasto (𝐶𝐷) e de sustentação (𝐶𝐿) em função do tempo para todos os métodos de marcha

no tempo, respectivamente. Como esperado, o coeficiente de sustentação (𝐶𝐿) para

Reynolds 40 é zero, enquanto que o coeficiente de arrasto obtido através dos métodos

de Euler e “The Midpoint” possuem uma grande concordância com os resultados

numéricos e experimentais apresentados na Tabela 11.

56

(a)

(c)

(b)

(d)

Figura 5.2 – Campo de pressão e streamlines para o número de Reynolds 40. (a) Método de Euler; (b) Método de

Heun; (c) Método “The Midpoint”; e (d) Método de Ralston.

Figura 5.3 – Variação no tempo do coeficiente de arrasto (𝐶𝐷) para os diversos métodos de marcha no tempo

(Reynolds 40).

57

Figura 5.4 – Variação no tempo do coeficiente de sustentação (𝐶𝐿) para os diversos métodos de marcha no tempo

(Reynolds 40).

5.2 Resultados para o Número de Reynolds 100

Aqui também são apresentados e comparados os valores do coeficiente de

arrasto (𝐶𝐷), sustentação (𝐶𝐿) e agora também o número de Strouhal (𝑆𝑡) obtidos

numericamente com os valores numéricos e experimentais encontrados na literatura

técnica para estudar tanto a influência quanto a precisão da solução de cada método de

marcha no tempo para este regime de escoamento, Reynolds 100. Estes resultados

podem ser visualizados na Tabela 12, logo abaixo.

58

Tabela 12 – Comparação dos coeficientes de arrasto (𝐶𝐷), sustentação (𝐶𝐿) e número de Strouhal (𝑆𝑡) obtidos através

dos diversos métodos de marcha no tempo para número de Reynolds 100.

Referência CD CL St Δt Tipo de Estudo

Herfjord (1995) 1.36 0.34 0.168 - Numérico (FEM)

Rengel (1999) 1.36 0.32 0.173 - Numérico (FVM)

Norberg (2003) - 0.32 0.164 - Experimental

Wieselsberger (1921) 1.41 - - - Experimental

Wanderley (2008a) 1.30 0.25 0.158 - Numérico (FDM)

Método de Euler

(Presente Estudo)

1.311 0.334 0.16 0.001 Numérico (FVM)

“Heun Method with a Single Corrector”

(Presente Estudo)

1.404 0.336 0.16 0.002 Numérico (FVM)

“The Midpoint Method”

(Presente Estudo)

1.345 0.333 0.16 0.002 Numérico (FVM)

“Ralston’s Method”

(Presente Estudo)

1.428 0.336 0.16 0.002 Numérico (FVM)

Pode-se observar que agora não existe uma diferença tão grande entre os

métodos de marcha no tempo de Ralston e de Heun quando comparados com os valores

do coeficiente de arrasto (𝐶𝐷) presentes na Tabela 12. Porém, estes valores ainda são

altos, mesmo apesar destes métodos representaram a física deste escoamento com

fidelidade, conforme Figura 5.5. Estes valores se distanciaram um pouco da média dos

resultados numéricos encontrados na literatura técnica. Mostrando-se desta vez, capazes

de serem utilizados como métodos de marcha no tempo para o estudo numérico do

escoamento ao redor de um cilindro para número de Reynolds 100. Todavia, os métodos

de Euler e “The Midpoint”, continuam mostrando-se mais eficientes do que estes

quando comparados com os demais resultados numéricos encontrados na literatura.

O campo de pressões com as linhas de corrente podem ser visualizadas na Figura

5.5 para os diferentes métodos de marcha no tempo. Estes indicam o comportamento do

fluido para o presente regime de escoamento. Como esperado, existe um padrão

oscilante das linhas de corrente causado pelo descolamento dos vórtices.

59

As Figuras Figura 5.6 e Figura 5.7 mostram a variação dos coeficientes de

arrasto (𝐶𝐷) e de sustentação (𝐶𝐿) em função do tempo para todos os métodos de marcha

no tempo, respectivamente. Como observado, o coeficiente de sustentação (𝐶𝐿) para

Reynolds 100 não encontrou maiores diferenças entre um método e outro, possuindo

uma satisfatória concordância com todos os resultados numéricos e experimentais

apresentados na Tabela 12. Enquanto que o coeficiente de arrasto obtido através dos

métodos de Euler e “The Midpoint” possuem uma grande concordância com os

resultados numéricos de Herfjord (1995), Rengel (1999) e Wanderley (2008a).

(a)

(c)

(b)

(d)

Figura 5.5 – Coeficiente de pressão e streamlines para o número de Reynolds 100. (a) Método de Euler; (b) Método

de Heun; (c) Método “The Midpoint”; e (d) Método de Ralston.

60

Figura 5.6 – Comparação do coeficiente de arrasto (𝐶𝐷) para os diversos métodos de marcha no tempo (Reynolds

100).

Figura 5.7 – Comparação do coeficiente de sustentação (𝐶𝐿) para os diversos métodos de marcha no tempo (Reynolds

100).

A Figura 5.8 é o resultado de uma análise espectral dos coeficientes de

sustentação (𝐶𝐿), Figura 5.7. O número de Strouhal é função da frequência de liberação

de vórtices (𝑓), do diâmetro da seção circular do cilindro (𝐷) e da velocidade do

escoamento (𝑈∞), 𝑆𝑡 =𝑓𝐷

𝑈∞. Na Figura 5.8, estão mostradas as oscilações

adimensionalizadas do coeficiente de sustentação no domínio da frequência. Onde

61

podemos observar, que todos os métodos de marcha no tempo capturam a mesma

frequência de desprendimento de vórtices. Estes também concordaram de forma

satisfatória com os resultados numéricos e experimentais presentes na Tabela 12.

(a)

(c)

(b)

(d)

Figura 5.8 – Número de Strouhal (𝑆𝑡) para o número de Reynolds 100. (a) Método de Euler; (b) Método de Heun; (c)

Método “The Midpoint”; e (d) Método de Ralston.

5.3 Resultados para o Número de Reynolds 200

De forma exatamente análoga à seção anterior, aqui também são apresentados e

comparados os valores do coeficiente de arrasto (𝐶𝐷), sustentação (𝐶𝐿) e o número de

Strouhal (𝑆𝑡) obtidos numericamente com os valores numéricos e experimentais obtidos

na literatura técnica. Esta comparação tem como objetivo, mais uma vez, estudar tanto a

influência quanto a precisão da solução de cada método de marcha no tempo para este

regime de escoamento, Reynolds 200. Estes resultados podem ser visualizados na

Tabela 13, logo abaixo.

62

Tabela 13 – Comparação dos coeficientes de arrasto (𝐶𝐷), sustentação (𝐶𝐿) e número de Strouhal (𝑆𝑡) obtidos através

dos diversos métodos de marcha no tempo para o número de Reynolds 200.

Referência CD CL St Δt Tipo de Estudo

Herfjord (1995) 1.35 0.70 0.196 0.002

Numérico

(FEM)

Rengel (1999) 1.35 0.67 0.203 - Numérico (FVM)

Norberg (2003) - 0.53 0.182 - Experimental

Wieselsberger (1921) 1.29 - - - Experimental

Wanderley (2008) 1.27 0.51 0.187 - Numérico (FDM)

Método de Euler

(Presente Estudo)

1.301 0.673 0.185 0.0005 Numérico (FVM)

“Heun Method with a Single Corrector”

(Presente Estudo)

1.340 0.681 0.185 0.0010 Numérico (FVM)

“The Midpoint Method”

(Presente Estudo)

1.325 0.680 0.185 0.0010 Numérico (FVM)

“Ralston’s Method”

(Presente Estudo)

1.345 0.681 0.185 0.0010 Numérico (FVM)

Pode-se observar que a diferença diminui ainda mais com o aumento do número

de Reynolds entre os métodos de marcha no tempo de Ralston e de Heun quando

comparados com os valores do coeficiente de arrasto (𝐶𝐷) presentes na literatura

técnica. Mostrando-se aqui também, capazes de serem utilizados como métodos de

marcha no tempo para o estudo numérico do escoamento ao redor de um cilindro para

número de Reynolds 200. E mais uma vez, os métodos de Euler e “The Midpoint”,

mostraram-se eficientes quando comparados com os demais resultados encontrados na

literatura.

O campo de pressões com as linhas de corrente podem ser visualizadas na Figura

5.9 para os diferentes métodos de marcha no tempo. Estes agora indicam não apenas o

comportamento do fluido para o presente regime de escoamento, mas também

evidenciam que quanto maior o número de Reynolds mais intenso ficam os núcleos de

baixa pressão e que estes são transportados pelo escoamento formando assim uma

63

esteira de vórtices de von Kárman. Também como esperado, existe um padrão oscilante

das linhas de corrente causado pelo descolamento destes vórtices.

As Figuras Figura 5.10 e Figura 5.11 mostram a variação dos coeficientes de

arrasto (𝐶𝐷) e de sustentação (𝐶𝐿) em função do tempo para todos os métodos de marcha

no tempo, respectivamente. Assim como o observado para o número de Reynolds 100, o

coeficiente de sustentação (𝐶𝐿) para Reynolds 200 não encontrou maiores diferenças

entre um método e outro, possuindo uma concordância com todos os resultados

numéricos e experimentais apresentados na Tabela 13. Tendo o coeficiente de arrasto o

mesmo comportamento, possuindo também uma grande concordância com os resultados

numéricos e experimentais também apresentados na Tabela 13.

(a)

(c)

(b)

(d)

Figura 5.9 – Coeficiente de pressão e streamlines para o número de Reynolds 200. (a) Método de Euler; (b) Método

de Heun; (c) Método “The Midpoint”; e (d) Método de Ralston.

64

Figura 5.10 – Comparação do coeficiente de arrasto para os diversos métodos de marcha no tempo (Reynolds 200).

Figura 5.11 – Comparação do coeficiente de sustentação para os diversos métodos de marcha no tempo (Reynolds

200).

Na Figura 5.12, estão mostradas as oscilações adimensionalizadas do coeficiente

de sustentação no domínio da frequência. Onde podemos observar, que todos os

métodos de marcha no tempo capturam a mesma frequência de desprendimento de

vórtices. Estes também concordaram de forma satisfatória com os resultados numéricos

e experimentais presentes na Tabela 13.

65

(a)

(c)

(b)

(d)

Figura 5.12 – Número de Strouhal (𝑆𝑡) para o número de Reynolds 200. (a) Método de Euler; (b) Método de Heun;

(c) Método “The Midpoint”; e (d) Método de Ralston.

6 Conclusão

Todos os resultados numéricos apresentados aqui concordaram satisfatoriamente

tanto com os resultados numéricos quanto com os resultados experimentais encontrados

na literatura técnica para os regimes de escoamento de Reynolds 100 e 200. Entretanto,

para o regime de escoamento de Reynolds 40 os métodos de marcha no tempo de

Ralston e Heun não foram satisfatórios. Estes obteram resultados para o coeficiente de

arrasto muito acima do esperado, mesmo representando o escoamento com certa

fidelidade. Porém os resultados obtidos pelos métodos explícitos de Euler de primeira

ordem e o “The Midpoint” (Runge-Kutta de 2ª ordem), mostraram-se capazes de obter

resultados satisfatórios para todos os regimes de escoamentos aqui estudados (Reynolds

40, 100 e 200).

66

Ressaltando, que estes resultados só foram possíveis porque foi desenvolvido um

estudo da influência da localização das fronteiras do domínio computacional, do

refinamento da malha e dos métodos de marcha no tempo no campo de escoamento.

Mostrando com isso, que estes afetam de forma considerável o campo de escoamento.

Tendo como resultado direto, um domínio computacional, uma malha e métodos de

marcha no tempo ditos ótimos para os regimes de escoamentos aqui estudados.

Com isso, os resultados apresentados pelos métodos explícitos de Euler e de

Runge-Kutta (“The Midpoint”) neste trabalho, demonstram que estes foram validados e

são possíveis obter resultados satisfatórios na investigação da vibração induzida por

vórtices em estruturas cilíndricas, como por exemplo risers. Ratificando assim, tanto a

qualidade quanto a eficiência destes códigos aqui implementados em volumes finitos

através destes métodos de marcha no tempo.

7 Trabalhos Futuros

Como trabalhos futuros fica aqui o desejo de transpor o código para a linguagem

C++, revisar o código com o intuito de otimizá-lo e paralelizá-lo de forma a torná-lo

mais rápido. Pois os métodos aqui implementados são explícitos e por isso estes são

bastante restritivos quanto ao número de CFL máximo possível de ser utilizado.

8 Referências Bibliográficas

Anderson Jr., D.A., 1990, MODERN COMPPRESSIBLE FLOW With Historical

Perspective. 2 ed. Sigapore, MacGraw-Hill.

Behr, M., Hastreiter, D., Mittal, S., et al., 1995, “Incompressible flow past a circular

cylinder: dependence of the computed flow field on the location of the lateral

boundaries”, Computer Methods in Applied Mechanics and Engineering, v. 123, pp.

209-316.

Behr, M., Liou, J., Shih, R., et al., 1991, “Vorticity-stream function formulation of

unsteady incompressible flow past a cylinder: Sensitivity of the computed flow field to

the location of the out flow boundary", International Journal for Numerical Methods

in Fluids, v. 12, pp. 323-342.

67

Chung, T.J., 2002, Computational Fluid Dynamics, 1 ed. Cambridge, Cambridge

University Press.

Constanceau, M., Bouard, R., 1977, “Experimental determination of the main features

of viscous flow in the wake of a circular cylinder in uniform translation steady flow”,

Journal of Fluid Mechanics, v. 79, pp. 231-256.

Fletcher, C.A.J., 1991, Computational Techniques for Fluid Dynamics 1

Fundamental and General Techniques. 2 ed. Germany, Springer-Verlag.

Fox, Robert W., McDonald, Alan T., 2001, INTRODUÇÃO À MECÂNICA DOS

FLUIDOS. 5 ed. Brasil, LTC – Livros Técnicos e Científicos.

Herfjord, K., 1995, A study of two-dimensional separated flow by a combination of the

finite element method and Navier-Stokes equations, Ph.D. dissertation, The Norwegian

Institute of Technology, Trondheim, Norway.

Kundu, Pijush K., Cohen, Ira M., 2001, FLUID MECHANICS. 2 ed. USA, Academic

Press.

Mavriplis, D., Jameson, A., 1987, “Multigrid solution of the two-dimensional Euler

equations on unstructured triangular meshes”. AIAA 25th Aerospace Sciences Meeting,

87-0353, Reno, Nevada, USA, 12-15 January.

Munson, Bruce R., Young, Donald F., Okiishi, Theodore H., 2004, FUNDAMENTOS

DA MECÂNICA DOS FLUIDOS. 4 ed. Brasil, Editora Edgard Blücher.

Norberg, C., 2003, “Fluctuating lift on a circular cylinder: review and new

measurements”, Journal of Fluids and Structures, v. 17, pp. 57-96.

Miliante, T.R.C.C, 2006, Simulação numérica de vibrações induzida por vórtice em

dois cilindros em linha em baixo Reynolds, Trabalho de Formatura, DENO/UFRJ, Rio

de Janeiro, RJ, Brasil.

Patankar, S.V., 1980, Numerical Heat Transfer and Fluid Flow, 1 ed. USA,

Hemisphere.

Rengel, J.E., Sphaier, S.H., 1999, “A projection method for unsteady Navier-Stokes

equation with finite volume method and collocated grid”, Hybrid Methods in Heat and

Mass Transfer, v. 1, no 4.

68

Sumer, B. Mutlu, FredsØe, JØrgen, 2006, HYDRODYNAMICS AROUND

CYLINDRICAL STRUCTURES. Revised Edition. Sigapore, World Scientific.

Tritton, D. J., 1959, “Experiments on the flow past a circular cylinder at low Reynolds

number”, Journal of Fluid Mechanics, v. 6, pp. 547-567.

Wanderley, J.B.V, 2001, “An algorithm for slightly compressible flows”. Proceedings

of the 20th International Offshore Mechanics and Arctic Engineering Conference, OFT-

1191, Rio de Janeiro, RJ, Brazil, 3-8 June.

Wanderley, Juan B. V., Souza, Gisele H.B., Sphaier, Sergio H., et al., 2008a, “Vortex-

induced vibration of an elastically mounted circular cylinder using an upwind TVD

two-dimensional numerical scheme”, Ocean Engineering, v. 35, pp. 1533-1544.

Wanderley, Juan B. V., Sphaier, Sergio H., Levi, Carlos, 2008b, “Uma investigação

numérica da vibração induzida por vórtices sobre um cilindro rígido elasticamente

apoiado”, 28th Congresso da Sobena, Rio de Janeiro, RJ, Brasil, .

Wieselsberger, C., 1921, “New data on the law of hydro and aerodynamic resistance”,

Physikalische Zeitschrift (in German), v. 22, pp. 321-382.

Versteeg, H.K., Malalasekera, W., 1995, An introduction to Computational Fluid

Dynamics The Finite Volume Method. 1 ed. Malaysia, Longman Scientific &

Technical.