Abstract

In this paper a fluid-structure algorithm for numerical simulation of wind tunnel aeroelastic analyses is proposed. The flow equations are resolved using the artificial compressibility approach and the explicit two-step Taylor-Galerkin scheme in the context of the Finite Element Method (FEM), where eight-point hexahedrons with one-point quadrature and hourglass control are adopted. Turbulence is simulated using Large Eddy Simulation (LES) with the classical Smagorinsky’s sub-grid scale model. A partitioned coupling model is adopted for fluid-structure interaction considering the three-dimensional rigid body motion approach on the structural field. Numerical results obtained here are compared with wind tunnel predictions on the CAARC tall building model and bridge sectional models to validate the present scheme.

Keywords: Computational Wind Engineering, Fluid-body interaction, Finite Element Method, Large Eddy Simulation (LES).

UM MODELO DE INTERAÇÃO FLUIDO-ESTRUTURA PARA SIMULAÇÃO NUMÉRICA DE EXPERIMENTOS AEROELÁSTICOS TRIDIMENSIONAIS EM TÚNEL DE VENTO
Resumo

Neste trabalho, um algoritmo fluido-estrutura para a simulação numérica de análises aeroelásticas de túnel de vento é proposto. As equações do escoamento são resolvidas usando a abordagem de compressibilidade artificial e o esquema explícito de dois passos de Taylor-Galerkin no contexto do Método dos Elementos Finitos (MEF), onde hexaedros de oito pontos com quadratura de um ponto e controle de modos espúrios são adotados. A turbulência é simulada usando o Large Eddy Simulation (LES) com o modelo de sub-malha clássico de Smagorinsky. Um modelo de acoplamento particionado é adotado para a interação fluido-estrutura considerando a abordagem tridimensional de movimento de corpo rígido no campo estrutural. Os resultados numéricos obtidos aqui são comparados com as previsões de túnel de vento no modelo de edifício alto do CAARC e uma seção de ponte para validar o esquema atual.

Palavras-chave: Engenharia do vento computacional, Interação fluido-estrutura, Método dos Elementos Finitos, Simulação de Grandes Escalas (LES).

1. Mestrando em Engenharia Civil, Programa de Pós-Graduação em Engenharia Civil, Universidade Federal do Rio Grande do Sul (PPGEC/UFRGS), Rio Grande do Sul-RS, Brasil. funezmarcelo@gmail.com

2 Dr. (alexandre.braun@ufrgs.br). Professor. Programa de Pós-Graduação em Engenharia Civil (PPGEC). Universidade Federal do Rio Grande do Sul (UFRGS). Av. Osvaldo Aranha, 99, 3º andar, CEP 90035-190. Porto Alegre-RS, Brasil

1 INTRODUÇÃO

A modelagem em túnel de vento foi extensivamente adotada para avaliar o comportamento aeroelástico de estruturas sob a ação do vento, especialmente pontes de grande extensão e edifícios altos. Neste contexto, técnicas experimentais são geralmente utilizadas para produzir condições ideais para as características do escoamento e estruturais do problema investigado e obter previsões confiáveis sobre a resposta estrutural. No entanto, observa-se que os resultados experimentais raramente são comparados com dados de campo, considerando que as análises de túnel de vento são normalmente realizadas durante a fase de projeto da estrutura. Com a constante melhoria na tecnologia de computadores, a simulação numérica pode ser empregada atualmente para verificar previsões experimentais obtidas em túneis de vento e simular as condições esperadas na estrutura real.

Review 645855537289-image1.png
Review 645855537289-image2.png
Review 645855537289-image3.png
Figura 1: Modelo aeroelástico de pontes com grandes vãos no túnel de vento Prof. Joaquim Blessmann (UFRGS).

No campo da aeroelasticidade de pontes, modelos seccionais são freqüentemente empregados para obter informações básicas sobre a resposta estrutural sob ação do vento, onde dispositivos especiais são usados para reproduzir a rigidez estrutural e as propriedades de amortecimento de forma adequada. Os modelos seccionais são feitos de um material rígido e simulam características físicas e geométricas da seção transversal usando apenas uma fração do seu comprimento total. Relações de similaridade são impostas entre o modelo e protótipo para definir os parâmetros físicos e geométricos. As propriedades de rigidez e amortecimento podem ser simuladas usando suportes de molas imersos em um fluido viscoso (ver Figura 1), onde modelos com um e dois graus de liberdade (movimento torsional e vertical) são considerados. Uma descrição detalhada de modelos aeroelásticos para pontes com grande vão em túneis de vento pode ser encontrada em Wardlaw [1], Scanlan [2] e Sarkar et al. [3].

Sabe-se que a resposta estrutural de edificios altos submetidos à ação do vento é principalmente governada pelos modos de vibração mais baixos (ver, por exemplo, Kareem [4]; Jeary e Ellis [5]), o que permite a adoção de modelos dinâmicos mais simples. A análise dinâmica de edifícios altos é geralmente realizada utilizando modelos equivalentes aeroelásticos, onde a geometria externa do edifício é reproduzida por uma casca com material rígido, a fim de obter a distribuição das forças aerodinâmicas sobre a superfície corretamente. Para edifícios altos, o modelo de edifício alto CAARC é geralmente utilizado para calibrar dispositivos de túnel de vento, onde uma geometria simples com seção transversal retangular é considerada. Uma descrição completa sobre o modelo de edifício alto CAARC pode ser encontrada em Melbourne [6].

Review 645855537289-image4.png
Review 645855537289-image5.png
Figura 2: Modelo aeroelástico stick para edifícios altos.

Uma das primeiras tentativas de modelagem aeroelástica de edifícios altos deve-se ao modelo conhecido como stick, onde a resposta estrutural induzida pela ação do vento é simulada considerando uma estrutura articulada na base. Dois modos translacionais fundamentais são levados em conta e os efeitos de torção são geralmente negligenciados (por exemplo, Isyumov [7]). O grau de liberdade de torção também pode ser modelado, embora apenas uma forma do modo uniforme possa ser representada com uma configuração básica. Um modelo aeroelástico de stick típico é mostrado na Figura 2, onde um bastão articulado com uma casca rígida é utilizado para reproduzir as propriedades mecânicas e corrigir relações entre os modos de vibração do modelo e proptótipo, que é ligado a molas para simular a rigidez estrutural e imerso em um fluido viscoso ou um campo eletromagnético para reproduzir o amortecimento estrutural. Thepmongkorn et al. [8] propuseram o modelo aeroelástico com base em dois graus de liberdade (BHA) para avaliações da resposta estrutural de edifícios altos sob a ação do vento. Além da resposta translacional ao longo de dois graus de liberdade ortogonais, o movimento torsional também pode ser obtido (ver Figura 3). O modelo aeroelástico BHA foi posteriormente utilizado para simular movimentos acoplados em edifícios altos devido a efeitos de interferência (Thepmongkorn et al. [9]). Uma revisão crítica sobre modelos aeroelásticos para predições experimentais de vibrações induzidas pelo vento em edifícios altos foi apresentada por Zhou e Kareem [10], onde são propostos ajustes na abordagem do modo linear adotado por modelos stick para edifícios com modos de vibração complexos e distribuição de massa complexa.

Review 645855537289-image6.png
Review 645855537289-image7.png
Figura 3: Modelo aeroelástico BHA.

Para reproduzir adequadamente as análises aerolelásticas de túnel de vento, o modelo numérico precisa de algumas características importantes: um esquema numérico para escoamentos tridimensionais com modelagem de turbulência é importante, considerando que os escoamentos de vento são inerentemente turbulentos e a turbulência é um fenômeno tridimensional. Embora modelos numéricos anteriores tenham adotado a hipótese de escoamento bidimensional (ver, por exemplo, Braun e Awruch [11]; Braun e Awruch [12]), observa-se que previsões numéricas recentes foram obtidas usando modelos de escoamentos tridimensionais para reproduzir escoamentos de vento com precisão. Por outro lado, a modelagem estrutural é freqüentemente realizada considerando a abordagem de corpo rígido, o que é uma hipótese válida, visto que os deslocamentos nos exemplos utilizados são maiores que as deformações que a estrutura poderia apresentar (Braun [13]). Além disso, os modos de vibração mais baixos associados aos edifícios altos são quase lineares. Por fim, um modelo de acoplamento fluido-estrutura deve ser considerado para aplicar cargas de vento na estrutura e impor condições cinemáticas no campo de escoamento. Embora os modelos monolíticos conduzam ao acoplamento exato em termos de conservação de energia e precisão, os modelos particionados são freqüentemente adotados devido às diferenças significativas observadas entre as propriedades mecânicas do fluido e da estrutura (ver Braun e Awruch [14]). Neste caso, o escoamento e os campos estruturais são resolvidos sequencialmente, onde a conservação de energia e a precisão só podem ser obtidas de forma assintótica, à medida que o intervalo de tempo é reduzido.

No presente trabalho, o modelo de interação fluido-estrutura apresentado por Braun e Awruch [11] [12] é estendido a fim de levar em conta o escoamento tridimensional turbulento e os movimentos tridimensionais do corpo rígido. As variáveis do escoamento são obtidas resolvendo as equações de Navier-Stokes e a equação de conservação de massa com a abordagem de compressibilidade artificial, que são discretizadas usando um esquema explícito de dois passos de Taylor-Galerkin no contexto do Método dos Elementos Finitos. Elementos hexaédricos de oito nós com quadratura de um ponto e controle de modos espúrios são adotados para discretizações espaciais. Simulação de Grandes Escalas (Large Eddy Simulation - LES) é empregada para escoamentos turbulentos e um esquema de acoplamento particionado para interação fluido-estrutura é adotado. Os movimentos estruturais seguem a suposição de corpo rígido, que são levados em conta no campo do escoamento usando uma descrição cinemática Arbitrária Lagrangeana-Euleriana (ALE). Para validar o modelo numérico proposto aqui, vibrações induzidas por vórtices em um cilindro circular e a atenuação do movimento rotacional de um cilindro retangular imerso em um fluido viscoso são simuladas. Além disso, análises de flutter em um modelo seccional da Great Belt East Bridge e testes aeroelásticos no modelo de edifício alto CAARC também são realizados.

2 ANÁLISE DO ESCOAMENTO

As expressões utilizadas na análise de escoamentos constituídos de um fluido qualquer são dadas pelas equações de balanço (conservação de massa, balanço de momentum e conservação de energia), relações constitutivas e relações termodinâmicas. Neste trabalho serão considerados apenas escoamentos isotérmicos, onde se despreza a variação de temperatura. Portanto, o sistema de equações fundamentais do escoamento é formado pelas equações e relações comentadas a seguir.

Equação de Navier-Stokes:

(1)


Na expressão acima, as propriedades físicas do fluido são a massa específica ( ), a viscosidade volumétrica ( ) e a viscosidade dinâmica ( ), sendo as variáveis do escoamento dadas pela pressão termodinâmica ( ) e as componentes do vetor de velocidades do escoamento ( ). As variáveis do escoamento (vi e p) são dadas em função de suas coordenadas cartesianas no espaço segundo a direção dos eixos e em função do tempo . Ainda tem-se a força de corpo e o delta de Kronecker ( ), onde se tem caso , e se .

Equação de conservação de massa usando a hipótese de pseudo-compressibilidade:

(2)


Equação da continuidade:

(3)


Em um fluido Newtoniano, as relações entre tensão e taxa de deformação em um ponto qualquer no espaço são dadas pelas expressões abaixo:

(4)
(5)


onde são as componentes do tensor de tensões totais e as componentes do tensor de tensões viscosas. Conforme White [15] Stokes propôs uma relação entre a viscosidade volumétrica ( ) e a viscosidade dinâmica ( ) do fluido compressível, reduzindo o número de propriedades que caracterizam o campo de tensões em um fluido compressível de dois para apenas um. A relação de Stokes, conhecida como Hipótese de Stokes, é dada por:

(6)


O tensor taxa de deformação do fluido é dado por:

(7)


Neste trabalho, para a análise numérica de escoamentos turbulentos, será utilizado a Simulação de Grandes Escalas (LES – Large Eddy Simulation) com modelos de submalha de fechamento. No modelo LES, os grandes vórtices (grandes escalas de turbulência) com frequências de flutuação mais baixas são resolvidos diretamente a partir da resolução de malha existente, já para os vórtices menores (escalas inferiores à malha) são utilizados modelos de turbulência que representarão os efeitos das pequenas escalas sobre o escoamento. Neste trabalho será utilizado o modelo clássico de Smagorinsky [16] e também o modelo dinâmico proposto inicialmente por Germano et. al [17] e corrigido por Lilly [18].

Para diferenciar as grandes das pequenas escalas no modelo LES, Leonard [19] propôs a filtragem de uma variável através da convolução da mesma por uma função filtro. Neste trabalho será utilizado uma função filtro do tipo box. Para a separação das grandes e pequenas escalas utiliza-se a Equação (8)]], na sequência a Equação (9)]] aplica o produto de convolução da variável pela função de filtragem, sendo a função de filtragem apresentada na Equação (10).

e (8)


sendo a componente do vetor velocidade na direção , a componente do vetor velocidade da parcela das grandes escalas, na direção , a componente do vetor velocidade da parcela das pequenas escalas na direção , a pressão, a pressão da parcela das grandes escalas e por final a pressão da parcela das pequenas escalas.

(9)
(10)


onde é a função de convolução, a função de filtragem, a dimensão do filtro na direção i e o número de dimensões do problema. Com o uso da filtragem e aplicando o operador de média nas variáveis da equação de Navier-Stokes, a expressão de balanço de momentum é apresentada na Equação (11) e a expressão de conservação de massa na Equação (12). Para mais detalhes ver Braun [20].

(11)
(12)
(13)


sendo o tensor de tensões submalha, a viscosidade turbulenta e as componentes do tensor taxa de deformação filtrado.

Pelo modelo clássico de Smagorinsky [16], a viscosidade turbulenta é obtida pela seguinte expressão:

(14)
(15)
(16)
(17)


sendo a constante de Smagorinsky que, conforme Murakami [21], normalmente assume valores entre 0.1 e 0.25, de acordo com o escoamento e a dimensão característica associada ao filtro utilizado.

Pelo modelo dinâmico, a viscosidade turbulenta é definida como:

(18)


onde é o coeficiente dinâmico, calculado automaticamente através da evolução do escoamento, variando em função da posição no espaço e do tempo. O coeficiente dinâmico é calculado da seguinte forma:

(19)
(20)
(21)


sendo o tensor de Leonard global e o tensor da parte anisotrópica das tensões de submalha. As variáveis estão associadas a um segundo filtro, nela utiliza-se um filtro de dimensão característica superior ao primeiro filtro. Neste trabalho a segunda filtragem é realizada conforme a Equação (22):

(22)


sendo uma variável genérica correspondente as grandes escalas do primeiro filtro e o seu valor associada ao segundo filtro. Ainda tem-se como o número de nós com conectividade direta ao nó , a distância euclidiana entre os nós e , e a variável associada ao primeiro filtro no nó . O tamanho do segundo filtro corresponde à raiz cúbica do volume total de elementos comuns ao nó i. Para maiores detalhes no modelo clássico e dinâmico, recomenda-se Braun [20].

3 ANÁLISE DA ESTRUTURA

No presente trabalho os corpos tratados serão considerados a partir de uma abordagem de corpo rígido, ou seja, sem apresentar deformações na sua forma durante todo o tempo de análise. Justifica-se essa hipótese visto que os deslocamentos nos exemplos utilizados são maiores que as deformações que a estrutura poderia apresentar (Braun [13]).Os corpos estudados serão tridimensionais, podendo apresentar deslocamento e rotação nas três direções (eixos e ), como é apresentado na Figura 4.

Review 645855537289-image8.png
(1) (b)
Figura 4: Modelos rígidos tridimensionais.

Na Figura 4 estão apresentados os dois casos que serão estudados no presente trabalho. O caso da direita será utilizado para os exemplos de edifícios altos, onde o centro de rotação está localizado na base, já o caso da esquerda será utilizado para os demais exemplos. A equação de equilíbrio dinâmico é apresentada como:

(23)


onde o subscrito indica que o termo se refere a estrutura e o sobrescrito indica que o termo se refere ao centro de massa da estrutura. Ainda temos como a matriz de massa, a matriz de amortecimento, a matriz de rigidez, o vetor de forças e , e como os vetores de deslocamento, velocidade e aceleração generalizados, respectivamente.

4 INTERAÇÃO FLUIDO-ESTRUTURA

Neste trabalho será usado um esquema de acoplamento particionado convencional, tratando os meios fluido e sólido separadamente e em sequência, como mostra a Figura 5.

Review 645855537289-image9.png
Figura 5: Algoritmo para análise de problemas de interação fluido-estrutura através de um modelo particionado convencional.

(Baseado em Braun [20]).

Na Figura 5, e representam, respectivamente, os vetores de deslocamento, velocidade e aceleração relativos a estrutura, e são, respectivamente, os vetores de pressão e velocidade referentes ao fluido e representa o vetor posição da malha no fluido. As expressões no topo da figura representam as equações de compatibilidade de deslocamentos, válidas apenas para a interface . O subíndice designa a posição no tempo. Neste esquema verifica-se que a posição da estrutura está sempre atrasada um passo de tempo em relação à posição do fluido, fazendo com que as condições de compatibilidade e de equilíbrio sejam impostas apenas de forma aproximada.

4.1 CONDIÇÕES DE COMPATIBILIDADE CINEMÁTICA E EQUILÍBRIO

Nas formulações de interação fluido-estrutura a condição de compatibilidade cinemática deve ser imposta na interface sólido-fluido, fazendo com que as velocidades e acelerações do fluido e da estrutura sejam iguais nos nós da malha pertencentes à interface. Devido à condição de não deslizamento adotada para fluidos sobre contornos sólidos, tem-se que:

e (24)


onde os subscritos E e F, representam quantidades relativas a estrutura e fluido, respectivamente, e o sobrescrito I refere-se a interface. Os valores de e são obtidos através dos vetores de velocidade e aceleração obtidos no centro de massa da estrutura e conforme a abordagem cinemática de corpos rígidos, usada neste trabalho, sendo dados pelas expressões abaixo:

(25)
(26)


onde e são os vetores de velocidade e aceleração angulares da estrutura e é o vetor posição relativa definido entre um ponto qualquer da interface e o centro de massa da estrutura. A matriz de translação e sua derivada temporal são definidas como:

(27)
(28)


As condições de equilíbrio de forças devem também ser impostas na interface fluido-estrutura, porém a carga produzida pelo fluido deve ser transladada para o centro de massa do corpo imerso no escoamento, onde as equações de movimento são escritas, ou seja:

(29)


sendo o vetor de forças do escoamento avaliadas em um dado ponto sobre a interface e é o mesmo vetor, apenas avaliado no centro de massa da estrutura. Sabendo que as forças que atuam nos contornos da estrutura são equivalentes aos vetores de tração do escoamento, avaliados sobre toda interface fluido-estrutura e com sinal invertido, tem-se a seguinte expressão:

(30)


onde é a transposta da matriz de translação avaliada em um ponto da interface.

4.2 DESCRIÇÃO ALE

As equações de Navier-Stokes e de conservação de massa comentadas anteriormente estão descritas com uma formulação Euleriana. Entretanto, em problemas de interação fluido-estrutura os deslocamentos ou deformações da estrutura exigem uma formulação arbitrária Lagrangeana-Euleriana (ALE) para a região próxima ao corpo. Para o movimento da malha deve-se usar o vetor velocidade de malha , definido conforme as condições de contorno dadas nas equações abaixo:

em
(31)
em
(32)


onde é a superfície de contorno fluido-estrutura, com velocidade de malha , já é a superfície de interface entre as regiões ALE e Euleriana da malha do fluido com velocidade de malha nula. Nos demais nós a velocidade é definida por um esquema arbitrário de movimento de malha.

Neste trabalho emprega-se o mesmo esquema adotado por Teixeira [22] e posteriormente por Braun [13], onde considera-se um procedimento de suavização das velocidades da malha a partir dos valores verificados na interface. Sua determinação é feita utilizando uma velocidade ponderada a partir da distância entre os nós do domínio do fluido e as distancias da interface, como é mostrado a seguir:

sendo (i=1,...,NALE; k=1,2,3)
(33)


onde NS é o número total de superfícies de fronteira e NALE o número de nós internos, ambos em relação ao domínio ALE. Os fatores são os coeficientes de influência entre os pontos do interior do domínio e os nós de fronteira , é a distância Euclidiana entre os nós e , já é um valor arbitrário definido pelo usuário com a finalidade de atenuar a ponderação da distância sobre os valores de velocidade de malha. Normalmente, é definido como 4, como pode-se ver no trabalho de Braun [13].

5 MODELO NUMÉRICO

Neste trabalho, as equações fundamentais do escoamento são discretizadas empregando-se o método explícito de dois passos de Taylor-Galerkin, onde aplica-se primeiramente uma discretização temporal sobre o sistema de equações através de expansões em série de Taylor até termos de segunda ordem seguida da aplicação do método de resíduos ponderados de Bubnov-Galerkin no contexto do método dos elementos finitos. Maiores detalhes sobre o método explícito de dois passos de Taylor-Galerkin podem ser obtidos em Braun [20].

O acoplamento na interação fluido-estrutura no contexto do método dos elementos finitos deve ser feito considerando um elemento do domínio do fluido localizado sobre a interface, como é apresentado na Figura 6, onde percebe-se que os nós 1, 3, 5 e 7 são os nós do elemento em contato com a estrutura.

Review 645855537289-image10.png
Figura 6: Elemento do domínio do fluido em contato com a estrutura.

As equações de Navier-Stokes discretizadas pelo método dos elementos finitos, em nível de elemento, podem ser representadas por:

(34)


onde , e são as matrizes de massa, advecção e difusão em nível de elemento. Para os elementos do fluido em contato com um corpo imerso, o sistema de equações deve ser rearranjado considerando-se separadamente os nós pertencentes a interface (com o sobrescrito ) e os nós que não pertencem a interface (com o sobrescrito ), onde apenas a parte com nós pertencentes a interface relevante para a análise de interação fluido-estrutura. Portanto:

(36)


Impõem-se então as condições de compatibilidade em nível de elemento sobre a equação acima. Estas condições são dadas abaixo:

(37)
(38)


onde e são os vetores velocidade e aceleração da estrutura em um ponto da interface, respectivamente, já e são os vetores velocidade e aceleração da estrutura no centro de massa, respectivamente. e são as matrizes de translação relativas a um elemento do fluido em contato com a estrutura imersa, dadas por:

(39)


onde as matrizes e são apresentadas nas equações (27) e (28), respectivamente. Vale ressaltar que as sub-matrizes e associadas a nós que não pertencem a interface, ou seja, com o sobrescrito , tem todos os seus elementos nulos. Para uma abordagem em corpo rígido, a equação de movimento está apresentada na Equação (23). Considerando as condições de compatibilidade e equilíbrio na interface, tem-se que:

(41)


onde é o tensor de tensões geradas pelo fluido e o vetor normal a superfície, ambos avaliados no vetor de coordenadas naturais.

Sendo que:

(42)


tem-se que:

(43)


Impondo as condições de compatibilidade na interface obtém-se:

(44)


Somando a contribuição dos elementos do fluido em contato com o corpo imerso, obtém-se a seguinte equação de movimento da estrutura:

(45)


onde NELI é o número de elementos de fluido sobre a interface.

6 RESULTADOS

6.1 ANÁLISE AEROELÁSTICA SOBRE UM CORPO RETANGULAR

O problema trata de um movimento de rotação livre do corpo quando submetido a uma rotação inicial de 5º, apresentado por Sarrate et al. [23], com o propósito de verificar a atenuação do movimento rotacional do prisma retangular quando submerso em diferentes fluidos viscosos. Este problema será analisado de duas formas, a primeira será analisada com apenas um elemento no sentido da profundidade da malha, ou seja, com um escoamento bidimensional, já a segunda análise é feita com 25 elementos no sentido da profundidade, ou seja, com um escoamento tridimensional, com o intuito de observar a influência da turbulência no exemplo analisado. A Figura 7 apresenta o domínio computacional e as condições de contorno para os exemplos com escoamento bidimensional e tridimensional, além do detalhe da malha próxima ao corpo retangular.

Review 645855537289-image11.png
Figura 7: Domínio computacional e condições de contorno para (a) escoamento bidimensional e (b) escoamento tridimensional. Em (c) detalhe da malha próximo ao corpo retangular.

A malha para o escoamento bidimensional apresenta um total de 38.720 nós e 19.200 elementos, já a malha para o escoamento tridimensional possui 503.360 nós e 480.000 elementos. A menor dimensão de um elemento na malha é de 1,25x10-2 m, levando a um incremento de tempo adimensional de Δt = 1,88x10-3 s para ambos os casos. Para as propriedades do fluido tem-se a velocidade do som igual a 0,29 m/s, velocidade característica de 0,029 m/s e a frequência torcional do corpo de 0,266 s-1, a viscosidade cinemática varia conforme o número de Reynolds, que nesse caso utilizaram-se três valores de Reynolds (3,628 ; 36,28 e 362,8). A velocidade característica do escoamento é considerada como a máxima velocidade na parte de trás do corpo, conforme Sarrate et al. [23], dado por , onde é a rotação inicial do corpo, igual a 5º, T é o período de vibração do corpo, L é o comprimento característico, considerado como metade do comprimento do corpo, ou seja, igual a 1,25 m, e ωr é a frequência torcional dada por , onde I é o momento de inércia de massa e é a rigidez rotacional. A Figura 8 apresenta os deslocamentos, velocidades e acelerações angulares para o caso com (a) escoamento bidimensional e (b) escoamento tridimensional. Nota-se que em todos os casos apresentados na Figura 8, a amplitude da rotação diminui com o tempo e a amplitude diminui de forma mais rápida quanto menor for o seu valor de Reynolds. Isto ocorre, pois, a maior viscosidade está associada ao menor número de Reynolds.

Review 645855537289-image12.png
Figura 8: Deslocamento, velocidade e aceleração angular para o caso com (a) escoamento bidimensional e (b) escoamento tridimensional.

Ainda calculou-se para os dois casos a taxa de incremento/decremento da seguinte forma: , onde e correspondem aos valores de pico entre um mesmo período de oscilação. A Tabela 1 apresenta os valores de decremento para ambos os casos, assim como a diferença percentual entre os resultados, onde pode-se notar que a diferença entre os decrementos aumenta percentualmente conforme o aumento do número de Reynolds. Esse fenômeno acontece visto que, com o aumento do Reynolds a turbulência começa a causar influência no resultado. Ou seja, a análise tridimensional, diferente da análise com apenas um elemento na profundidade, trará resultados com algumas diferenças, causadas pela influência da tridimensionalidade do escoamento, e essas diferenças são acentuadas com o aumento do número de Reynolds. Os resultados estão de acordo com os apresentados por Sarrate et al. [23].

Tabela 1: Decremento médio para os casos bidimensional e tridimensional, assim como a diferença percentual entre os resultados.
Re Decremento médio (caso bidimensional) Decremento médio (caso tridimensional)
3,628 0,20567478 0,20064923 2,44%
36,28 0,03373386 0,03052540 9,51%
362,8 0,01342618 0,01012951 24,55%


6.2 ANÁLISE AEROELÁSTICA SOBRE UM CILINDRO

Para a validação do programa estudou-se o exemplo apresentado em Ahn e Kallinderis [24] e Borazjani et al. [25], onde se deseja observar o fenômeno de lock-in relacionado ao desprendimento de vórtices em torno de um cilindro. As condições de contorno, características geométricas do exemplo e a malha utilizada estão apresentadas na Figura 9.

Review 645855537289-image13.png
Figura 9: (a) Condições de contorno, características geométricas e (b) malha utilizada para o cilindro.

A malha em elementos finitos possui 13.790 nós e 6.750 elementos, sendo a menor dimensão da malha encontrada de 0,02 m, levando a um incremento de tempo de Δt = 7x10-5 s. Para as propriedades físicas do cilindro utilizou-se um valor de Reynolds de 150, massa específica do fluido de 1,0 kg/m3, comprimento característico de 1,0 m, viscosidade cinemática de 0,2 m2/s, velocidade do som de 150 m/s e velocidade do escoamento não perturbado de 30 m/s. Já para a estrutura tem-se o cilindro com uma massa de 2 kg e a rigidez vertical dada por , onde utilizou-se velocidades reduzidas entre 3 e 8, com incrementos de 1. Ressalta-se que não há amortecimento estrutural. Os deslocamentos máximos do cilindro obtidos no presente trabalho, assim como o das referências estão apresentados na Figura 10 junto com o histórico de deslocamento vertical para a velocidade reduzida 4. Observa-se que a região do lock-in - onde a frequência de desprendimento de vórtices coincide com a frequência da estrutura, resultando numa amplitude de vibração alta - ocorre entre as velocidades reduzidas 4 e 7, visto que nessa região o deslocamento do cilindro ultrapassa 10% do seu diâmetro, fora dessa região o deslocamento é reduzido drasticamente. Os resultados obtidos da região de lock-in são similares aos de Ahn e Kallinderis [24] e Borazjani et al. [25].

Review 645855537289-image14.png

Figura 10: (a) Deslocamentos máximos do cilindro e (b) histórico de deslocamento vertical para a Vred 4.

6.3 ANÁLISE AEROELÁSTICA SOBRE UMA SEÇÃO DE PONTE

Neste trabalho será analisada a seção da Great Belt East Bridge, onde será feita uma análise do fenômeno de instabilidade por flutter, através do método direto com base na observação do comportamento da estrutura para diferentes valores de velocidades reduzidas. Esse método é mostrado em Selvam et al. [26], onde determina-se a taxa de incremento ou decremento da resposta da estrutura para as diversas velocidades reduzidas. Com esses dados constrói-se um gráfico em função das velocidades reduzidas, obtendo a velocidade crítica quando a curva gerada cruza o eixo das abscissas. As características da seção da ponte e da elevação estão apresentadas na Figura 11. O domínio computacional, condições de contorno e a malha de elementos finitos estão apresentados na Figura 12. A malha de elementos finitos possui 32.490 nós e 15.940 elementos, sendo a menor dimensão da malha de 0,11 m.

Review 645855537289-image15.png
Figura 11: Características gerais da ponte Great Belt East (a) seção e (b) elevação. (Braun [13].
Review 645855537289-image16.png
Figura 12: (a) Domínio computacional, condições de contorno e (b) malha de elementos finitos.
Review 645855537289-image17.png
Figura 13: Histórico do deslocamento angular para as diferentes velocidades reduzidas.

Para o presente exemplo utilizou-se um valor de Reynolds de 1x105, massa específica do fluido de 1,32 kg/m3, dimensão característica de 31 m e velocidade do som de 337 m/s, a viscosidade cinemática, incremento de tempo e a velocidade de referência do escoamento variam conforme a velocidade reduzida utilizada, neste exemplo utiliza-se as velocidades reduzidas de 2, 4, 6 e 10. Para a estrutura tem-se o momento de inércia de massa com o valor de 2,47x106 kgm2/m, a rigidez rotacional 7,21x106 Nm/rad e a frequência natural angular de 0,272 Hz. Vale ressaltar que não há amortecimento rotacional. Os históricos do deslocamento angular para as diferentes velocidades reduzidas estão apresentados na Figura 13, já a taxa de incremento/decremento está apresentado na Figura 14.

Review 645855537289-image18.png
Figura 14: Taxa de incremento/decremento em função da velocidade reduzida.

Nota-se na Figura 13 que o deslocamento angular amortece para as velocidades reduzidas 2, 4 e 6, aumentando o amortecimento conforme o aumento da velocidade reduzida, já para a velocidade reduzida 10 ocorre o aumento do deslocamento angular, ou seja, divergindo o resultado. Observando a curva obtida na Figura 14, conclui-se que ocorre a instabilidade por flutter no momento que a curva corta o eixo das abscissas, ou seja, no momento em que a velocidade reduzida possui um valor de 8,3, correspondendo a uma velocidade dimensional de 70 m/s.

Ainda se analisou a Great Belt East Bridge com escoamento tridimensional, onde o domínio computacional, condições de contorno e detalhe para a malha no contorno da ponte apresentadas na Figura 15.

Review 645855537289-image19.png
Figura 15: (a) Domínio computacional, condições de contorno e (b) malha no contorno da ponte.

As constantes do fluido são as mesmas apresentadas para o exemplo anterior, mudando apenas as constantes da estrutura, onde tem-se o momento de inércia de massa com o valor de 7,66x107 kgm2, a rigidez rotacional 2,24x108 Nm2/rad e a frequência natural angular de 0,272 Hz. Calcularam-se as taxas de incremento/decremento dos históricos de deslocamentos angulares para o caso com escoamento tridimensional, resultando no gráfico apresentado na Figura 16, onde pode-se observar que a curva corta o eixo das abcissas num valor de velocidade reduzida de 8,42, o que leva a uma velocidade dimensional de 71 m/s. A Tabela 2 apresenta os resultados de outros trabalhos para a velocidade crítica, onde pode-se notar que a variação entre os resultados com escoamento bidimensional e tridimensional foi muito pequena, em torno de 1,4%, isso pode ser explicado pelo fato de o escoamento ser predominantemente bidimensional. Pode-se concluir que os resultados do presente trabalho estão em acordo com os apresentados pelos outros autores, com uma pequena diferença na velocidade crítica em relação as referências, principalmente aos obtidos numericamente e experimentalmente por Larsen e Walther [27], onde se tem uma diferença máxima do resultado em túnel de vento para o do presente trabalho em torno de 4%.

Review 645855537289-image20.png
Figura 16: Taxa de incremento/decremento em função das velocidades reduzidas para o caso com escoamento tridimensional.
Tabela 2: Resultados de velocidade crítica de drapejamento.
Referência Vcrit (m/s)
Presente trabalho (escoamento bidimensional) 70
Presente trabalho (escoamento tridimensional) 71
Larsen e Walther [27] (Numérico) 74
Larsen e Walther [27] (Testes em túnel de vento) 73
Braun [13] 69-73
Frandsen [28] (Numérico) 65-75
Waterson e Baker [29] (Numérico) 70
Lee et al. [30] (Numérico) 72-73


A Figura 17 apresenta dois casos de linhas de corrente para a velocidade reduzida 2, no instante t = 10s para Y = 45 m e Y = 47.5 m. Nota-se pelas figuras abaixo, que o escoamento é predominantemente bidimensional. Futuras análises devem ser feitas para verificar a influência da extensão do trecho longitudinal tomado para a simulação numérica tridimensional sobre os resultados, além do uso de condições de contorno periódicas ao invés de condições de simetria.

Review 645855537289-image21.png
Figura 17: Linhas de corrente para escoamento tridimensional para a velocidade reduzida 2, no instante t = 10s, para Y = 45 m e Y = 47.5 m.

6.4 ANÁLISE AEROELÁSTICA SOBRE UM EDIFÍCIO ALTO

Para a análise aeroelástica sobre um edifício alto escolheu-se um edifício modelo, amplamente estudado, conhecido como modelo de edifício alto CAARC (Commonwealth Advisory Aeronautical Council). Para isso, analisou-se a resposta da estrutura para 3 valores de velocidades reduzidas diferentes (VRED = 2,6 e 8). Para as dimensões do domínio de análise, utiliza-se um domínio semelhante ao apresentado por Braun e Awruch [14], apresentado na Figura 18. Detalhes da malha do domínio e no contorno do edifício estão apresentados na Figura 19.

Review 645855537289-image22.png
Figura 18: Domínio computacional e condições de contorno do CAARC. (baseado em Braun e Awruch [14]).
Review 645855537289-image23.jpeg
Figura 19: Malha em elementos finitos do (a) domínio computacional e (b) próxima ao CAARC.

A malha possui 919.959 nós e 891.250 elementos. O menor elemento da malha possui 0,55 m. Para as constantes do fluido tem-se o número de Reynolds com o valor de 8x104, massa específica de , já a viscosidade cinemática, velocidade do som, velocidade de referência e incremento de tempo variam conforme o valor da velocidade reduzida utilizada. Os dados da estrutura são os seguintes: massa do edifício de 3,888x107 kg, inércia de massa de 4,265x1011 kgm2, inércia de massa de 4,228x1011 kgm2, rigidez rotacional de 6,735x1011 Nm/rad, rigidez rotacional de 6,677x1011 Nm/rad, coeficiente de amortecimento de 0,01 e frequência (para ambas as direções) de 0,2 Hz.

Review 645855537289-image24.png
Figura 20: Histórico de deslocamentos no topo do edifício para a direção longitudinal e transversal ao escoamento para a velocidade reduzida 2.
Review 645855537289-image25.png
Figura 21: (a) Média das amplitudes obtidas normalizadas pela seção transversal (x/L) e (b) desvio padrão normalizado pela seção transversal ( )

Na Figura 20 é apresentada a resposta da estrutura nas direções longitudinal e transversal para uma velocidade reduzida Vr = 2. Na Figura 21, os resultados obtidos na presente simulação são comparados com resultados numéricos e experimentais obtidos por outros autores. Nota-se que tanto as amplitudes para a resposta longitudinal como para a resposta transversal aumentam conforme aumenta a velocidade reduzida. Ressalta-se que os resultados para os deslocamentos transversais ao escoamento das referências são das respostas r.m.s. (root mean square) da estrutura. Para a Figura 21 (a) obteve-se resultados muito parecidos aos resultados numéricos de Braun e Awruch [14] e aos resultados experimentais de Thepmongkorn et al. [8], não muito distantes da curva sugerida por Melbourne [6]. Já para os resultados apresentados na Figura 21 (b)obteve-se resultados entre os de Braun e Awruch [14] e os demais autores, ainda mais próximos aos de Braun e Awruch [14]. Portanto, pode-se concluir que os resultados obtidos no presente trabalho estão de acordo com as referências tanto para os deslocamentos longitudinais e transversais ao escoamento.

Review 645855537289-image26.png
Figura 22: Campo de pressão instantâneo para a velocidade reduzida 2 para (a) a face frontal do CAARC e (b) a face traseira do CAARC.

Na Figura 22 pode-se notar maiores resultados de pressão na face frontal, com uma distribuição de pressões de acordo com outros trabalhos, além de valores de pressão negativos na face traseira do CAARC, visto que a parte da frente é influenciada pelo escoamento se chocando com o edifício, já a parte de trás do edifício sofre influência da turbulência gerada pelo escoamento, sendo uma zona de sucção. Na Figura 23 nota-se a influência da turbulência no escoamento, com a criação de diversos vórtices em diferentes posições. Nota-se em ambas figuras um campo de pressão alto na frente do edifício e uma zona de sucção atrás do edifício. Na Figura 24 observa-se área de recirculação do escoamento próximo a base a frente do edifício, já na parte de trás nota-se a zona de sucção, assim como o efeito da turbulência, criando diversas zonas de recirculação do escoamento. Em todas figuras anteriores pode-se notar que a estrutura do escoamento prevista pelo LES é irregular e complexa. Além disso, Huang et al. [31] comenta que nem todos vórtices obtidos com o modelo LES conseguiram ser reproduzidos com outros modelos como e LK.

Review 645855537289-image27.png
Figura 23: Campo de pressão instantâneo e linhas de corrente para a velocidade reduzida 2 no instante t = 70s em diferentes valores de Z.
Review 645855537289-image28.png
Figura 24: Campo de pressão instantâneo e linhas de corrente para a velocidade reduzida 2 no instante t = 70s em diferentes valores de Y.

Padrões de circulação típicos de escoamento em torno de edifícios altos com características de camada limite são reproduzidas numericamente conforme pode-se observar nas figuras acima (ver Peterka et al. [32]), onde, por exemplo, a separação do escoamento é facilmente visualizada nas paredes laterais e no topo do edifício, além do surgimento do vórtice em formato de ferradura (horseshoe vortex) próximo ao chão. Nota-se na Figura 25 uma zona de recirculação a frente do edifício próximo ao chão, causadas por parte do escoamento (até uma altura em torno de 0,75H sendo H a altura do edifício) que ao se chocar com o edifício, é direcionada ao chão e na sequência na direção oposta ao escoamento principal, contribuindo para a criação de vórtices em formato de ferradura, semelhantes ao apresentado por Braun e Awruch [14] e comentados por Peterka et al. [32]). Já a parte do escoamento mais próximo ao topo do edifício é direcionada para cima e após passar o edifício volta a diminuir sua altitude drasticamente, como também pode ser observado nas linhas de corrente da Figura 24]]. O mesmo acontece para a parte do escoamento que escapa pelas laterais do edifício, ascendendo em direção a esteira.

Review 645855537289-image29.png
Figura 25: (a) Linhas de corrente verticais para a velocidade reduzida 2 no instante de tempo t = 100s e em (b) tem-se as linhas de corrente horizontais para a velocidade reduzida 2 e instante de tempo t = 100s para diferentes valores de Z.

7 CONCLUSÕES

Com o presente trabalho foi possível a atualização do modelo bidimensional proposto por Braun [13] para um modelo tridimensional, a fim de obter um correto tratamento dos efeitos da turbulência, que é um fenômeno intrinsicamente tridimensional. Os exemplos analisados obtiveram resultados semelhantes aos apresentados nas referências, podendo mostrar a boa precisão do código numérico desenvolvido. Para a validação realizou-se a análise aroelástica sobre um cilindro, obtendo a região de lock-in semelhante aos das referências. No estudo da seção da Great Belt East Bridge pode-se obter resultados de velocidade crítica de drapejamento semelhantes aos das referências. Ainda se analisou o escoamento bidimensional e tridimensional, obtendo resultados semelhantes, pois o escoamento simulado possui características bidimensionais, podendo ser observado nas linhas de corrente apresentadas. Para o exemplo do edifício modelo CAARC, obteve-se resultados de deslocamentos transversais e longitudinais semelhantes aos das referências, principalmente aos apresentados por Braun e Awruch [14]. Ainda pode-se observar valores altos de pressão na face frontal do edifício e valores negativos de pressão na face traseira do edifício. Observou-se, também, a recirculação do escoamento na parte frontal do edifício próximo ao chão e a criação de diversos vórtices no escoamento.

REFERÊNCIAS

[1] Wardlaw R.L. Sectional versus full model wind tunnel testing of bridge road decks. Proceedings of Indian Academy of Sciences (Engineering Sciences), 3, 177-198, 1980.

[2] Scanlan R.H. Interpreting aeroelastic models of cable-stayed bridges. Journal of Engineering Mechanics, 114:555-575, 1987.

[3] Sarkar P.P., Jones N.P., Scanlan R.H. Identification of aeroelastic parameter of flexible bridges. Journal of Engineering Mechanics, 120:1718-1742, 1994.

[4] Kareem A. Acrosswind response of buildings. Journal of Structural Engineering, 108:869-887, 1982.

[5] Jeary A., Ellis B.R. On predicting the response of tall buildings to wind excitation. Journal of Wind Engineering and Industrial Aerodynamics, 13:173-182, 1983.

[6] Melbourne W.H. Comparison of measurements on the CAARC standard tall building model in simulated model wind flows. Journal of Wind Engineering and Industrial Aerodynamics, 6:73-88, 1980.

[7] Isyumov N. The aeroelastic modeling of tall buildings. Proceedings of the International Workshop on Wind Tunnel Modelling Criteria and Techniques in Civil Engineering Applications, 373-407, 1982.

[8] Thepmongkorn S., Kwok K.C.S., Lakshmanan N. A two-degree-of-freedom base hinged aeroelastic (BHA) model for response predictions. Journal of Wind Engineering and Industrial Aerodynamics, 83:171-181, 1999.

[9] Thepmongkorn S., Wood G.S., Kwok K.C.S. Interference effects on wind-induced coupled motion of a tall building. Journal of Wind Engineering and Industrial Aerodynamics, 90:1807-1815, 2002.

[10] Zhou Y., Kareem A. Aeroelastic balance. Journal of Engineering Mechanics, 129:283-292, 2003.

[11] Braun A.L., Awruch A.M. Aerodynamic and aeroelastic analysis of bundled cables by numerical simulation. Journal of Sound and Vibration, 284:51-73, 2005.

[12] Braun A.L., Awruch A.M. Finite element simulation of the wind action over bridge sectional models: Application to the Guamá River Bridge (Pará State, Brazil). Finite Elements in Analysis and Design, 44:105-122, 2007.

[13] Braun, A.L. Um modelo para simulação numérica da ação do vento sobre seções de ponte. Tese de Mestrado – Programa de Pós-Graduação em Engenharia Civil, Universidade Federal de Rio Grande do Sul, Porto Alegre, 2002.

[14] Braun A.L., Awruch A.M. A partitioned model for fluid-structure interaction problems using hexahedral finite elements with one-point quadrature. International Journal for Numerical methods in Engineering, 79:505-549, 2009.

[15] White, F. M. Fluid Mechanics, 8th edition. McGraw Hill, New York, 2015.

[16] Smagorinksy J. General circulation experiments with the primitive equations, I, the basic experiment. Monthly Weather Review, 91:99-135, 1963.

[17] Germano M., Piomelli U., Moin P., Cabot W.H. A dynamic subgrid-scale eddy viscosity model. Physics of Fluids, A3:1760-1765, 1991.

[18] Lilly D.K. A proposed modification of the Germano subgrid-scale closure method. Physics of Fluids, A4:633-635, 1992.

[19] Leonard A. Energy cascade in large-eddy simulations of turbulent flows. Advances in Geophysics, 18A:237-248, 1975.

[20] Braun, A. L. Simulação numérica na Engenharia do Vento incluindo efeitos de interação fluido-estrutura. 2007. 300 f. Tese (Doutorado em Engenharia) – Programa de Pós-Graduação em Engenharia Civil, Universidade Federal de Rio Grande do Sul, Porto Alegre, 2007.

[21] Murakami S. Current status and future trends in computational wind engineering. Journal of Wind Engineering and Industrial Aerodynamics, 67-68:3-34, 1997.

[22] Teixeira, P. R. F. Simulação numérica da interação de escoamentos tridimensionais de fluídos compressíveis e incompressíveis e estruturas deformáveis usando o Método dos Elementos Finitos. 2001, 237. Tese (Doutorado em Engenharia) – Programa de Pós-Graduação em Engenharia Civil, Universidade Federal do Rio Grande do Sul, Porto Alegre, 2001.

[23] Sarrate J., Huerta A., Donea J. Arbitrary lagrangian-eulerian formulation for fluid-rigid body interaction. Computer Methods in Applied Mechanics and Engineering, 190:3171-3188, 2001.

[24] Ahn H.T., Kallinderis, Y. Strongly coupled flow/structure interactions with a geometrically conservative ALE scheme on general hybrid meshes. Journal of Computational Physics, 219(2):671-696, 2006.

[25] Borazjani I, GE L, Sotiropoulos F. Curvilinear immersed boundary method for simulating fluid structure interaction with complex 3D rigid bodies. Journal of Computational Physics, 227:7587-7620, 2008.

[26] Selvam R.P., Govindaswamy S., Bosch H. Aeroelastic analysis of bridges using FEM and moving grids. Wind and Structures, 5:257-266, 2002.

[27] Larsen A., Walther J.H. Aeroelastic analysis of bridge girder sections based on discrete vortex simulation. Journal of Wind Engineering and Industrial Aerodynamics, 67&68:253-265, 1997.

[28] Frandsen J.B. Numerical bridge deck studies using finite elements. Part I: flutter. Journal of Fluids Structure, 19:171–191, 2004.

[29] Waterson N.P., Baker N. Numerical prediction of flutter behaviour for long-span bridge decks. Proceedings of the Fifth International Symposium on Computational Wind Engineering, 2010.

[30] Lee N., Lee H. Baek C., Lee S. Aeroelastic analysis of bridge deck flutter with modified implicit coupling method. Journal of Wind Engineering and Industrial Aerodynamics, 155:11-22, 2016.

[31] Huang S., Li Q.S., Xu S. Numerical evaluation of wind effects on a tall steel building by CFD. Journal of Constructional Steel Research, 63(5):612-627, 2007.

[32] Peterka J.A., Meroney R.N., Kothari K.M. Wind flow patterns about buildings. Journal of Wind Engineering and Industrial Aerodynamics, 21(1):21-38, 1985.

Back to Top