ESTUDO NUMÉRICO DO ESCOAMENTO AO REDOR …monografias.poli.ufrj.br/monografias/monopoli10029663.pdfcomparados com os resultados apresentados na norma IEC 61400-12-1. vii Graduation

ESTUDO NUMÉRICO DO ESCOAMENTO AO REDOR DE UMA TORRE

ANEMOMÉTRICA UTILIZANDO O MÉTODO DE VÓRTICES: COMPARAÇÃO

COM A IEC 61400-12-1

Rodrigo Rodrigues de Souza Pereira

Projeto de Graduação apresentado ao Curso de

Engenharia Mecânica da Escola Politécnica da

Universidade Federal do Rio de Janeiro como

parte dos requisitos necessários à obtenção do

título de Engenheiro Mecânico.

Professor: D. Sc., Roney Leon Thompson

Rio de Janeiro

Setembro de 2019



COM A IEC 61400-12-1


PROJETO DE GRADUAÇÃO APRESENTADO AO CURSO DE ENGENHARIA

MECÂNICA DA ESCOLA POLITÉCNICA DA UNIVERSIDADE FEDERAL DO

RIO DE JANEIRO COMO PARTE DOS REQUISITOS NECESSÁRIOS À

OBTENÇÃO DO TÍTULO DE ENGENHEIRO MECÂNICO.

Examinada por:

_______________________________________________

Prof., Roney Leon Thompson, D. Sc.

_______________________________________________

Vanessa Gonçalves Guedes , D. Sc.

_______________________________________________

Prof., Gustavo Cesar Rachid Bodstein, D. Sc.

RIO DE JANEIRO, RJ - BRASIL

SETEMBRO DE 2019

iii

Pereira, Rodrigo Rodrigues de Souza

Estudo Numérico do Escoamento ao Redor de uma Torre

Anemométrica Utilizando o Método de Vórtices:

Comparação com a IEC 61400-12-1/ Rodrigo Rodrigues de

Souza Pereira – Rio de Janeiro: UFRJ/Escola Politécnica,

2019.

XII, 65 p.: il.; 29,7 cm.

Orientador: Roney Leon Thompson

Projeto de Graduação – UFRJ/ Escola Politécnica/ Curso

de Engenharia Mecânica, 2019.

Referências Bibliográficas: p. 61-65.

1.Escoamento Atmosférico 2.Torre Anemométrica

3.Método de Vórtices Discretos 4.Método dos Painéis 5.IEC

61400-12-1. I. Thompson, Roney Leon. II. Universidade

Federal do Rio de Janeiro, Escola Politécnica, Programa de

Engenharia Mecânica. III. Estudo Numérico do Escoamento

ao Redor de uma Torre Anemométrica Utilizando o Método

de Vórtices: Comparação com a IEC 61400-12-1.

iv

AGRADECIMENTOS

Fui criado em uma família que trabalhava na indústria de óleo e gás – meu avô

paterno trabalhou na Shell e na Petrobras e, em seguida, ambos os meus pais

trabalharam para a Petrobras. Por esse motivo, não só em encontros familiares, mas

também em eventos com amigos de meus parentes, os assuntos normalmente

perpassavam petróleo e energia. Foi somente no início dos anos 2000, quando minha

mãe passou a se envolver em projetos de arquitetura na empresa com foco no padrão

LEED (Leadership in Energy and Environmental Design), que tive meu primeiro

contato com sustentabilidade. À medida que aprofundava meus conhecimentos em

energia e sustentabilidade, entendi que esses temas são interdependentes e que poderia

usar meus conhecimentos em engenharia para ajudar a encontrar soluções sustentáveis

para atendimento da demanda energética por meio do uso de energia renováveis.

Em julho de 2016, iniciei um estágio em pesquisa no CEPEL (Centro de

Pesquisas de Energia Elétrica da Eletrobras), onde trabalhei por 2 anos ajudando a

equipe do atual Departamento de Materiais, Eficiência Energética e Geração

Complementara (DME) a otimizar parques eólicos, resultando numa melhor

distribuição das turbinas que levaram a ganhos significativos para as empresas

controladoras e para a indústria como um todo.

Esses anos foram formativos no meu desenvolvimento profissional, adicionando a

pesquisa no meu rol de paixões: nesse período, ajudei na produção de 6 artigos, muitos

dos quais foram apresentados em conferências de energia pelo Brasil. Esse trabalho foi

a inspiração para usar meus anos estudando, pesquisando e trabalhando com energias

renováveis para analisar como otimizar a geração de energia eólica, tema sobre o qual

me debrucei para desenvolver o trabalho a seguir.

Portanto, agradeço a minha família, a meus pais – Fernando e Christina – e meus

irmãos – Guilherme, Leonardo, Michelle e Bianca – pela educação e apoio dado ao

longo deste tempo. Agradeço aos meus colegas de faculdade, em especial ao Daniel

Agnese Ramos, pelo esclarecimento de dúvidas sobre o mercado energia eólica.

Agradeço aos colegas do Cepel e à D.Sc. Vanessa Guedes Gonçalves pela orientação e

conhecimentos ao longo do estágio e do desenvolvimento deste trabalho. Por fim,

agradeço aos professores do Programa de Engenharia Mecânica da UFRJ pelos

ensinamentos ao longo do curso e ao Prof. D.Sc. Roney Leon Thompson pela

v

preocupação em passar esses ensinamentos através de uma didática excepcional em sala

de aula.

vi

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

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



COM A IEC 61400-12-1


Setembro/2019

Orientador: Roney Leon Thompson

Curso: Engenharia Mecânica

Este trabalho utilizou o método de vórtices discretos associado ao método dos painéis

com singularidade constante, para simular numericamente o escoamento incompressível

bidimensional ao redor de um cilindro triangular. Após a determinação de hipóteses

simplificadoras e de condições de não-escorregamento e impenetrabilidade no corpo, o

algoritmo lagrangiano de convecção e difusão dos vórtices com avanço de 2ª ordem de

Adam-Bashfotth foi comparado a resultados experimentais presentes na literatura.

Utilizou-se de um cilindro retangular para validação pela existência de um maior

número de trabalhos na literatura acadêmica e pela presença de cantos vivos similares

ao cilindro triangular. A motivação para o estudo do escoamento ao redor do cilindro

triangular é o fato de que esta geometria se aproxima de torres treliçadas utilizadas em

campanhas de medição de vento de parques eólicos. Os resultados obtidos são

comparados com os resultados apresentados na norma IEC 61400-12-1.

vii

Graduation Project’s abstract presented to Polytechnic School/UFRJ as a mandatory

requirement for obtaining the title of Mechanical Engineer.

NUMERICAL STUDY OF THE ATMOSFERIC FLOW AROUND ANEMOMETRIC

TOWERS BY USING THE VORTEX METHOD: A COMPARISON WITH IEC

61400-12-1


September/2019

Advisor: Roney Leon Thompson

Course: Mechanical Engineering

This study used the discrete vortex method associated with the panel method to with

constant singularity to simulate an incompressible, bidimensional flow past a triangular

cylinder. Firstly, simplifying hypothesis, no slip and impenetrability conditions were

considered to form the lagrangean algorithm with vortex convection and diffusion and a

two-step Adam-Bashfotth method. The results were compared to experimental academic

results. The validation geometry chosen was the rectangular cylinder due to its sharp

corners and its frequency in the academic literature. This study was motivated by the

resemblance of an anemometric tower used to measure wind speed campaigns and a

triangular cylinder. The final results were then compared to the IEC 61400-12-1.

viii

SUMÁRIO

1 Introdução .................................................................................................... 13

1.1 MOTIVAÇÃO ................................................................................................................ 13

2 Revisão Bibliográfica .................................................................................. 16

3 Metodologia ................................................................................................. 18

3.1 MODELAGEM MATEMÁTICA ................................................................................... 19

3.1.1 Domínio ....................................................................................................................................... 19

3.1.2 Hipóteses simplificadoras ........................................................................................................... 20

3.2 PROBLEMA DE VALOR DE CONTORNO ................................................................ 21

3.2.1 Adimensionalização do problema de valor de Contorno ........................................................... 22

3.2.2 Equação do Transporte da Vorticidade ...................................................................................... 23

3.2.3 Método de Vórtices Discretos (MVD) ......................................................................................... 24

3.2.4 Método dos Painéis .................................................................................................................... 27

3.2.4.1 Escoamento Potencial ............................................................................................................ 29

3.2.4.2 Velocidades Induzidas............................................................................................................. 31

3.2.4.3 Somatório do Escoamento Incidente com as Velocidades Induzidas..................................... 34

3.2.4.4 Condição de impenetrabilidade e de não-escorregamento ................................................... 36

3.2.5 Implementação Numérica .......................................................................................................... 36

3.2.5.1 Discretização do Corpo ........................................................................................................... 37

3.2.5.2 Cálculo da Matriz de Influência .............................................................................................. 38

3.2.5.3 Geração de Vorticidade .......................................................................................................... 38

3.2.5.4 Convecção dos Vórtices .......................................................................................................... 38

3.2.5.5 Reflexão de Vórtices ............................................................................................................... 39

3.2.5.6 Difusão dos Vórtices ............................................................................................................... 39

3.2.5.7 Cálculo do Coeficiente de Pressão sobre o Corpo .................................................................. 40

3.2.5.8 Cálculo de Forças sobre o Corpo ............................................................................................ 40

3.2.5.9 Avanço no Tempo ................................................................................................................... 40

3.2.5.10 Conferência de Critérios de Parada ........................................................................................ 41

4 Resultados .................................................................................................... 42

ix

4.1 VALIDAÇÃO COM O CILINDRO DE SEÇÃO RETANGULAR .............................. 42

4.2 RESULTADOS PARA O CILINDRO DE SEÇÃO TRIANGULAR ............................ 47

4.2.1 CASO I – Seção a 0° ..................................................................................................................... 48

4.2.2 CASO II – Seção a 30° .................................................................................................................. 51

4.2.3 CASO III – Seção a 60° ................................................................................................................. 54

5 Comparação e Discussão ............................................................................ 57

6 Conclusão ..................................................................................................... 59

7 Referências ................................................................................................... 61

x

LISTA DE FIGURAS

Figura 1 – Geometria e região fluida .............................................................................. 20 Figura 2 – Vórtice circular ............................................................................................. 25

Figura 3 – Vórtices de Lamb .......................................................................................... 26 Figura 4 – Pontos para a formação de painéis ................................................................ 27 Figura 5 – Definição dos pontos de controle e dos pontos de geração de vórtices ........ 28 Figura 6 – Rotação do corpo em torno da origem .......................................................... 28

Figura 7 – Grandezas para o calculo de 𝐼𝑓𝑛 e 𝐼𝑓𝜏 ......................................................... 33 Figura 8 – Fluxograma do algoritmo implementado ...................................................... 37 Figura 9 – Discretização do cilindro retangular e seus pontos de controle .................... 42 Figura 10 – Valores de Cp ao redor do corpo para um cilindro circular: (a) simulação

numérica realizada (b) experimental (LEE, 1975) ......................................................... 43

Figura 11 – Coeficientes de arrasto e sustentação ao longo do tempo de simulação para

o cilindro retangular........................................................................................................ 44

Figura 12 – Linhas de corrente no entorno de um cilindro retangular à Re = 105 obtida

em (TAMURA, OHTA e KUWAHARA, 1990) ............................................................ 45 Figura 13 – Campo de velocidade para simulação do cilindro retangular ..................... 45 Figura 14 – Campo de velocidade com as linhas de corrente para simulação do cilindro

retangular ........................................................................................................................ 46 Figura 15 – Dispersão de vórtices na esteira do cilindro retangular .............................. 46

Figura 16 – Ângulos de incidência dos cilindro triangulares ......................................... 47 Figura 17 – Resultados experimentais do escoamento ao redor de cilindros triangulares

com variação de ângulos obtido em (IUNGO e BURESTI, 2009) ................................ 48

Figura 18 – Valores de Cp ao redor do corpo para a simulação numérica para o cilindro

triangular a 00º ................................................................................................................ 48


o cilindro triangular a 00º ............................................................................................... 49

Figura 20 – Campo de velocidade para simulação do cilindro triangular a 00º ............. 49 Figura 21 – Campo de velocidade com as linhas de corrente para simulação do cilindro

triangular a 00º ................................................................................................................ 50 Figura 22 - Dispersão de vórtices na esteira do cilindro triangular a 00º ....................... 50


o cilindro triangular a 30º ............................................................................................... 51 Figura 24 – Valores de Cp ao redor do corpo para a simulação numérica para o cilindro

triangular a 30º ................................................................................................................ 51 Figura 25 – Campo de velocidade para simulação do cilindro triangular a 30º ............. 52

Figura 26– Campo de velocidade com as linhas de corrente para simulação do cilindro

triangular a 30º ................................................................................................................ 52 Figura 27 – Dispersão de vórtices na esteira do cilindro triangular a 30º ...................... 53 Figura 28 – Coeficientes de arrasto e sustentação ao longo do tempo de simulação para

o cilindro triangular a 60º ............................................................................................... 54 Figura 29 – Valores de Cp ao redor do corpo para a simulação numérica para o cilindro

triangular a 60º ................................................................................................................ 54

Figura 30 – Campo de velocidade para simulação do cilindro triangular a 60º ............. 55 Figura 31 – Campo de velocidade com as linhas de corrente para simulação do cilindro

triangular a 60º ................................................................................................................ 55 Figura 32 – Dispersão de vórtices na esteira do cilindro triangular a 60º ...................... 56

xi

Figura 33 – Gráfico de comparação entre os Cd’s obtidos e os resultados apresentados

em (IUNGO e BURESTI, 2009) .................................................................................... 57 Figura 34 – Gráfico de comparação entre os Cl’s obtidos e os resultados apresentados

em (IUNGO e BURESTI, 2009) .................................................................................... 57

Figura 35 – Comparação entre simulações: (a) Média temporal da série de simulações a

00º (b) Imagem apresentada na Norma IEC 61400-12-1 .............................................. 59

xii

LISTA DE TABELAS

Tabela 1 – Comparação de resultados para o cilindro retangular................................... 44 Tabela 2 – Compilado dos resultados do cilindro retangular para diferentes ângulos ... 58

13

1 INTRODUÇÃO

1.1 MOTIVAÇÃO

A geração de energia em larga escala se desenvolveu, primordialmente, em torno

do uso de combustíveis não renováveis, sendo esse o meio mais usual de suprir nossa

demanda energética. Isso remonta da Revolução Industrial, quando a fonte energética

que impulsionou a automação dos meios de produção era proveniente de um

combustível fóssil: o carvão mineral, que alimentava as primeiras máquinas a vapor

(FOUQUET e PEARSON, 2012). Nos séculos XIX, com o querosene substituindo o

óleo de baleia como a principal fonte de iluminação pública (YERGIN, 2010), e XX,

com a gradual massificação dos motores a combustão e a ascensão da gasolina, elevou-

se o status dos derivados do petróleo, os quais se tornaram as mais importante fontes de

energia (GRUBLER, 2012). Essa dependência ao petróleo impactou todos os setores da

sociedade, sendo responsável por novas dinâmicas geopolíticas, com crises econômicas

nacionais e internacionais (SILVÉRIO e SKLO, 2012), e especulação nos mercados de

commodities (TOKIC, 2010). O uso de combustíveis fósseis também foi responsável

por problemas ambientais ao longo de sua cadeia produtiva, bem como o aumento de

incidência de doenças pulmonares e oftalmológicas (DOCKERY, SCHWARTZ e

SPENGLER, 1992) decorrentes da emissão de poluentes (FINLAYSON-PITTS e

PITTS JR., 1997). Além dessas discussões, outra questão relevante é a previsão do

momento da escassez do petróleo, a chamada previsão do pico de produção mundial do

petróleo – o ponto a partir do qual a produção em novos campos não vai mais

compensar a queda natural de produção nos campos ao redor do mundo ( (HÖÖK,

HIRSCH e ALEKLETT, 2009) e (OWEN, INDERWILDI e KING, 2010)).

Esses debates chamaram a atenção da população, governos e empresas, que viram

o desenvolvimento de novas tecnologias que permitissem a conservação do meio-

ambiente e da nossa saúde como uma necessidade e uma oportunidade de diversificar

seu portfólio energético, sem perder de vista questões fundamentais como viabilidade

econômica e escalabilidade. Apesar do aumento da quantidade de fontes de energia, a

matriz energética mundial ainda tem grande dependência da queima de combustíveis

fósseis. As maiores dificuldades para sua substituição são sua densidade energética, a

14

facilidade logística de transporte e estoque dos derivados de petróleo em forma líquida

ou gasosa, o baixo preço em relação às outras formas não convencionais de energia, e a

infraestrutura desenvolvida ao longo de mais de um século com combustíveis fósseis

em mente.

Porém, graças à demanda pública, à perspectiva de escassez de energias não

renováveis, às preocupações com meio-ambiente e saúde, e ao barateamento de fontes

não convencionais de energia, as energias renováveis vêm se mostrando mais presentes

tanto no meio urbano quanto rural, sendo cada vez mais visadas, como indicado pelos

futuros investimentos mundiais em energias renováveis (IEA, 2017).

Uma forte evidência dessa mudança de paradigma pode ser observada pelo

aumento do interesse por construções sustentáveis, como o de zero-energybuildings

(CELLURA, GUARINO, et al., 2015), e fachadas solares (CAO, HASAN e SIRÉN,

2013), para sistemas off grid. Temos também um aumento da participação de energias

renováveis na matriz energética do mundo e em especial da energia eólica no Brasil,

que por tradição teve um grande foco em geração de energia hidroelétrica e que começa

a explorar seu potencial eólico e solar. O que vemos, então, é o aumento dos

investimentos em energia eólica e, com isso, nos estudos em relação aos melhores locais

para sua produção e maneiras de como otimizar a geração de energia por essa fonte.

Para a definição adequada e desenvolvimento de parques eólicos, medições

anemométricas – após a prospecção de uma área de interesse – e cálculos feitos a partir

dessas medições embasam a escolha adequada do tipo de aerogerador e o melhor layout

do site para produção de energia eólica – com foco na otimização do posicionamento

dessas turbinas para a maximização do retorno sobre o investimento, resultado do valor

do projeto subtraído do seu custo de investimento e, então, dividido pelo próprio custo

de investimento.

Escoamentos ao redor de corpos rombudos são caracterizados pela separação da

camada limite e pela formação de uma generosa esteira atrás do corpo. Esses fenômenos

tornam a previsão do escoamento muito difíceis, e muitas vezes métodos empíricos são

utilizados para obtenção de forças aerodinâmicas. A maioria dos estudos numéricos e

experimentais do escoamento em corpos rombudos se dão ao redor de um cilindro

circular (WISSINK e RODI, 2008). Uma variedade de métodos com e sem malhas

também são usados ( (CHORIN, 1973) e (LEONARD, 1980)). No geral, resultados para

escoamentos incompressíveis ao redor de uma seção transversal de um cilindro

15

triangular com diferentes ângulos de ataque, para um número de Reynolds maior que 1

x 104, são raros na literatura (IUNGO e BURESTI, 2009).

Este trabalho tem como um dos objetivos explorar o potencial de simulações

lagrangeanas bidimensionais incompressíveis ao redor de corpos rombudos, as quais

dispensam o uso de geração de complexas malhas e conseguem simular com sucesso

esse tipo de escoamento turbulento, concentrando a região de cálculo na camada limite e

na esteira viscosa.

O estudo a seguir se baseia no estudo de Guedes (2003), adaptando o modelo

utilizado para descrever os escoamentos ao redor de corpos rombudos, com seções

retangulares, para um modelo que descreva como esses escoamentos se comportam

quando a seção é um triângulo equilátero em diferentes ângulos em relação ao fluxo de

vento. O interesse na geometria escolhida é devido ao fato de o uso da torre treliçada

para medições anemométricas está prevista na norma atual, com padrões de montagem

exigidos estabelecidos pela IEC 61400-12-1. Como a metodologia aplicada não possui

alto tempo computacional devido à suas vantagens citadas, o objetivo é poder realizar

simulações mais detalhadas do escoamento (ao invés de utilizar fórmula empírica

sugerida por norma) em função da direção preferencial do escoamento no local de

interesse, a fim de melhor quantificar as incertezas nas medições de velocidade de vento

em função do posicionamento da torre.

16

2 REVISÃO BIBLIOGRÁFICA

O trabalho apresentado nos próximos capítulos foi constituído a partir de uma

simulação numérica de adaptação própria. Para aferir se tal adaptação condiz com a

realidade física, se faz necessário um melhor entendimento dos casos já estudados e

disponíveis na literatura acadêmica.

Ao iniciar o estudo, deve-se observar o fenômeno físico em um ambiente com

condições controladas para validação de modelos matemáticos que descrevam o

problema. Por isso, de forma a tornar este trabalho inteligível, aborda-se primariamente

trabalhos experimentais.

Aqui estão listados alguns dos principais trabalhos para os quatro casos que serão

apresentados no Capítulo 4 – Resultados. A começar pelo caso utilizado para validação

do modelo numérico – cilindro retangular de razão de aspecto igual à 1 – e pelo caso

que foco é deste trabalho – cilindro triangular.

Rosenhead (1963) trata do desenvolvimento e da estabilidade da camada limite

laminar para escoamentos incompressíveis. Em Rosenhead (1931) tentou-se, pela

primeira vez, a simulação do campo de vorticidade por meio da implementação do

método de vórtices discretos. O problema desta implementação era a divergência que

ocorria após alguns passos de tempo.

Vickery (1966) apresentou um estudo da variação de arrasto e sustentação para

um cilindro retangular de razão de aspecto igual à 1. Neste mesmo trabalho, foram

aplicadas pequenas variações angulares, ideia similar à utilizada neste trabalho para a

seção triangular.

Bearman e Trueman (1972) estudou a possibilidade da redução do coeficiente de

arrasto em cilindro retangulares por meio de um modelo físico bidimensional.

Chorin (1973) traz em seu trabalho um método numérico para a solução das

equações de Navier-Stokes no espaço bidimensional. Este método serve de base para

este trabalho, com a ideia principal da geração e dispersão de vórtices. Neste trabalho

foi utilizado o método do avanço randômico.

Lee (1975) mediu o campo de pressão atuante em um cilindro retangular imerso

em um escoamento uniforme e turbulento e como o aumento da turbulência impacta o

arrasto no corpo.

17

Laneville e Yong (1983) mostraram a visualização de um escoamento médio ao

redor de um cilindro retangular bidimensional. Nele é possível ver uma bolha de

separação na parede do cilindro, bem como a distribuição de vorticidade.

Blevins (1984) dispõe das aplicações do estudo de fluidos para engenharia e

técnicas computacionais. Neste trabalho são demonstradas abordagens para solução das

equações de Navier-Stokes e da camada limite. Ainda neste, Blevins traz medidas de

turbulência e escoamentos potenciais.

O livro de Lewis (1991) aborda algoritmos com as programações de geração de

vorticidade na superfície do corpo e o desprendimento da esteira. Em Lewis (1981), o

Método dos Painéis é utilizado para o Método de Vórtices. Isso facilitou o estudo de

geometrias com maior complexidade.

O trabalho de Tamura, Ohta e Kuwahara (1990) foi usado como referência na

validação do problema. Nele é gerada uma malha radial para visualização das linhas de

corrente no entorno de um cilindro retangular de razão de aspecto 1 imerso em um

fluido com Re = 105.

Norberg (1993) fez diversos experimentos para observar o escoamento que passa

por um cilindro retangular com diferentes razões de aspecto e diferentes ângulos de

incidência.

Chen e Liu (1999) realizou um trabalho similar ao de (NORBERG, 1993), com a

maior diferença entre eles sendo a faixa de número de Re estudada e a exclusividade da

razão de aspecto de 1.

A tese de doutorado de Alonso (2005) trata do estudo de diferentes perfis de

cilindros triangulares com variações angulares. Nele são levantados os coeficiente de

arrasto, de pressão e de sustentação para números de Reynolds entre 104 e 2,6 x 105.

Iungo e Buresti (2009) fizeram um trabalho no qual o foco é aumentar a

quantidade de dados para esta geometria, escassa na literatura. Este estudo experimental

feito com anemometria de fio quente é a comparação utilizada para validar os resultados

numéricos obtidos pelo algoritmo adaptado de Guedes (2003).

18

3 METODOLOGIA

A modelagem matemática utilizada, bem como a implementação númerica-

computacional deste trabalho, foram adaptadas a partir do trabalho original de Guedes

(2003), para escoamentos ao redor de cilindros circulares e cilindros retangulares.

O estudo de escomentos se baseia em duas principais metodologias, experimental

e numérica-computacional. A primeira delas consiste em um ambiente – geralmente em

escala – com condições controladas, onde o estudo do escoamento e a interação entre o

corpo e o fluido é observado para coleta de resultados. Após essa coleta, as conclusões

são realizadas e finalmente as leis são verificadas.

A outra metodologia consiste na implementação numérica de modelos para

representar experimentos. Devido à complexidade de acesso a laboratórios com

condições controladas, à complexidade de corpos estudados e aos avanços

computacionais, essa metodologia se estabeleceu como a mais comum nos dias atuais e

será aplicada desenvolvida neste estudo.

Pode-se apontar quatro critérios para caracterização do método aplicado em

estudos por simulações numéricas, são eles:

Causa x Consequência

Abordagem Euleriana x Lagrangeana

Invariância Galileana x Euclidiana

Dependência x Independência do Usuário

Primeiramente, a diferença entre “Causa” e “Consequência”, como o próprio

nome infere, diz respeito a vorticidade estar ligada a um conjunto de forças ou ser

consequência de um conjunto de forças.

A abordagem Euleriana atenta a região onde ocorrem as vorticidades. Geralmente

é feita com implementação de um sistema de geração de malhas para determinar os nós

de interesse. Já abordagem Lagrangiana acompanha a partícula gerada e percebe seu

desenvolvimento ao longo do tempo.

A invariância Galileana trata do princípio de relatividade, pelo qual há

independência da movimentação de um observador para descrição das leis da natureza,

enquanto a invariância Euclidiana se baseia na dependência do espaço e leva em conta

apenas ângulos e distâncias em sua descrição.

19

Por fim, a dependência do usuário impacta na forma com que o método numérico

será implementado, sendo ajustada a sensibilidade de acordo com os parâmetros

escolhidos de forma empírica pelo usuário.

Esse trabalho tem uma abordagem lagrangeana, com parâmetros inputados pelo

usuário. Nele há dependência do espaço e tem a vorticidade como consequência de um

conjunto de forças provenientes dos pontos gerados a cada iteração.

3.1 MODELAGEM MATEMÁTICA

Primeiramente, apresenta-se a modelagem matemática escolhida para abordar o

problema. Esta modelagem começa pela região de domínio em que se encontra a

geometria e pelas hipóteses simplificadoras e adimensionalizações utilizadas para o

escoamento em questão. Com essas hipóteses estabelecidas, é tratado o Problema de

Valor de Contorno, de onde se obtém a equação de transporte de vorticidade.

Posteriormente, apresenta-se a modelagem numérica-computacional utilizada,

com suas simplificações e adaptações para resolução do problema específico deste

trabalho, cilindros triangulares. Neste trabalho especificamente, escolheu-se por utilizar

o Método dos painéis para modelagem pelo vasto uso em estudos aerodinâmicos de

aerofólios e grande validação na literatura acadêmica (ROY e BANDYOPADHYAY,

2006).

3.1.1 Domínio

Considerou-se um escoamento uniforme para o problema estudado. O sistema de

coordenadas polares foi escolhido como a melhor forma de se analisar as fronteiras da

região fluida: Sc (superfície do corpo) e S∞ (superfície no infinito). Suas fronteiras são

definidas como:

𝑆𝑐 ∶ 𝑟′ = 𝑟′𝑜 (3.1)

S∞: 𝑟′ → ∞ (3.2)

𝑆 = 𝑆𝑐 ∪ S∞ (3.3)

20

Figura 1 – Geometria e região fluida

3.1.2 Hipóteses simplificadoras

Foram consideradas algumas hipóteses simplificadoras para este estudo. São elas:

1 – Hipótese do Contínuo

Essa hipótese estabelece que as propriedades do fluido são continuas, por menor

que seja a resolução estudada. Ou seja, a composição molecular, e seus átomos

discretos, não impacta o movimento do fluido, permitindo este ser estudado através do

cálculo de funções contínuas.

As propriedades atômicas e moleculares são consideradas em outras variáveis

como massa específica e viscosidade.

2 – Fluido newtoniano com viscosidade constante

O fluído tem tensão de cisalhamento diretamente proporcional à taxa de

deformação e viscosidade aparente constante (FOX, PRITCHARD e MCDONALD,

2013).

3 – Fluido incompressível

21

Para velocidades de escoamento muito altas, a massa específica do fluido é

afetada. Isso caracteriza fluidos compressíveis, onde o numero de Mach (Ma = U/c,

onde c é a velocidade do som) é maior que 0,3.

O estudo em questão lida com escoamentos atmosféricos e, portanto, com

velocidades muito abaixo da velocidade do som. Por isso, a hipótese de não alteração da

massa específica pela incompressibilidade do fluido foi utilizada.

4 – Velocidade incidente U constante sobre o corpo

Considerou-se a incidência de um escoamente constante U como modelagem do

escoamento atmosférico na torre anemométrica modelada.

5 – Escoamento bidimensional

Por simplicidade de implementação e pelo Método de Painéis escolhido para

modelagem do corpo, foi considerada apenas uma seção do cilindro estudado e

consequente bidimensionalização do problema.

6 – Região fluida infinita

A região fluida é expansível de acordo com o andamento dos vórtices calculados

para geração de resultados, onde o fator limitante para o tamanho do domínio é a

capacidade de processamento e o tempo da máquina em que a simulação é rodada.

7 – Não há interferência de fatores externos ao corpo

O fenômeno físico foi considerado como sendo afetado única e exclusivamente

pelos métodos elencados neste trabalho, sem que houvesse considerações de fatores do

entorno do corpo que pudessem impactar de outras formas na vorticidade.

3.2 PROBLEMA DE VALOR DE CONTORNO

As equações de Navier-Stokes regem o escoamento ao redor de um corpo.

Consistem dos princípios de conservação de massa e da conservação de quantidade de

movimento.

∇′. 𝒖′ = 0 , na região fluida (3.4)

𝐷𝒖′

𝐷𝑡′ = 𝜕𝒖′

𝜕𝑡′ + 𝒖′. ∇′𝒖′ = −1

𝜌∇𝑝′ + 𝜈∇′2

𝒖 (3.5)

O escoamento na superfície do corpo atende a premissa de não escorregamento e

de impenetrabilidade.

22

𝒖′ = 𝒒′ , em 𝑟′ = 𝑟′0 (3.6)

onde q’ é a velocidade da superfície do corpo. Ao se separar a equação 3.6 nas

componentes normal e tangencial, obtém-se:

𝑢′𝒏 = 𝑞′

𝒏 , em 𝑟′ = 𝑟′0 → Condição de Impenetrabilidade (3.7)

𝑢′𝒕 = 𝑞′

𝒕 , em 𝑟′ = 𝑟′0 → Condição de Escorregamento nulo (3.8)

Estando o corpo em repouso com o fluido se movimento a uma velocidade U,

tem-se 𝒒′ = 0 e, portanto:

𝑢𝒏 = 𝒖′. 𝒏 = 0 , em 𝑟′ = 𝑟′0 (3.9)

𝑢𝒕 = 𝒖′. 𝒕 = 0 , em 𝑟′ = 𝑟′0 (3.10)

Com a condição de contorno em 𝑆∞:

|𝒖′| → 𝑈 , em 𝑟′ → ∞ (3.11)

3.2.1 Adimensionalização do problema de valor de Contorno

O escoamento bidimensional é influenciado por grandezas como o comprimento

característico, a massa específica, a viscosidade do fluido e a velocidade do escoamento

incidente. Essas grandezas podem ser adimensionalizadas por:

𝑥 ≡𝑥′

𝑙, comprimento adimensionalizado na direção x (3.12a)

𝑦 ≡ 𝑦′

𝑙, comprimento adimensionalizado na direção y (3.132b)

𝑟 ≡ 𝑟′

𝑙, comprimento adimensionalizado na direção radial (3.12c)

𝑟0(𝜃) ≡ 𝑟0′(𝜃)

𝑙, forma do corpo adimensionalizada (3.12d)

𝑡 ≡ 𝑡′𝑈

𝑙, instante de tempo adimensionalizado (3.12e)

𝑢 ≡𝑢′

𝑈 , componente da velocidade adimensionalizada na direção x (3.12f)

𝑣 ≡ 𝑣′

𝑈, componente da velocidade adimensionalizada na direção (3.12g)

23

𝑝 ≡ 𝑝′

1

2𝜌𝑈2

, instante de tempo adimensionalizado (3.12h)

E tem-se, por fim o problema na região fluida na forma:

∇. 𝒖 = 0 (3.13)

𝐷𝒖

𝐷𝑡=

𝜕𝒖

𝜕𝑡+ 𝒖. ∇𝒖 = −∇𝑝 +

1

𝑅𝑒∇2𝒖 (3.14)

Com condições de impenetrabilidade, de escorregamento nulo e uma imposição

no infinito descritas anteriormente, respectivemente, se apresentam como:

𝑢𝒏 = 𝒖. 𝒏 = 0 (3.15)

𝑢𝒕 = 𝒖. 𝒕 = 0 (3.16)

|𝒖| → 1 , em 𝑟 → ∞ (3.17)

Calcula-se então os coeficientes de arrasto (Cd), de sustentação (Cl) e de Strouhal

(St).

𝐶𝑑 ≡ 𝐷

1

2𝜌𝑈2𝐴

(3.18)

𝐶𝑙 ≡ 𝐿

1

2𝜌𝑈2𝐴

(3.19)

𝑆𝑡 ≡𝑓𝐷

𝑈 (3.20)

Para validação da precisão do método numérico aplicado, os valores de Cd, Cl e St

são calculados ao fim de cada rodada.

3.2.2 Equação do Transporte da Vorticidade

O vetor vorticidade ω’ é definido como:

𝝎′ ≡ ∇′ × 𝐮′ (3.21)

Ao aplicar essa definição nas equações 3.4 e 3.5, o termo da pressão é cortado e a

equação resultante é mostrada a seguir:

24

𝐷𝝎′

𝐷𝑡′ = 𝜕𝝎′

𝜕𝑡′ + 𝒖′. ∇′𝝎′ = 𝝎′. ∇′𝒖′ + 𝜇∇′2𝝎′

(3.22)

Adimensionaliza-se essa equação com as variáveis das equações 3.12a, 3.12f e

𝝎 =𝑙𝝎′

𝑈, e obtem-se a equação de transporte da vorticidade adimensionalizada:

𝐷𝝎

𝐷𝑡=

𝜕𝝎

𝜕𝑡+ 𝒖. ∇𝝎 = 𝝎. ∇𝒖 + 𝜇∇2𝝎 (3.23)

Ao analisar a equação, entende-se o sentido físico das parcelas, tais quais:

𝒖. ∇𝝎 é a variação da vorticidade causada pela convecção.

1

𝑅𝑒∇2𝝎 é a taxa de difusão de vorticidade por ação da viscosidade

𝝎. ∇𝒖 é a taxa de deformação das linha de vorticidade.

Por se tratar de um escoamento bidimensional, tem-se que o valor de 𝝎. ∇𝒖 é

nulo. Portanto, a vorticidade é representada como uma grandeza escalar por estar

sempre perpendicular ao plano de escoamento, conforme a equação 3.24 abaixo:

𝐷𝜔

𝐷𝑡=

𝜕𝜔

𝜕𝑡+ 𝒖. ∇𝜔 = 𝜔. ∇𝒖 +

1

𝑅𝑒∇2𝜔 (3.24)

3.2.3 Método de Vórtices Discretos (MVD)

Para simulação numérica-computacional de vórtices bidimensionais, a serem

comparados com o método utilizado na IEC 61400-12-1, foi utilizado o Método de

Vórtices Discretos (MVD) pela maior simplicidade de implementação quando

comparado a outras metodologias, tais quais o Método de Elementos Finitos ou o

Método de Volumes Finitos, que necessitam de geração de malhas para o cálculo das

grandezas de forma euleriana.

Pode-se então separar a vorticidade em duas regiões. A primeira, no interior do

círculo amarelo de raio 𝜎0 da Figura 2 com valor igual a ω e a área externa ao círculo

com valor zero

25

Figura 2 – Vórtice circular

Pelo Teorema de Stokes, tem-se que:

∫ (∇ × 𝒒). 𝒏𝑑𝑆𝑆

= ∮ 𝒒𝑑𝑠𝐶

(3.25)

Ao resolvermos o Teorema de Stokes, obtém-se duas equações que variam de

acordo com o raio 𝑟∗ < 𝜎0 < 𝑟, com q e q* constantes no entorno da circunferência

para um mesmo r.

𝑞∗ =1

2𝝎𝒓∗, 𝑟∗ ≤ 𝜎0 (3.26)

𝑞 =𝝎𝜎0

2

2𝒓, 𝑟 > 𝜎0 (3.27)

Neste trabalho opta-se pelo modelo do vórtice de Lamb. Seu modelo elimina

singularidades e tornar as funções velocidade e vorticidade contínuas conforme

representado na Figura 3.

q*

r*

26

Figura 3 – Vórtices de Lamb

O potencial complexo adimensional Ѱ pode ser definido como a soma da parte

real composta pela função corrente ψ com a parte imaginária composta pelo potencial

de velocidade ϕ.

A velocidade complexa externa ao vórtice e a velocidade induzida por ele em um

ponto 𝑧 = 𝑟𝑒𝑖𝜃 pode ser escrita como:

𝑄 = 𝑞𝑒−𝑖(𝜃+𝜋

2) (3.28)

Define-se, também, as grandezas adimensionais: a intensidade do vórtice como

𝜅 =𝝎𝜎0

2

2 e a circulação como 𝛤 = 2𝜋𝜅.

Com essas definições, resolve-se a equação:

𝑄 =𝑑Ѱ

𝑑𝑧 (3.29)

E obtém-se a velocidade complexa induzida por um sistema de 𝑁𝑣 vórtices de

intensidade 𝛤𝜅, onde k=1, 2, 3, ..., 𝑁𝑣:

𝑄 = 𝑢 − 𝑖𝑣 = −𝑖

2𝜋∑

𝛤𝜅

𝑧−𝑧𝜅

𝑁𝑣𝜅=1 (3.30)

Na metodologia implementada neste trabalho, os vórtices são gerados de acordo

com a geometria e a discretização do corpo rombudo estudado por meio do Método dos

Painéis. Os vórtices gerados na superfície do corpo se deslocam por convecção e

difusão – caracterizando o método como puramente lagrangiano – e modelam o campo

de velocidade.

27

3.2.4 Método dos Painéis

Esta seção descreve a forma de desratização escolhida para o corpo imerso no

fluido. O Método dos Painéis empregado neste trabalho divide a fronteira onde ocorre a

interação fluido-corpo em painéis geradores de vórtices. Esta divisão é feita de forma

não uniforme a partir do número par de vórtices Nvgfc por iteração para cada face da

seção estudada. Com as arestas da seção definidas, geram-se pontos com distâncias

inversamente proporcionais a 2n a partir do ponto anterior (que serão os limites do

painéis) – a começar pelo centro da face e terminar em uma de suas arestas –, onde n é

uma lista dos números naturais não nulos até (𝑁𝑣𝑔𝑓𝑐 2) − 1⁄ .

Para exemplificação, na Figura 4 abaixo foi considerada uma geometria com face

de tamanho 2𝑥. Observa-se a formação dos pontos que definem a dimensão do paínel a

partir do centro da face com a proporcionalidade mencionada.

Figura 4 – Pontos para a formação de painéis

Após a formação de todos os pontos e painéis na face, é atribuído um ponto de

controle ao centro de cada painel gerado. Desses pontos de controle é dada uma

pequena distância ε, no sentido do vetor normal ao painel, onde serão gerados os

vórtices discretos a cada passo de tempo, conforme mencionado em Guedes (2003) e

demonstrado em um polígono arbitrário na Figura 5 abaixo.

28

Figura 5 – Definição dos pontos de controle e dos pontos de geração de vórtices

Como o foco deste trabalho é o estudo do escoamento ao redor de um corpo

rombudo, pode-se rotacionar a incidência do escoamento ou o corpo em si. Pelo fato das

duas maneiras apresentarem a mesma diferença angular relativa, escolheu-se por manter

a direção do escoamento de intensidade U constante, e rotacionar o corpo imerso no

fluido por um ângulo α em torno da origem, conforme Figura 6.

Figura 6 – Rotação do corpo em torno da origem

29

Dentre os diversos métodos para o cálculo de escoamento potencial, o escolhido

para este trabalho foi o Método dos Painéis. Algumas das vantagens deste método,

proveniente do Método de Elementos de Contorno (BREBBIA, TELLES e WROBEL,

1984), são a facilidade na aplicação para diversas geometrias e a comprovação da

representatividade dos resultados indicados na literatura acadêmica. Esta

representatividade, vinculada ao Método de Vórtices Discretos (MVD) descrita na seção

anterior, pode ser percebida pelos trabalhos de Lewis (1991), Guedes (2003), Alcântara

Pereira (1999) e Silva (2005).

O corpo discretizado por painéis de diferentes dimensões e imerso no fluido é

representado por uma série de singularidades com distribuição constante de fontes. Para

desenvolver a metodologia empregada, esta capítulo foi dividido em: cálculo do

escoamento potencial; cálculo das velocidades induzidas; somatório do escoamento

incidente com as velocidades induzidas; imposição das condições de não-

escorregamento e impenetrabilidade.

3.2.4.1 Escoamento Potencial

Conforme descrito em Batchelor (1967), o escoamento rotacional possui um

campo de velocidade associado que pode ser descrito como o gradiente de uma função

escalar 𝜙, denominado Potencial de Velocidade. Para melhor entendimento deste

potencial de velocidade, considere dois pontos a e b distantes um comprimento s um do

outro em um plano bidimensional. Pode-se definir uma intensidade 𝜆 variante de

acordo com o comprimento. Ao se pensar em um comprimento infinitesimal ds tem-se

que o potencial de velocidade induzido em um ponto por uma fonte de comprimento ds

é descrito por Anderson (1985):

𝑑𝜙 =𝜆𝑑𝑠

2𝜋ln 𝑟 (3.31)

Portanto, para todo o comprimento de corpo que se pretende encontrar o potencial

de velocidade induzido pela distribuição de fontes, deve-se integrar a formulação acima

(3.51) nos limites a e b que compõem o comprimento s.

𝜙(𝑥, 𝑦) = ∫𝜆𝑑𝑠

2𝜋ln 𝑟

𝑏

𝑎 (3.32)

30

Estende-se este conceito para um corpo arbitrário e considera-se a seção de tal

corpo como sendo o plano bidimensional no qual se calcula o potencial de velocidade.

Para este corpo arbitrário existem np painéis que discretizam o perímetro formado, onde

cada painel possui uma intensidade de fonte diferente (𝜆1, 𝜆2, 𝜆3, … , 𝜆𝑗 , … , 𝜆𝑛𝑝). Como o

corpo está imerso em um fluido com escoamento uniforme U (conforme descrito nas

hipóteses simplificadoras) e a superfície do corpo é coincidente com uma linha de

corrente do escoamento resultante, sem que haja escoamento no interior do corpo, o

problema se resume a encontrar valores apropriados para todos 𝜆𝑗 , de modo a atender a

as condições de impenetrabilidade e não escorregamento nos pontos de controle.

Para um ponto P(x,y) dentro do domínio estudado, tem-se 𝑟𝑝𝑗 =

√(𝑥 − 𝑥𝑗)2

+ (𝑦 − 𝑦𝑗)2 como a distância de qualquer ponto de controle do painel j

para este ponto P.

O potencial de velocidade induzido pelo painel j neste ponto P é dado por:

𝛥𝜙𝑗 =𝜆𝑗

2𝜋∫ 𝑙𝑛 𝑟𝑝𝑗𝑗

𝑑𝑠𝑗 (3.33)

Como a singularidade escolhida para o estudo tem distribuição constante de

fontes, o potencial é impactado apenas pelo painel j em questão, e o potencial em P

gerado por todos os painéis é o somatório dos potenciais de velocidade induzidos por

cada um dos painéis, conforme equação (3.34) abaixo:

𝜙(𝑃) = ∑ 𝛥𝜙𝑗𝑛𝑝𝑗=1 = ∑

𝜆𝑗

2𝜋∫ 𝑙𝑛 𝑟𝑝 𝑗𝑗

𝑑𝑠𝑗𝑛𝑝𝑗=1 (3.34)

Ao se considerar em um ponto P coincidente com o ponto de controle do painel i

nas coordenadas 𝑥𝑐𝑖 e 𝑦𝑐𝑖, obtem-se uma fórmla similar à supracitada, porém com 𝑟𝑖𝑗 =

√(𝑥𝑐𝑖 − 𝑥𝑗)2

+ (𝑦𝑐𝑖 − 𝑦𝑗)2. A equação (3.35) trata, mais especificamente, dos

potenciais de velocidade induzidos em cada um dos pontos de controle definidos.

𝜙(𝑥𝑐𝑖, 𝑦𝑐𝑖) = ∑𝜆𝑗

2𝜋∫ 𝑙𝑛 𝑟𝑖 𝑗𝑗


31

3.2.4.2 Velocidades Induzidas

Após o cálculo do escoamento potencial mencionado na subseção acima, utiliza-

se de uma operação de derivada para aferição das velocidades induzidas. Para a

componente induzida da velocidade normal ao painel 𝑢𝑓𝑛 gerada nos pontos de

controle, tem-se que:

𝑢𝑓𝑛(𝑥𝑐𝑖, 𝑦𝑐𝑖) =𝜕

𝜕𝑛𝑖[𝜙(𝑥𝑐𝑖, 𝑦𝑐𝑖)] (3.36)

Ao substituir o valor de 𝜙 da equação (3.35) e realizar a operação, a velocidade

terá como parâmetros 𝜆, x e y. A excepcionalidade deste caso é em i=j em que sua

contribuição no próprio painel é 𝜆𝑖

2 (ANDERSON, 1985).

𝑢𝑓𝑛(𝑥𝑐𝑖, 𝑦𝑐𝑖) =𝜆𝑖

2+ ∑

𝜆𝑗

2𝜋∫

𝜕

𝜕𝑛𝑖𝑙𝑛 𝑟𝑖 𝑗𝑗

𝑑𝑠𝑗𝑛𝑝𝑗=1𝑗≠𝑖

(3.37)

Considera-se, então, a hipótese simplificadora 4 – Velocidade incidente U

constante sobre o corpo – e a equação (3.6) – condição de impenetrabilidade – para

aferição de 𝑢𝑓𝑛 no ponto de controle do painel, como a componente da velocidade U

projetada pelo ângulo 𝛽𝑖 formado entre a direção do escoamento e a normal 𝑛𝑖 ao

painel.

𝜆𝑖

2+ ∑

𝜆𝑗

2𝜋∫

𝜕

𝜕𝑛𝑖𝑙𝑛 𝑟𝑖 𝑗𝑗

𝑑𝑠𝑗𝑛𝑝𝑗=1𝑗≠𝑖

+ 𝑈 cos 𝛽𝑖 = 0 (3.38)

Analogamente, para a componente tangencial faz-se a derivada em função do

vetor tangencial e posterior substituição do parâmetro 𝜙 pela equação (3.35), e obtém-

se:

𝑢𝑓𝜏(𝑥𝑐𝑖, 𝑦𝑐𝑖) =𝜕

𝜕𝜏𝑖[𝜙(𝑥𝑐𝑖, 𝑦𝑐𝑖)] = ∑

𝜆𝑗

2𝜋∫

𝜕

𝜕𝜏𝑖𝑙𝑛 𝑟𝑖 𝑗𝑗


Outra vez, sob a condições das hipóteses simplificadoras, juntas à equação (3.7) –

condição de não escorregamento –, projeta-se a velocidade U para se encontrar 𝑢𝑓𝜏 no

ponto de controle do painel.

32

∑𝜆𝑗

2𝜋∫

𝜕

𝜕𝜏𝑖𝑙𝑛 𝑟𝑖 𝑗𝑗

𝑑𝑠𝑗𝑛𝑝𝑗=1 + 𝑈 sin 𝛽𝑖 = 0 (3.40)

Conforme mencionado anteriormente, o problema pode ser resolvido ao se

encontrar valores adequados de 𝜆𝑗. Para tal, isola-se parcelas que contenham o

parâmetro 𝜆 nas equações (3.38) e (3.40), e obtém-se, respectivamente:

𝜋𝜆𝑗 + ∑ 𝜆𝑗𝐼𝑓𝑛(𝑖, 𝑗)𝑛𝑝𝑗=1𝑗≠𝑖

= −2𝜋𝑈 cos 𝛽𝑖 (3.41a)

∑ 𝜆𝑗𝐼𝑓𝜏(𝑖, 𝑗)𝑛𝑝𝑗=1 = −2𝜋𝑈 sin 𝛽𝑖 (3.41b)

onde 𝐼𝑓𝜏(𝑖, 𝑗) é a integral, multiplicada por 2𝜋, que compõe uma matriz para o

calculo da velocidade normal induzida pelo painel j em um ponto de controle i. A

integral 𝐼𝑓𝜏(𝑖, 𝑗) é análoga para a velocidade tangencial induzida em um ponto de

controle.

Pelo princípio de conservação de massa, a soma das intensidades das fontes deve

ser zero. Portanto, os somatórios dos valores de 𝜆𝑗 multiplicados pelo tamanho do painel

𝑆𝑝𝑗 deve ser zero.

∑ 𝜆𝑗𝑆𝑝𝑗𝑛𝑝𝑗=1 = 0 (3.42)

As grandezas relevantes para o cálculo das integrais 𝐼𝑓𝑛(𝑖, 𝑗) e 𝐼𝑓𝜏(𝑖, 𝑗) são

expostas na Figura 7 abaixo e sua solução pode ser descrita como o sistema de equações

(3.43).

33

Figura 7 – Grandezas para o calculo de 𝐼𝑓𝑛 e 𝐼𝑓𝜏

𝐼𝑓𝑛(𝑖, 𝑗) = ∫𝐶𝑠𝑗+𝐷

𝑠𝑗2+2𝐴𝑠𝑗+𝐵

𝑑𝑠𝑗𝑆𝑝𝑗

0 (3.43)

𝐴 = −(𝑥𝑐𝑖 − 𝑋𝑝𝑗) cos 𝜑𝑗 − (𝑦𝑐𝑖 − 𝑌𝑝𝑗) sin 𝜑𝑗 (3.43a)

𝐵 = (𝑥𝑐𝑖 − 𝑋𝑝𝑗)2

+ (𝑦𝑐𝑖 − 𝑌𝑝𝑗)2 (3.43b)

𝐶 = sin(𝜑𝑖 − 𝜑𝑗) (3.43c)

𝐷 = (𝑦𝑐𝑖 − 𝑌𝑝𝑗) cos 𝜑𝑖 − (𝑥𝑐𝑖 − 𝑋𝑝𝑗) sin 𝜑𝑖 (3.43d)

𝑆𝑝𝑗 = √(𝑋𝑝𝑗+1 − 𝑋𝑝𝑗)2

+ (𝑌𝑝𝑗+1 − 𝑌𝑝𝑗)2 (3.43e)

𝐸 = √𝐵 − 𝐴2 = (𝑥𝑐𝑖 − 𝑋𝑝𝑗) sin 𝜑𝑗 − (𝑦𝑐𝑖 − 𝑌𝑝𝑗) cos 𝜑𝑗 (3.43f)

Com este sistema de equações, obtém-se o resultado para a integral da parcela

normal induzida como:

𝐼𝑓𝑛(𝑖, 𝑗) =𝐶

2𝑙𝑛 (

𝑆𝑝𝑗2+2𝐴𝑆𝑝𝑗+𝐵

𝐵) +

𝐷−𝐴𝐶

𝐸(tan−1 𝑆𝑝𝑗+𝐴

𝐸− tan−1 𝐴

𝐸) (3.44)

De forma análoga, para parcela tangencial tem-se:

34

𝐼𝑓𝜏(𝑖, 𝑗) =𝐷−𝐴𝐶

2𝐸𝑙𝑛 (


𝐵) − 𝐶 (tan−1 𝑆𝑝𝑗+𝐴


𝐸) (3.45)

Pode-se descrever, também, as velocidades induzias normal 𝑢𝑓(𝑥, 𝑦) e tangencial

𝑣𝑓(𝑥, 𝑦) através de um sistema de equações conforme equações (3.46) e (3.47) abaixo:

𝑢𝑓(𝑥, 𝑦) = 1 + ∑𝜆𝑗

2𝜋(

𝐶

2𝑙𝑛 (


𝐵) +

𝐷−𝐴𝐶



𝐸))

𝑛𝑝𝑗=1 (3.46)

𝐴 = −(𝑥 − 𝑋𝑝𝑗) cos 𝜑𝑗 − (𝑦 − 𝑌𝑝𝑗) sin 𝜑𝑗 (3.46a)

𝐵 = (𝑥 − 𝑋𝑝𝑗)2

+ (𝑦 − 𝑌𝑝𝑗)2 (3.46b)

𝐶 = −cos(𝜑𝑗) (3.46c)

𝐷 = 𝑥 − 𝑋𝑝𝑗 (3.46d)


+ (𝑌𝑝𝑗+1 − 𝑌𝑝𝑗)2 (3.46e)

𝐸 = √𝐵 − 𝐴2 = (𝑥 − 𝑋𝑝𝑗) sin 𝜑𝑗 − (𝑦 − 𝑌𝑝𝑗) cos 𝜑𝑗 (3.46f)

𝑣𝑓(𝑥, 𝑦) = ∑𝜆𝑗

2𝜋(

𝐶

2𝑙𝑛 (


𝐵) +

𝐷−𝐴𝐶



𝐸))

𝑛𝑝𝑗=1 (3.47)

𝐴 = −(𝑥 − 𝑋𝑝𝑗𝑑𝑔) cos 𝜑𝑗 − (𝑦 − 𝑌𝑝𝑗) sin 𝜑𝑗 (3.47a)

𝐵 = (𝑥 − 𝑋𝑝𝑗)2

+ (𝑦 − 𝑌𝑝𝑗)2 (3.47b)

𝐶 = −sin(𝜑𝑗) (3.47c)

𝐷 = 𝑦 − 𝑌𝑝𝑗 (3.47d)


+ (𝑌𝑝𝑗+1 − 𝑌𝑝𝑗)2 (3.47e)


3.2.4.3 Somatório do Escoamento Incidente com as Velocidades Induzidas

Por fim, calcula-se as componentes horizontal e vertical da velocidade como a

soma do escoamento incidente com a velocidade induzida dos painéis fontes e a

35

contribuição dos vórtices, e podem ser descritas, respectivamente, pelos sistemas de

equações (3.48) e (3.49):

𝑢(𝑥, 𝑦) = 1 + ∑𝜆𝑗

2𝜋(

𝐶

2𝑙𝑛 (


𝐵) +

𝐷−𝐴𝐶



𝐸))

𝑛𝑝𝑗=1 −

1

2𝜋∑ 𝛤𝑗

(𝑦−𝑦𝑗)

(𝑥−𝑥𝑗)2

+(𝑦−𝑦𝑗)2

𝑁𝑣𝑗=1

(3.48)


𝐵 = (𝑥 − 𝑋𝑝𝑗)2

+ (𝑦 − 𝑌𝑝𝑗)2 (3.48b)

𝐶 = −cos(𝜑𝑗) (3.48c)

𝐷 = 𝑥 − 𝑋𝑝𝑗 (3.48d)


+ (𝑌𝑝𝑗+1 − 𝑌𝑝𝑗)2 (3.48e)


𝑣(𝑥, 𝑦) = ∑𝜆𝑗

2𝜋(

𝐶

2𝑙𝑛 (


𝐵) +

𝐷−𝐴𝐶



𝐸))

𝑛𝑝𝑗=1 +

1

2𝜋∑ 𝛤𝑗

(𝑥−𝑥𝑗)

(𝑥−𝑥𝑗)2

+(𝑦−𝑦𝑗)2

𝑁𝑣𝑗=1

(3.49)


𝐵 = (𝑥 − 𝑋𝑝𝑗)2

+ (𝑦 − 𝑌𝑝𝑗)2 (3.49b)

𝐶 = −sin(𝜑𝑗) (3.49c)

𝐷 = 𝑦 − 𝑌𝑝𝑗 (3.49d)

36


+ (𝑌𝑝𝑗+1 − 𝑌𝑝𝑗)2 (3.49e)


3.2.4.4 Condição de impenetrabilidade e de não-escorregamento

A condição de impenetrabilidade do painel deve levar em conta a influência dos

demais vórtices gerados, np, e espalhados ao longo do escoamento. Para tal, adapta-se a

equação (3.41a), adicionando mais uma parcela e obtém-se:

𝜋𝜆𝑖 + ∑ 𝜆𝑗𝐼𝑓𝑛(𝑖, 𝑗)𝑛𝑝𝑗=1𝑗≠𝑖

+ ∑ 𝛤𝑗𝐼𝑣𝑛(𝑖, 𝑗)𝑛𝑝𝑗=1𝑗≠𝑖

= −2𝜋𝑈 cos 𝛽𝑖 (3.50)

Analogamente, o mesmo é feito para a equação (3.41b), de onde obtém-se a

equação de não escorregamento:

∑ 𝜆𝑗𝐼𝑓𝜏(𝑖, 𝑗)𝑛𝑝𝑗=1 + ∑ 𝛤𝑗𝐼𝑣𝜏(𝑖, 𝑗)

𝑛𝑝𝑗=1 = −2𝜋𝑈 sin 𝛽𝑖 (3.51)

3.2.5 Implementação Numérica

Para colocar em prática o Método de Vórtice Discretos e solucionarmos a

Equação de Transporte de Vorticidade – equação (3.24) –, foi feito um ajuste do

algoritmo utilizado em Guedes (2003). Este algoritmo utiliza o Método dos Painéis para

geração da vorticidade e pode ser entendido de forma resumida como o fluxograma

apresentado na Figura 8 abaixo.

37

Figura 8 – Fluxograma do algoritmo implementado

Este processo foi aplicado para em uma seção retangular de razão de aspecto igual

à 1 para validação do algoritmo. Posteriormente, o mesmo foi aplicado para uma seção

de triângulo equilátero com rotações de 00˚, 30˚, 60˚, para melhor representação de uma

torre anemométrica exposta às condições de escoamento atmosférico.

3.2.5.1 Discretização do Corpo

A discretização de todos os corpos – 1ª etapa do fluxograma acima – segue

explícito para um polígono arbitrário no subcapítulo 3.3.4 – Método de Vórtices. A

particularidade entre o retângulo e o triângulo se dá pelo número de lados e pela

convergência do algoritmo com o número de painéis gerados.

Para o retângulo, os melhores resultados foram obtidos com 56 painéis, enquanto

para um triângulo números abaixo de 18 painéis careciam de detalhes e números acima

de 18 painéis não apresentaram bons resultados.

Discretização do Corpo

Cálculo da Matriz de Influência Geração de

Vorticidade

Convecção dos Vórtices

Difusão dos Vórtices

Cálculo do Coeficiente de Pressão

sobre o Corpo

Cálculo das Forças sobre

o CorpoAvanço no

Tempo

Conferência de Critèrios de Parada

38

3.2.5.2 Cálculo da Matriz de Influência

A segunda etapa do fluxograma é o cálculo da matriz de influência. Esta etapa é

resolvida pela solução de um sistema que visa o atendimento às condições de não-

escorregamento e de impenetrabilidade – equações (3.50) e (3.51). O esquema

apresentado foi feito a partir de Guedes, Bodstein e Hirata (2002) e solucionado por um

algoritmo de eliminação de Gauss.

∑ (𝐼𝑓𝑛(𝑖, 𝑗)𝜆𝑗(𝑡) + 𝐼𝑣𝑛(𝑖, 𝑗)𝛤𝑗(𝑡)) = 𝑏𝑖(𝑡),𝑁𝑗=1 1 ≤ 𝑖 ≤ 𝑁 (3.52a)

∑ (𝐼𝑓𝜏(𝑖, 𝑗)𝜆𝑗(𝑡) + 𝐼𝑣𝜏(𝑖, 𝑗)𝛤𝑗(𝑡)) = 𝑏𝑖(𝑡),𝑁𝑗=1 𝑁 + 1 ≤ 𝑖 ≤ 2𝑁 − 2 (3.52b)

∑ 𝜆𝑗(𝑡) = 𝑏2𝑁−1(𝑡)𝑁𝑗=1 (3.52c)

∑ 𝛤𝑗(𝑡) = 𝑏2𝑁(𝑡)𝑁𝑗=1 (3.52d)

3.2.5.3 Geração de Vorticidade

Os vórtices são gerados a uma distância 𝜎0 do corpo de acordo com a Figura 2

com uma intensidade 𝛤𝑗. Considera-se também o sistema apresentado na Figura 7 com

coordenadas 𝑥𝑐𝑗 e 𝑦𝑐𝑗. Então, para o cálculo da velocidade induzida no ponto de

controle do painel j (𝑥𝑐𝑗 , 𝑦𝑐𝑗) por um vórtice gerado a uma distancia 𝜎0 do painel, tem-

se a equação:

|𝑢𝑗(𝑥𝑐𝑗, 𝑦𝑐𝑗)| =𝛤𝑗

2𝜋𝜎0 (3.53)

3.2.5.4 Convecção dos Vórtices

Considera-se a convecção sem a presença de difusão, logo, ela pode ser descrita

por:

𝜕𝜔

𝜕𝑡+ 𝒖. ∇𝝎 = 0 (3.54)

39

Têm-se que, conforme mencionado no subcapítulo 3.3.4.2 – Velocidades

Induzidas – e posteriormente desenvolvido no subcapítulo 3.3.4.3 – Cálculo das

Velocidades –, a velocidade pode ser representada pela soma do escoamento incidente

com a velocidade induzida pelos painéis fontes e a contribuição dos vórtices. Como o

escoamento incidente tem apenas componente horizontal de modulo U, o sistema que

descreve as velocidades escalares u e v é:

(𝑢𝑣

) = (𝑈0

) + (𝑢𝑓

𝑣𝑓) + (

𝑢𝑣

𝑣𝑣) (3.55)

Foi escolhido realizar o esquema de Adams-Bashforth de 2ª ordem, pois o mesmo

apresenta erros menores por ser de ordem 𝛥𝑡2.

A discretização de 𝑑𝑥

𝑑𝑡= 𝒖(𝑥(𝑡), 𝑡) pelo esquema de Adam-Bashforth de 2ª ordem

é dada, em coordenadas cartesianas, por:

(𝑥(𝑡 + 𝛥𝑡)𝑦(𝑡 + 𝛥𝑡)

) = (𝑥(𝑡)𝑦(𝑡)

) + 𝛥𝑡 ([1,5 𝑢(𝑡) − 0,5𝑢(𝑡 − 𝛥𝑡)]

[1,5 𝑣(𝑡) − 0,5𝑣(𝑡 − 𝛥𝑡)]) (3.56)

3.2.5.5 Reflexão de Vórtices

Para solucionar o problema numérico dos vórtices que entram no corpo, é

utilizado o Método de Reflexão apresentado em Lewis (1991). Neste método, a

velocidade induzida em um painel é redirecionada para uma imagem do painel com uma

pequena distância perpendicular.

3.2.5.6 Difusão dos Vórtices

Considera-se a difusão sem a presença de convecção, logo, ela pode ser descrita

por:

𝜕𝝎

𝜕𝑡=

1

𝑅𝑒∇2𝝎 (3.57)

𝜎2(𝑡 + 𝛥𝑡) = 𝜎2(𝑡) +4𝛥𝑡

𝑅𝑒 (3.58)

40

3.2.5.7 Cálculo do Coeficiente de Pressão sobre o Corpo

Para o cálculo do coeficiente de pressão sobre o corpo, utilizamos das

adimensionalizações mencionadas nas equações (3.12) e a decomposição utilizada em

Mustto C., Hirata e Bodstein (1998), Guedes, Hirata e Bodstein (1999) e Guedes (2003)

para obtenção de uma formulação geral do gradiente de pressão.

−1

2∇𝑝 =

𝜕𝒖

𝜕𝑡 (3.59)

O cálculo de 𝐶𝑝 tem seu desenvolvimento apresentado por Lewis (1991). O 𝐶𝑝 no ponto

m é dado por:

𝐶𝑝𝑚 = 1 + 2 (∑𝛤𝑣𝑛

∆𝑡𝑚𝑛=1 − ∑

𝛤𝑣𝑛

∆𝑡

𝑛𝑚𝑎𝑥𝑛=1 ) (3.60)

Foi utilizada a implementação do trabalho Mustto C., Hirata e Bodstein (1998),

onde a pressão máxima encontrada se torna a pressão de referência das demais pressões.

Para tal, criou-se uma variável com valor nulo e realizou-se uma varredura em todo o

corpo para encontrar 𝑝𝑠 = 𝑝𝑚𝑎𝑥 em 𝑛 = 𝑛𝑚𝑎𝑥 a cada instante de tempo.

3.2.5.8 Cálculo de Forças sobre o Corpo

O cálculo das forças é feito utilizando a integral da pressão ao longo do corpo.

3.2.5.9 Avanço no Tempo

Para o caso do cilindro de seção retangular, tem-se o avanço de tempo 𝛥𝑡 –

adimensionalizado pelas equações (3.12) – para um deslocamento puramente

convectivo 𝛥𝑙 de mesma ordem que a distância 𝛥𝑠 entre dois vórtices gerados na

proximidade do cilindro.

𝛥𝑡 = 𝑘𝛥𝑠 (3.61)

Para o cilindro retangular, tem-se que:

𝛥𝑠 =2𝑑+2𝑙

𝑁𝑣𝑔=

2

𝑁𝑣𝑔(1 +

𝑙

𝑑) (3.62)

41

onde l e d são, respectivamente, o cumprimento e a altura da seção estudada.

Como a razão de aspecto (RA) é igual à 𝑙

𝑑 , a formula do avanço temporal pode ser

escrita como:

𝛥𝑡 =2𝑘

𝑁𝑣𝑔(1 + 𝑅𝐴) (3.63)

3.2.5.10 Conferência de Critérios de Parada

O critério de parada neste trabalho é dado como input pelo usuário. O usuário

entra o número de iterações que deseja realizar, e, ao atingir este número, o programa

sai da rotina apresentada na Figura 8, gera os últimos arquivos exportados para pós-

processamento e se encerra.

42

4 RESULTADOS

4.1 VALIDAÇÃO COM O CILINDRO DE SEÇÃO RETANGULAR

A validação da metodologia e implementação explicitadas nos capítulos anteriores

foram dadas através de simulações para gerar os resultados do escoamento ao redor da

seção retangular e da seção triangular. O objetivo dessa etapa é encontrar os resultados

de modo a validar e, por conseguinte, assegurar o correto uso do algoritmo

desenvolvido para os três estudos de casos apresentados a seguir que utilizam o avanço

de 2ª ordem de Adam-Bashfort.

Pelo fato de o objeto foco do trabalho – cilindro de seção triangular – apresentar

uma menor quantidade de estudos na literatura, a forma de validação escolhida foi o

cilindro de seção retangular com razão de aspecto igual a 1. Uma similaridade entre as

formas que motivou a escolha de tal geometria para validação é a presença de cantos

vivos. Escolheu-se um número de Reynolds de 1,20 x 105 devido à disponibilidade de

resultado experimental para comparação com o cilindro triangular. O número total de

vórtices gerados por iteração foi 56, distribuídos igualmente entre os lados da seção

retangular (14 cada), com o maior refinamento da discretização nos vértices da seção,

conforme ilustração abaixo:

Figura 9 – Discretização do cilindro retangular e seus pontos de controle

‘

U

43

Para o caso de validação, foram realizadas 500 iterações, dividida em duas

análises. A primeira análise, quantitativa, compara o coeficiente de pressão (Cp) em

torno do corpo e a média dos coeficientes de arrasto e de sustentação (Cd e Cl) com

resultados apresentados na literatura acadêmica. A segunda análise, com cunho mais

qualitativo, mira a esteira gerada na última iteração por meio do julgamento do

comportamento físico apresentado.

Para o coeficiente de sustentação, espera-se um comportamento oscilatório em

torno do zero pela simetria da seção estudada. Já para o coeficiente de arrasto, uma

comparação com resultados de outros experimentos numéricos e empíricos se faz

necessária para aferir a precisão do modelo desenvolvido.

Primeiramente, comparam-se as médias do coeficiente de pressão, de arrasto e de

sustentação, obtidos através da simulação, com o presente na literatura acadêmica para

escoamentos de mesmo tipo. O gráfico presente na Figura 10 mostra a distribuição do

coeficiente de pressão ao redor da superfície da seção e compara com o resultado

experimental de Lee (1975).

Figura 10 – Valores de Cp ao redor do corpo para um cilindro circular: (a) simulação numérica

realizada (b) experimental (LEE, 1975)

Os gráficos que detalham os demais coeficientes ao longo das iterações são

apresentados na Figura 11, onde podemos perceber uma correta oscilação do coeficiente

de arrasto em torno do zero. Outro ponto interessante a ser notado é a amplitude da

oscilação do coeficiente de arrasto com o passar do tempo de simulação.

(a) (b)

44

Figura 11 – Coeficientes de arrasto e sustentação ao longo do tempo de simulação para o

cilindro retangular

Em seguida, obtém-se os coeficientes médios para o arrasto e a sustenção. Os

valores encontrados para a simulação Cd = 1,966 e Cl = 0,044. Ambos são comparados

com valores obtidos na literatura na tabela abaixo.

Tabela 1 – Comparação de resultados para o cilindro retangular

Fonte Cd Std

Simulação Adam-Bashfort 1,966 0,127

Guedes (2003) [Num.] 1,99 0,129

Vickery (1966) [Exp.] 2,05 0,118

Bearman e Trueman (1972) [Exp.] 2,17 0,124

Vê-se, pelos resultados obtidos, que a simulação de Adam-Bashfort se aproximara

de valores numéricos e experimentais obtidos em outros estudos e, portanto, pode ser

utilizada para a comparação qualitativa do escoamento no entorno do corpo.

Para a comparação da física do problema, utilizou-se como base o trabalho

numérico de Tamura, Ohta e Kuwahara (1990), apresentado na Figura 12. Isso foi feito,

pois no estudo foi utilizada a mesma razão de aspecto e o mesmo número de Reynolds.

Através desse experimento, podemos aferir, visual e qualitativamente, se o

comportamento da esteira condiz com o fenômeno físico esperado para o cilindro de

seção retangular.

45

Figura 12 – Linhas de corrente no entorno de um cilindro retangular à Re = 105 (TAMURA,

OHTA e KUWAHARA, 1990)

As simulações de 2ª ordem gerada conforme metodologia apresentada é

apresentada nas Figura 15 e Figura 13.

Figura 13 – Campo de velocidade para simulação do cilindro retangular

46

Figura 14 – Campo de velocidade com as linhas de corrente para simulação do cilindro

retangular

Figura 15 – Dispersão de vórtices na esteira do cilindro retangular

47

Nela, pode-se perceber um comportamento similar das linhas de corrente ao

trabalho utilizado como referência. Percebe-se também a presença de recirculações ao

longo da esteira.

Do exposto, concluímos que o método apresentado é válido e pode ser utilizado

para a continuação do trabalho – estudo do escoamento em cilindros de seção triangular.

Pode-se verificar também que o número de iterações foi superdimensionado e que para

obtenção de resultados satisfatórios, sem que erros numéricos interfiram no campo de

vórtices, este número pode ser reduzido.

4.2 RESULTADOS PARA O CILINDRO DE SEÇÃO TRIANGULAR

Aqui são apresentadas as simulações realizadas para escoamentos incompressíveis

com Re = 1,20 x 105 ao redor de um cilindro de seção triangular utilizando, a

modelagem e a implementação descritas no Capítulo 3. O objetivo deste capítulo é

testar se a metodologia aplicada ao cilindro retangular pode ser utilizada no cilindro

triangular com sucesso. Ambas as seções apresentam cantos vivos e, portanto, utilizou-

se a mesma implementação.

Para o estudo do escoamento no entorno da seção triangular, foram feitas

simulações para três posições diferentes do cilindro triangular, conforme explicado no

capítulo anterior, o que gerou 3 casos diferentes – seções a 0°, 30° e 60°. Dentro de cada

um desses casos, estudou-se a implementação do avanço de 2ª ordem de Adam-

Bashford. A Figura 16 demonstra de forma visual o ângulo de incidência de cada caso.

Figura 16 – Ângulos de incidência dos cilindro triangulares

Por fim, implementação numérica do cilindro de seção triangular foi sumarizada e

comparada aos resultados experimentais de anemometria de fio quente expostos em

Iungo e Buresti (2009).

‘

UU

00 30

U

60

I II III

48

Figura 17 – Resultados experimentais do escoamento ao redor de cilindros triangulares com

variação de ângulos (IUNGO e BURESTI, 2009)

4.2.1 CASO I – Seção a 0°

Nesse caso foi utilizado uma discretização da seção triangular em 18 vórtices – 6

em cada lado. Após 310 iterações, foram obtidos os seguintes gráficos dos coeficientes

de arrasto e de sustentação para o avanço de 2ª ordem.


triangular a 00º

49


cilindro triangular a 00º

Os valores da média dos coeficientes, a partir da 50ª iteração, foram: Cd = 1,66 e

Cl = -0,27 e o número de St = 0,25 na simulação com avanço de 2ª ordem de Adam-

Bashford. Pela Figura 19, é possível observar um comportamento oscilatório adequado

dos coeficientes em relação ao tempo da simulação.

É valido notar que o campo de velocidade gerado com 310 simulações condiz

com o comportamento físico esperado, principalmente nas regiões mais próximas à

seção triangular. As Figura 20 e Figura 21 demonstram a presença de recirculações na

esteira simulada para o modelo com avanço de 2ª ordem. Nessas figuras, é possível

observar as linhas de corrente, que ajudam na compreensão do fenômeno físico.

Figura 20 – Campo de velocidade para simulação do cilindro triangular a 00º

50


triangular a 00º

Figura 22 - Dispersão de vórtices na esteira do cilindro triangular a 00º

51

Conforme vemos

Figura 22, os vórtices de Von-Kármán foram bem capturados quando a metodologia do

avanço de 2ª ordem de Adam-Bashford foi utilizada, pois percebe-se uma concentração

das recirculações com zonas delimitadas.

52

4.2.2 CASO II – Seção a 30°


em cada lado. Esse caso apresenta ângulo de 30o em relação ao Caso I apresentado.

Após 300 iterações foram obtidos os seguintes gráficos dos coeficientes de arrasto e de

sustentação para o avanço de 2ª ordem.




triangular a 30º

Os valores da média dos coeficientes a partir da 50ª iteração foram: Cd = 1,44 e


53

Bashford. Pela Figura 23, é possível observar um comportamento oscilatório adequado


É válido notar que o campo de velocidade gerado com 310 simulações condiz

com o comportamento físico esperado, principalmente nas regiões mais próximas à

seção triangular. A Figura 25 e abaixo apresentam a presença de recirculações na esteira

simulada para o método com avanço 2ª ordem. Nessas figuras é possível ver também as

linhas de corrente, que ajudam na compreensão do fenômeno físico. Observa-se

também, na Figura 26, a variação da magnitude de velocidade e linhas de corrente que

oscilam. Isso condiz com o comportamento esperado devido aos vórtices de von-

Kármán.


Figura 26– Campo de velocidade com as linhas de corrente para simulação do cilindro

triangular a 30º

54

Figura 27 – Dispersão de vórtices na esteira do cilindro triangular a 30º

Mais uma vez, conforme vemos Figura 27, os vórtices de Von-Kármán foram bem

capturados quando a metodologia do avanço de 2ª ordem de Adam-Bashford foi

utilizada. Esse caso mostrou com clareza a formação de recirculações atrás de corpo.

Vê-se também que, para o avanço de 2ª ordem, os vórtices gerados não se limitam a

esteira logo atrás do corpo.

55

4.2.3 CASO III – Seção a 60°


em cada lado. Esse caso apresenta ângulo de 60º em relação ao Caso I apresentado.

Após 310 iterações foram obtidos os seguintes gráficos dos coeficientes de arrasto e de

sustentação para o avanço de 2ª ordem.




triangular a 60º

Os valores da média dos coeficientes a partir da 50ª iteração foram: Cd = 1,35 e


56

Bashford. Pela Figura 28 é possível observar um comportamento oscilatório adequado


É valido notar que o campo de velocidade gerado com 310 iterações condiz com o

comportamento físico esperado, principalmente nas regiões mais próximas à seção

triangular. As Figura 30 e Figura 31 abaixo apresentam a presença de recirculações na

esteira simulada tanto para o método com avanço de 2ª ordem. Nessas figuras é possível

ver também as linhas de corrente, que ajudam na compreensão do fenômeno físico.

Observa-se também, na Figura 31, a formação de recirculações de sentidos opostos

intercaladas e linhas de corrente que oscilam. Isso condiz com o comportamento

esperado devido aos vórtices de von-Kármán.



triangular a 60º

57

Figura 32 – Dispersão de vórtices na esteira do cilindro triangular a 60º

Novamente, para este caso, vemos com clareza a formação dos vórtices de Von-

Kármán pela Figura 32. Porém, era de se esperar que o valor de Cl fosse próximo de

zero, uma vez que o corpo apresenta uma simetria em relação ao escoamento incidente.

58

5 COMPARAÇÃO E DISCUSSÃO

Os valores de Cd, Cl e St obtidos nas simulações acima são contrapostos com os

valores experimentais obtidos em Iungo e Buresti (2009) e Alonso (2005) para uma

melhor análise dos casos individuais e da variação dos coeficientes conforme o corpo é

girado.

Figura 33 – Gráfico de comparação entre os Cd’s obtidos e os resultados apresentados em Iungo

e Buresti (2009)

Figura 34 – Gráfico de comparação entre os Cl’s obtidos e os resultados apresentados em Iungo

e Buresti (2009)

0,6

0,8

1

1,2

1,4

1,6

1,8

2

0 10 20 30 40 50 60

Cd

Graus rotacionados [o]

Cd x Graus rotacionados

Alonso (2005) h/w=3 h/w=2 h/w=1.5 h/w=1 Resultados Obtidos ESDU (1979)

-1,2

-1

-0,8

-0,6

-0,4

-0,2

0

0,2

0 10 20 30 40 50 60

Cl


Cl x Graus rotacionados

Alonso (2005) h/w=3 h/w=2 h/w=1,5 h/w=1 Resultados Obtidos

59

onde ℎ𝑤⁄ é a razão de aspecto do corpo do experimento de Iungo e Buresti

(2009).

Os valores utilizados para a plotagem do gráfico foram as médias das simulações

com 310 iterações, caso a caso, conforme indicado na apresentação dos resultados. Os

valores foram compilados e exibidos na tabela abaixo para melhor quantificar a análise

proposta.

Tabela 2 – Compilado dos resultados do cilindro retangular para diferentes ângulos

Graus [o] Cd Cl St

0 1,66 -0,27 0,25

30 1,44 -0,87 0,38

60 1,35 -0,70 0,13

Pela tabela acima e pelas

Figura 33 e Figura 34, os coeficientes resultaram em médias diferentes dos

esperados para a vasta maioria dos casos estudados. Ao analisar o gráfico comparativo,

percebe-se que na situação em que o corpo se encontra a 0o os resultados obtidos

ficaram próximos do resultado experimental, em que o Cl deveria ser bastante próximo

de zero pela simetria em relação ao escoamento. A Figura 22 nos mostra uma esteira

muito coerente com regiões de maior intensidade de vórtices bem definidas.

0,6

0,8

1

1,2

1,4

1,6

1,8

2

0 10 20 30 40 50 60

Cd


Cd x Graus rotacionados

Alonso (2005) h/w=3 h/w=2 h/w=1.5 h/w=1 Resultados Obtidos ESDU (1979)

60

No cilindro triangular rotacionado em 30o, o avanço de 2ª ordem apresenta um

comportamento um pouco mais divergente para ambos os coeficientes, conforme

percebido da Tabela 2, isso se deve ao comportamento não natural da esteira formada.

Por fim, para o caso com o corpo posicionado a 60o, vê-se um comportamento dos

valores de Cd e Cl similares aos ocorrido no caso do corpo posicionado a 30º e

comparados aos resultados obtidos em experimentos. Era de se esperar que em 60o o

corpo se comportasse de forma similar ao caso em que o corpo se situava a 0o pela

simetria em relação ao escoamento, mas houve um maior coeficiente de sustentação de

-0,70.

Vê-se que os resultados obtidos foram muito próximos ao apresentado em Alonso

(2005), com exceção do Cl para o Caso III. Isto valida as demais simulações e permite

uma comparação com a norma IEC 61400.

61

6 CONCLUSÃO

Figura 35 – Comparação entre simulações: (a) Média temporal da série de simulações a 00º

(b) Imagem apresentada na Norma IEC 61400-12-1

A Figura 35 tem como objetivo comparar as simulações realizadas e a norma IEC-

61400-12-1: Medições do desempenho de potência de aerogeradores. Para tanto, foi

feita a média da série temporal de 100 instantes de tempo diferentes simulados por meio

da implementação do avanço de Adam-Bashfort para o Caso I – cilindro de seção

triangular a 0° –, por ter apresentado resultados mais próximos da literatura acadêmica.

Com isso, obteve-se a imagem apresentada na Figura 35a, que, quando comparada ao

exposto na norma, mostra um comportamento similar.

Conforme discutido no Capítulo 1, a correta mensuração da velocidade de uma

região de prospecção de produção de energia de parques eólicos é de suma importância

uma vez que essa variável é elevada ao cubo no cálculo da potência disponível do vento.

Qualquer erro em tal medição traz consigo um impacto direto no cálculo da produção de

energia esperada e na taxa interna de retorno de um empreendimento comercial.

Este, que se baseia em Guedes (2003), adaptando o modelo utilizado para

descrever os escoamentos ao redor com seções retangulares, para um modelo que

descreva como os escoamentos se comportam quando a seção é um triângulo equilátero

em diferentes ângulos em relação ao fluxo de vento, apresentou resultados satisfatórios

mas precisa de ajustes para a obtenção de melhores de resultados de coeficientes de

arrasto, sustentação e pressão. Como a metodologia aplicada não possui alto tempo

computacional devido às suas vantagens, o objetivo é, em trabalhos futuros,

62

implementar uma ferramenta computacional para realizar simulações mais detalhadas

do escoamento (ao invés de utilizar apenas fórmula empírica sugerida pela norma) em

função da direção preferencial do escoamento no local de interesse, a fim de melhor

quantificar as incertezas nas medições de velocidade de vento em função do

posicionamento da torre.

63

7 REFERÊNCIAS

ALCÂNTARA PEREIRA, L. A. Simulação Numérica do Escoamento em Torno de

um Corpo de Forma Arbitrária Utilizando o Método de Vórtices Discretos.

EFEI, Itajubá, MG, Brasil. 1999.

ALONSO, G. Fenomenos de galope en obstaculos de section no rectangular.

Universidad Politecnica de Madrid, Madrid, Espanha. 2005.

ANDERSON, J. D. Fundamentals of Aerodynamics. New York: McGraw-Hill Book

Company, 1985.

BATCHELOR, G. K. An Introduction to Fluid Dynamics. Cambridge: Cambridge

University Press, 1967.

BEARMAN, P. W.; TRUEMAN, D. M. An investigation of the flow around retangular

cylinders. Aeronautical Quarterly, p. 229-237, 1972.

BLEVINS, R. D. Applied Fluid Dynamics Handbook. [S.l.]: Van Nostrand Reinhold

Co., 1984.

CAO, S.; HASAN, A.; SIRÉN, K. On-site energy matching indices for buildings with

energy conversion, storage and hybrid grid connections. Energy and Buildings, p.

423-438, 2013.

CELLURA, M. et al. Different energy balances for the redesign of nearly net zero

energy buildings: An Italian case study. Renewable and Sustainable Energy

Reviews, v. 45, p. 100-112, 2015.

64

CHEN, J. M.; LIU, C. Vortex shedding and surface pressures on a square cylinder at

incidence to a uniform air stream. International Journal of Heat and Fluid Flow,

v. 20, p. 592-597, 1999.

CHORIN, A. J. Numerical study of slightly viscous flow. Journal of Fluid Mechanics,

v. 57, p. 785-796, 1973.

DOCKERY, D. W.; SCHWARTZ, J.; SPENGLER, J. D. Air pollution and daily

mortality: Associations with particulates and acid aerosols. Environmental

Research, v. 59, p. 362-373, 1992.

FINLAYSON-PITTS, B. J.; PITTS JR., J. N. Tropospheric air pollution: ozone,

airborne toxics, polycyclic aromatic hydrocarbons, and particles. Science, v. 276, p.

1045-1051, 1997.

FOUQUET, R.; PEARSON, P. The Long Run Demand For Lightning: Elasticities And

Rebound Effects In Different Phases Of Economic Development. Economics of

Energy and Environmental Policy, v. 1, 2012.

FOX, R. W.; PRITCHARD, P. J.; MCDONALD, A. T. Introdução à mecânica dos

fluidos. Rio de Janeiro: LTC, 2013.

GRUBLER, A. Energy transitions research: Insights and cautionary tales. Energy

Policy, v. 50, p. 8-16, 2012.

GUEDES, V. G. Estudo Númerico do Escoamento ao Redor de Cilindros

Circulares e Retangulares Utilizando o Método de Vórtices. COPPE/UFRJ, Rio

de Janeiro, RJ, Brasil. 2003.

65

GUEDES, V. G.; BODSTEIN, G. C. R.; HIRATA, M. H. Numerical Simulation of the

Flow Around a Square Cylinder Using the Vortex Method. Encit 2002.

Caxambu: Brazilian Society of Mechanical Sciences. 2002. p. 1-9.

GUEDES, V. G.; HIRATA, M. H.; BODSTEIN, G. C. R. Vortex Method Simulation

of the Flow Around a Circular Cylinder Using the Multipole Expansion

Algorithm. International Conference on Computational Heat and Mass Transfer.

Gazimagusa: TRNC. 1999.

HÖÖK, M.; HIRSCH, R.; ALEKLETT, K. Giant oil field decline rates and their

influence on world oil production. Energy Policy, v. 37, p. 2262-2272, 2009.

IUNGO, G. V.; BURESTI, G. Experimental investigation on the aerodynamic loads and

wake flow features of low aspect-ratio triangular prisms at different wind

directions. Jounal of Fluids and Structures, v. 25, p. 1119-1135, 2009.

LANEVILLE, A.; YONG, L. Z. Mean flow patterns around two-dimensional

rectangular cylinders and their interpretation. Journal of Wind Engineering and

Industrial Aerodynamics, v. 14, p. 387-398, 1983.

LEE, B. E. The effect of turbulence on the surface pressure field of a square prism.

Journal of Fluid Mechanics, v. 69, p. 263-282, 1975.

LEONARD, B. P. Vortez Method for Flow Simulation. Journal of Computational

Physics, v. 37, p. 289-335, 1980.

LEWIS, R. I. The Prediction of Separated Flows From Bluff Bodies. Journal

Mechanical Engineering Science, v. 23, p. 1-12, 1981.

66

LEWIS, R. I. Vortex Elements Methods for Fluid Dynamic analysis of engineering

systems. Cambridge: Cambridge University Press, 1991.

MUSTTO C., A. A.; HIRATA, M. H.; BODSTEIN, G. C. R. Numerical Simulation

oof the flow around a rotating circular cylinder using vortex method. VII Encit

98. Rio de Janeiro: Brazilian Society of Mechanical Sciences. 1998. p. 31-36.

NORBERG, C. Flow around rectangular cylinders: pressure forces and wake

frequencies. Journal of Wind Engineering and Industrial Aerodynamics, v. 49,

p. 187-196, 1993.

OWEN, N.; INDERWILDI, O.; KING, D. A. The status of conventional world oil

reserves—Hype or cause for concern? Energy Policy, v. 38, p. 4743-4749, 2010.

ROSENHEAD, L. The Formation of Vortices from a Surface of Discontinuity. Proc.

Royal Soc. London Ser. A, v. 134, p. 170-192, 1931.

ROSENHEAD, L. Laminar Boundary Layers. New York: Dover Publications, 1963.

SILVA, D. F. C. Simulção numérica do Escoamento ao redor de Aerofólios Via

Método de Vórtice Associado ao Método dos Painéis. COPPE/UFRJ, Rio de

Janeiro, RJ, Brasil. 2005.

SILVÉRIO, R.; SKLO, A. The effect of the financial sector on the evolution of oil

prices: Analysis of the contribution of the futures market to the price discovery

process in the WTI spot market. Energy Economics, v. 34, p. 1799-1808, 2012.

TAMURA, T.; OHTA, I.; KUWAHARA, K. On the reliability of two dimensional

simulation for unsteady flows around a cylinder-type structure. Journal of Wind

Engineering and Industrial Aerodynamics, v. 35, p. 275-298, 1990.

67

TOKIC, D. The 2008 oil bubble: Causes and consequences. Energy Policy, v. 38, p.

6009-6015, 2010.

VICKERY, B. J. Fluctuating lift and drag on a long cylinder of square cross-section.

International Journal for Numerical methods in Fluids, v. 5, p. 485-494, 1966.

WISSINK, J. G.; RODI, W. Numerical study of the near wake of a circular cylinder.

International Journal of Heat and Fluid Flow, p. 1060-1070, 2008.

YERGIN, D. O Petróleo: uma história mundial de conquistas, poder e dinheiro. São

Paulo: Editora Paz e Terra, 2010.

Documents

ESTUDO NUMÉRICO DO ESCOAMENTO AO REDOR …monografias.poli.ufrj.br/monografias/monopoli10029663.pdfcomparados com os resultados apresentados na norma IEC 61400-12-1. vii Graduation