UNIVERSIDADE ESTADUAL PAULISTA "JÚLIO DE MESQUITA FILHO" CAMPUS DE SÃO JOÃO DA BOA VISTA VINÍCIUS ROSOLEN ALTARNI Análise Aeroelástica de uma Seção Típica Considerando o Modelo Aerodinâmico Beddoes-Leishman São João da Boa Vista 2024 Vinícius Rosolen Altarni Análise Aeroelástica de uma Seção Típica Considerando o Modelo Aerodinâmico Beddoes-Leishman Trabalho de Graduação apresentado ao Conselho de Curso de Graduação em Engenharia Aeronáutica do Campus de São João da Boa Vista, Universidade Estatual Paulista, como parte dos requisitos para obtenção do diploma de Graduação em Engenharia Aeronáutica . Orientador: Profº Dr. Vagner Candido de Sousa São João da Boa Vista 2024 A465a Altarni, Vinícius Rosolen Análise aeroelástica de uma seção típica considerando o modelo aerodinâmico Beddoes-Leishman / Vinícius Rosolen Altarni. -- São João da Boa Vista, 2024 69 p. : il., tabs. Trabalho de conclusão de curso (Bacharelado - Engenharia Aeronáutica) - Universidade Estadual Paulista (UNESP), Faculdade de Engenharia, São João da Boa Vista Orientador: Vagner Candido de Sousa 1. Aeroelasticidade. 2. Engenharia. 3. Aeronáutica. 4. Oscilações. 5. Simulação por computador. I. Título. Sistema de geração automática de fichas catalográficas da Unesp. Biblioteca da Universidade Estadual Paulista (UNESP), Faculdade de Engenharia, São João da Boa Vista. Dados fornecidos pelo autor(a). Essa ficha não pode ser modificada. Dedico esse trabalho aos meus pais, Débora e Junior AGRADECIMENTOS Em primeiro lugar, quero agradecer a Deus por ter me possibilitado trilhar esse caminho acadêmico até aqui. Agradeço também à minha família, que me proporcionou todo o suporte necessário durante a graduação para que eu pudesse me dedicar e desenvolvê-la com excelência. Ainda, deixo meu agradecimento especial aos professores que de uma forma ou de outra contri- buíram para minha formação durante a faculdade. Um agradecimento particular ao meu orientador, Vagner, que me apresentou o tema de estudo e me ajudou a desenvolver este trabalho. Por fim, registro aqui meu mais sincero obrigado aos amigos que estiveram comigo compartilhando essa jornada: Guilherme Lellis, Lucas Tazitu, Yago Zampolli, Filipe Valentim e Felipe Mendes. A companhia e amizade de vocês foi um grato presente que sem dúvidas levarei comigo para o resto da vida. Agradeço todo o apoio mútuo, cumplicidade, alegrias e lembranças colecionadas nesse período. “A [...] [plane] is safe [...] [on the ground], but that’s not what [...] [planes] are for.” (William G. T. Shedd) RESUMO Compreender como é o comportamento de uma estrutura imersa em um escoamento, e sua resposta dinâmica quanto aos esforços inerciais, elásticos e aerodinâmicos é de suma importância e a aeroelasti- cidade é o ramo da engenharia aeronáutica responsável por isso. Asas de aviões, rotores de helicópteros e pás de geradores eólicos são exemplos de estruturas presentes no universo aeronáutico, e sujeitas a um efeito chamado flutter, de natureza particularmente catastrófica. Diversos estudos tentam representar e entender as características desse efeito, tentando prever e mitigá-lo. Um desses estudos é a metodologia semi-empírica proposta por Beddoes-Leishman, que reúne diversas características atraentes, como a utilização de aerodinâmica não estacionária e aerodinâmica não linear. Esta, por sua vez, é responsável por descrever um fenômeno conhecido como estol dinâmico, e levar em consideração na resposta os efeitos gerados por essa não linearidade, como por exemplo o desprendimento de vórtice e como ele afeta o carregamento aerodinâmico. Isso garante uma resposta bastante realista ao modelo simulado. Ainda, o modelo de Beddoes-Leishman é aplicável para elevados ângulos de ataque, é facilmente representado por equações em espaço de estado, e tem um custo computacional baixo para solução. A integração numérica das equações diferenciais ordinárias é realizada via MATLAB pelo método de Runge-Kutta de quarta ordem. Neste trabalho, é feito um estudo mais aprofundado especificamente na influência dos parâmetros não lineares na contribuição dos coeficientes totais - sustentação e momento de arfagem - em situação de flutter sob estol dinâmico num perfil NACA 0012. PALAVRAS-CHAVE: aeroelasticidade; flutter; Beddoes-Leishman; estol dinâmico; aerodinâmica não linear; aerodinâmica não estacionária; espaço de estados; Runge-Kutta de 4ª ordem. ABSTRACT Understanding the behavior of a structure immersed in a flow, and its dynamic response to inertial, elastic, and aerodynamic forces is of paramount importance. Aeroelasticity is the branch of aeronautical engineering responsible for this. Airplane wings, helicopter rotors, and wind turbine blades are examples of structures in the aeronautical realm subject to a particularly catastrophic effect called flutter. Various studies attempt to represent and understand the characteristics of this effect, trying to predict and mitigate it. One of these studies is the semi-empirical methodology proposed by Beddoes- Leishman, which includes several attractive features, such as the use of unsteady and nonlinear aerodynamics. This, in turn, is responsible for describing a phenomenon known as dynamic stall, and takes into account the effects generated by this nonlinearity, such as vortex shedding and how it affects aerodynamic loading. This ensures a very realistic response to the simulated model. Furthermore, the Beddoes-Leishman model is applicable for high angles of attack, is easily represented by state-space equations, and has a low computational cost for solving. The numerical integration of the ordinary differential equations is carried out via MATLAB using the fourth-order Runge-Kutta method. In this work, a more in-depth study is conducted specifically on the influence of nonlinear parameters on the contribution of total coefficients - lift and pitching moment - in a flutter situation under dynamic stall on a NACA 0012 airfoil. KEYWORDS: aeroelasticity; flutter; Beddoes-Leishman; dynamic stall; nonlinear aerodynamics; unsteady aerodynamics; state-space; Runge-Kutta 4th order. LISTA DE ILUSTRAÇÕES Figura 1 Diagrama de Collar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 Figura 2 Desenvolvimento Estol Dinâmico . . . . . . . . . . . . . . . . . . . . . . . . . 31 Figura 3 Dinâmica Típica do Comportamento do Estol . . . . . . . . . . . . . . . . . . 32 Figura 4 Representação da Seção Aeroelástica . . . . . . . . . . . . . . . . . . . . . . . 33 Figura 5 Fluxograma Resolução Simulação . . . . . . . . . . . . . . . . . . . . . . . . 40 Figura 6 Respostas Estados Estruturais, U < Vcrit . . . . . . . . . . . . . . . . . . . . . 43 Figura 7 Respostas Estados Estruturais, U = Vcrit . . . . . . . . . . . . . . . . . . . . . 43 Figura 8 Respostas Estados Estruturais, U > Vcrit . . . . . . . . . . . . . . . . . . . . . 44 Figura 9 Resposta h - Efeito Estol Dinâmico . . . . . . . . . . . . . . . . . . . . . . . . 44 Figura 10 Resposta θ - Aerodinâmica Linear × Não Linear . . . . . . . . . . . . . . . . 45 Figura 11 Amplitude Ângulo de Ataque Aerodinâmico × Velocidade do Escoamento . . . 45 Figura 12 Ângulo de Ataque × Tempo . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 Figura 13 Ângulo de Ataque Aerodinâmico e Ângulo de Pitch Estrutural × Tempo . . . . 46 Figura 14 Coeficientes de Força Normal e Momento de Arfagem, U = Vcrit . . . . . . . . 47 Figura 15 Ângulo de Ataque Aerodinâmico × Tempo, U = 14, 48m/s . . . . . . . . . . 48 Figura 16 Coeficientes de Força Normal e Momento de Arfagem, U = 17m/s . . . . . . 48 Figura 17 Coeficientes de Força Normal e Momento de Arfagem - Regime Permanente, U = 17m/s . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 Figura 18 Ângulo de Ataque Aerodinâmico × Tempo, U = 17m/s . . . . . . . . . . . . 50 Figura 19 Coeficientes de Força Normal e Momento de Arfagem, Induzidos por Desprendi- mento de Vórtices, U = 17m/s . . . . . . . . . . . . . . . . . . . . . . . . . 50 Figura 20 Coeficientes Totais de Força Normal e Momento de Arfagem × α, Regime Permanente, U = 17m/s . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 Figura 21 Parcelas dos Coeficientes de Força Normal e Momento de Arfagem × α, Regime Permanente, U = 17m/s . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 Figura 22 Cn × α: Regime Transiente + Permanente, U = 17m/s . . . . . . . . . . . . 53 Figura 23 Cm × α: Regime Transiente + Permanente, U = 17m/s . . . . . . . . . . . . 53 Figura 24 Coeficientes de Sustentação, U = 17m/s . . . . . . . . . . . . . . . . . . . . 54 Figura 25 Coeficiente de Sustentação × α, U = 17m/s . . . . . . . . . . . . . . . . . . 54 Figura 26 Estados Aerodinâmicos de Beddoes-Leishman, U = 17m/s . . . . . . . . . . 55 Figura 27 Contador de Tempo Efeitos Não Lineares, U = 21m/s . . . . . . . . . . . . . 56 Figura 28 x9 e α, U = 21m/s . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 Figura 29 Estados Aerodinâmicos Não Lineares BL, U = 21m/s . . . . . . . . . . . . . 57 Figura 30 Coeficientes de Força Normal e Momento de Arfagem, U = 21m/s . . . . . . 58 Figura 31 Coeficientes Totais de Força Normal e Momento de Arfagem × α, Regime Permanente . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 Figura 32 Início e Fim Fase Não Linear . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 Figura 33 Estados Aerodinâmicos Não Lineares BL, U = 21m/s, f = 0.5 . . . . . . . . 60 Figura 34 Coeficientes de Força Normal e Momento de Arfagem, U = 21m/s, f = 0.5 . 60 Figura 35 Ângulo de Ataque Aerodinâmico × ponto de separação, U = 21m/s . . . . . 61 Figura 36 Coeficiente de Sustentação × ponto de separação, U = 21m/s . . . . . . . . 61 Figura 37 Coeficiente de Sustentação e Momento NACA 0012 . . . . . . . . . . . . . . . 71 LISTA DE TABELAS Tabela 1 – Parâmetros Estruturais do Modelo . . . . . . . . . . . . . . . . . . . . . . . . . 41 Tabela 2 – Parâmetros Aerodinâmicos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 Tabela 3 – Constantes Experimentais: Modelagem Aerodinâmica Não-Linear de Beddoes- Leishman . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 Tabela 4 – Descrição Escoamento e Efeitos Carregamento Aerodinâmico . . . . . . . . . . . 51 LISTA DE ABREVIATURAS E SIGLAS BA Bordo de ataque BF Bordo de fuga BL Beddoes-Leishman LISTA DE SÍMBOLOS h Deslocamento translacional θ Ângulo de incidência do perfil c Corda b Semicorda xea Posição relativa do eixo elástico xac Posição relativa do centro aerodinâmico xθ Distância relativa adimensional entre o eixo elástico e CG CG Centro de Gravidade do perfil L Força aerodinâmica de sustentação Mα Momento de arfagem U Velocidade do escoamento K Energia cinética P Energia potencial Wnc Trabalho virtual das componentes não conservativas mt Massa total do aparato, em translação (aerofólio + fixações) m Massa do aerofólio Iθ Momento de inércia ao redor do eixo elástico kθ Rigidez torcional kh Rigidez translacional µe Razão de massas translacional e rotacional rθ Raio de giração di,j Elementos da matriz de amortecimento viscoso de Rayleigh ωθ Frequência natural rotacional ωh Frequência natural translacional ρ Densidade do ar Cl Coeficiente de sustentação bidimensional Cmea Coeficiente de momento de arfagem no eixo elástico α Ângulo de ataque q Taxa de arfagem β Fator de compressibilidade de Prandtl-Glauert M Número de Mach Vs Velocidade do som TI Razão entre corda e velocidade do som Kα Termo de modelagem Kq Termo de modelagem KαM Termo de modelagem KqM Termo de modelagem Cnα Coeficiente angular da força normal por ângulo de ataque α1 Ângulo de estol estático s1 Constante experimental s2 Constante experimental fa Ponto de separação do escoamento fb Ponto de separação do escoamento K0 Constante experimental K1 Constante experimental K2 Constante experimental Cm0 Coeficiente de momento de arfagem sob força normal nula η Percentual da força tangencial medida em relação ao esperado αE Ângulo de ataque efetivo em termos de parâmetros de circulação Cn1 Coeficiente de força normal na iminência do estol Tp Constante de tempo modelagem aerodinâmica não-linear Tf Constante de tempo modelagem aerodinâmica não-linear Tυ Constante de tempo modelagem aerodinâmica não-linear Tυl Intervalo de tempo viagem do vórtice BA-BF τυ Contador δα1 Variação máxima de α1 Cp n Coeficiente de força normal potencial Cf n Coeficiente de força normal (parcela efeito da separação) Cυ n Coeficiente de força normal (parcela efeito do desprendimento de vórtices) Cp m Coeficiente de momento de arfagem potencial Cf m Coeficiente de momento de arfagem (parcela efeito da separação) Cυ m Coeficiente de momento de arfagem (parcela efeito do desprendimento de vórtices) Cf t Coeficiente de força tangencial (parcela efeito separação) ts Passo de tempo φ Estimativas da inclinação método RK4 SUMÁRIO 1 INTRODUÇÃO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 1.1 OBJETIVOS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 1.2 ESTRUTURA DO TRABALHO . . . . . . . . . . . . . . . . . . . . . . . . . . 32 2 MODELO MATEMÁTICO . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 2.1 EQUACIONAMENTO ESTRUTURAL . . . . . . . . . . . . . . . . . . . . . . 33 2.2 EQUACIONAMENTO AERODINÂMICO - MODELAGEM DE BEDDOES- LEISHMAN . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 2.2.1 Equacionamento Aerodinâmica Linear . . . . . . . . . . . . . . . . . . . . . . 35 2.2.2 Equacionamento Aerodinâmica Não-linear . . . . . . . . . . . . . . . . . . . 36 3 METODOLOGIA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.1 MÉTODO COMPUTACIONAL . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.2 FLUXOGRAMA DE CÁLCULO . . . . . . . . . . . . . . . . . . . . . . . . . . 39 4 PARÂMETROS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 5 RESULTADOS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 5.1 VERIFICAÇÃO PRELIMINAR . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 5.2 COMPORTAMENTOS AERODINÂMICOS . . . . . . . . . . . . . . . . . . . . 45 5.2.1 Bifurcação α vs U . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 5.2.2 Ângulo de Ataque (α) vs Tempo . . . . . . . . . . . . . . . . . . . . . . . . . . 46 5.2.3 Ângulo de Ataque Aerodinâmico (α) e Ângulo de Pitch Estrutural (θ) vs Tempo 46 5.2.4 Coeficientes de Força Normal e Momento de Arfagem . . . . . . . . . . . . . 47 5.2.4.1 Resposta no Tempo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 5.2.4.2 Envelope de α . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 5.2.5 Coeficiente de Sustentação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 5.2.5.1 Resposta no Tempo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 5.2.5.2 Envelope α . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 5.2.6 Estados Aerodinâmicos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 5.2.7 Contador τv e Estados Aerodinâmicos Não Lineares . . . . . . . . . . . . . . 55 5.2.8 Coeficientes de Força Normal e Momento de Arfagem . . . . . . . . . . . . . 58 5.2.8.1 Resposta no Tempo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 5.2.8.2 Envelope de α . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 5.2.9 Fase Não Linear . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 5.2.10 Ponto de separação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 6 CONCLUSÃO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 7 CONSIDERAÇÕES FINAIS . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 REFERÊNCIAS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 APÊNDICE A – RELAÇÕES UTILIZADAS NA ADIMENSIONALIZAÇÃO DA MODELAGEM ESTRUTURAL . . . . . . . . . . . . . 69 ANEXO A – DADOS EXPERIMENTAIS NA LITERATURA NACA 0012 . 71 29 1 INTRODUÇÃO Aeroelasticidade é o ramo da engenharia aeronáutica responsável por estudar como forças aerodinâ- micas, inerciais e elásticas interagem entre si e são responsáveis por diversos fenômenos de diferentes naturezas. Essas interações são, ainda, dependente de diversos fatores geométricos e estruturais, portanto, dimensionais. Entretanto, não é incomum a utilização da adimensionalização para garantir uma abordagem mais ampla ao tema. Problemas aeroelásticos não existiriam se as estruturas fossem perfeitamente rígidas. Asas de aviões modernos são cada vez mais flexíveis, e esse é o motivo principal da existência dos diversos fenômenos aerelásticos (BISPLINGHOF; ASHLEY; HALFMAN, 1955). Por mais conveniente que seja, em diversas aplicações de engenharia a modelagem de estruturas rígidas na prática nem sempre se observa. E quando submetidas a um escoamento sobre elas, é observado que um corpo flexível está sujeito à interação da tríplice força, como é o caso das asas de aeronaves. Isso ainda é aplicável para rotores de helicópteros, geradores de energia eólica, e até mesmo presente em pontes. Figura 1 – Diagrama de Collar fonte: (COLLAR, 1978) Ainda, é de interesse da aeroelasticidade compreender o limite da estabilidade de tais fenômenos. A estabilidade estática é uma propriedade do estado de equilíbrio. A definição de equilíbrio por sua vez, no caso de um avião, por exemplo, equivale a um voo nivelado em que a resultante das forças é nula. Dessa forma, pode-se definir estabilidade estática como sendo a tendência inicial de um corpo retornar ao equilíbrio após uma perturbação (NELSON, 1997). A estabilidade estática pode ser classificada como i) estável, quando o equilíbrio é restaurado, ii) instável, quando ele não é restaurado, ou iii) neutra, quando um novo ponto de equilíbrio é encontrado. Por sua vez, a estabilidade dinâmica está relacionada com a resposta temporal de uma perturbação. Também é classificada em estável, instável ou neutra. Um corpo dinamicamente estável é aquele cujas 30 oscilações tendem a reduzir com o tempo, pois tem dissipação de energia, logo um amortecimento positivo. Já aquele dinamicamente instável, as oscilações aumentam - adição de energia no sistema. Por fim, no dinamicamente neutro, elas permanecem constantes. Vale ressaltar que um corpo pode ser estaticamente estável, mas dinamicamente instável. Estabili- dade estática não garante estabilidade dinâmica, porém, para haver estabilidade dinâmica, é preciso ter estabilidade estática. Além disso, no caso de movimentos oscilatórios, como é o caso desse trabalho, frequência e período são bastante relevantes (NELSON, 1997). Aprofundando um pouco mais as interações das forças de diferentes naturezas envolvidas nos fenômenos aeroelásticos, como apresentado na Figura 1, tem-se que para forças aerodinâmicas associadas com forças elásticas (estruturais) o que define-se por aeroelasticidade estática, como casos de distribuição de carga, eficiência de controle, divergência ou reversão do sistema de controle. A divergência por exemplo, é um caso em que as forças aerodinâmicas estáticas são elevadas de forma que a rigidez torcional da asa não é mais capaz de suportar a deformação, levando ao aumento do ângulo de incidência, que também aumenta os esforços aerodinâmicos até a falha da estrutura. Adicionando as forças inerciais nesse sistema, elevamos a complexidade dele e com isso, novos efeitos podem resultar da tríplice interação, chamada de aeroelasticidade dinâmica, sendo o flutter o principal deles. O flutter é um fenômeno de instabilidade autossustentado caracterizado por oscilações do tipo flexo- torcionais, com amplitudes crescentes. Dentre os tipos de flutter, pode-se citar o clássico, envolvendo torção e flexão, mas também o stall flutter, relacionado ainda com o estol dinâmico, em que a separação do escoamento durante o movimento oscilatório e consequentemente os vórtices desprendidos estão envolvidos na análise. Devido à natureza catastrófica do flutter, é um fenômeno relevante e com diversos estudos na área de aeroelasticidade. Stall flutter é um fenômeno aeroelástico inerentemente não linear. Como o nome sugere, isso ocorre com a separação parcial ou completa do escoamento enquanto a oscilação acontece. Diferentemente do flutter clássico (escoamento colado durante o tempo todo), o mecanismo de transferência de energia do escoamento para o aerofólio em oscilação não está no acoplamento elástico e/ou aerodinâmico entre os modos oscilação ou ainda um atraso de fase do deslocamento e reação aerodinâmica. A característica essencial do stall flutter é a reação aerodinâmica não linear. Além disso, a instabilidade básica e o comportamento principal são explicados em termos não lineares da força normal e momento de arfagem (DOWELL, 2015). Não obstante, quando o campo do escoamento é visualizado durante o efeito do stall flutter, é observado que vórtices livres são gerados nas vizinhanças dos pontos de separação. Esses vórtices se desprendem periodicamente criando regiões de redução ou até reversão da velocidade do escoamento na proximidade do aerofólio (DOWELL, 2015). O estol dinâmico tem esse nome devido ao perfil aerodinâmico sofrer uma perda na capacidade de gerar sustentação a partir do momento que o ângulo de ataque ultrapassa o ângulo de estol estático, até que retorna para um valor abaixo deste. Esse fenômeno ocorre para valores positivos quanto negativos de α. Contudo, em perfis assimétricos, o módulo dele é diferente. Modelar o comportamento de estol dinâmico e seus efeitos derivados é um trabalho bastante 31 Figura 2 – Desenvolvimento Estol Dinâmico fonte: (MULLENERS K., 2013) desafiador, e métodos de fluidodinâmica computacional podem falhar em descrever as características do carregamento aerodinâmico, dependendo do modelo de turbulência adotado. Simulações de alta fidelidade têm uma resposta melhorada, porém com custo computacional mais elevado ainda. Uma solução para isso, foi então desenvolvida pela indústria de helicópteros, que são modelos semi- empíricos. Essas abordagens modelam a física principal do problema com o auxílio de parâmetros avaliados empiricamente, e a principal vantagem disso é o baixo custo computacional (SANTOS, 2021). Outro ponto relevante que tange a aeroelasticidade dinâmica é modelagem aerodinâmica utilizada para descrever e prever tais fenômenos. Por se tratar de movimentos oscilatórios, a aerodinâmica estacionária não é apropriada para bem representar e capturar os resultados com precisão. Para baixas frequências, modelos estacionários ou quase estacionários podem ser satisfatórios, mas em elevadas frequências, especialmente para o caso de flutter, um modelo não estacionário é necessário para garantir respostas mais realistas. Além disso, o tipo de aerodinâmica envolvida na modelagem é de suma importância. Na aero- dinâmica linear, o carregamento é descrito por relações mais simplificadas, enquanto que na prática existem efeitos aerodinâmicos relevantes a serem considerados, como desprendimento de vórtices e o estol dinâmico. Isso causa, por exemplo, redução na capacidade de produção de sustentação e também estão associados com uma defasagem entre o movimento estrutural e o carregamento aerodinâmico. Dessa forma, a aerodinâmica não linear é mais apropriada e proporciona resultados mais coerentes. Entre os modelos de aerodinâmica não linear e não estacionária, o semi-empírico proposto por Beddoes-Leishman (LEISHMAN; BEDDOES, 1989) é um método que inclui correções para contabili- zar os efeitos causados pelo estol dinâmico em aerofólios com elevados ângulos de ataque (SANTOS; PEREIRA; MARQUES, 2017), de fundamental relevância para o presente trabalho. Não obstante, o equacionamento pode ser convenientemente escrito em representação do tipo espaço de estados, sendo técnicas de integração numérica facilmente aplicáveis. Por fim, da literatura, a perturbação causada pelo desprendimento de vórtices é percorrida ao longo da corda e induz uma onda de pressão na superfície do aerofólio. Essa variação na pressão provoca aumentos significativos na sustentação e na tendência de nose-down para o momento de arfagem, como apresentado na Figura 3. Para oscilações senoidais sob condições de escoamento colado, loops elípticos do coeficientes de sustentação e momento de arfagem deveriam ser obtidos. Entretanto, é possível notar distorções, especialmente para o coeficiente de momento. Esse comportamento foi mapeado como um terceiro harmônico que tem o efeito de aumentar a contribuição de pitch próximos do mínimo e máximo ângulo de ataque (LEISHMAN; BEDDOES, 1989). 32 Figura 3 – Dinâmica Típica do Comportamento do Estol fonte: (LEISHMAN; BEDDOES, 1989) 1.1 OBJETIVOS Como objetivo principal do trabalho, está a identificação da influência dos parâmetros não lineares na contribuição dos coeficientes totais - sustentação e momento de arfagem - em situação de flutter sob estol dinâmico. O objeto de estudo é um perfil NACA 0012, com dois graus de liberdade. Como objetivos secundários, estão a verificação e validação do estol dinâmico na simulação computacional em consonância com a literatura. 1.2 ESTRUTURA DO TRABALHO O trabalho está estruturado em modelagem estrutural, aerodinâmicas linear e não linear, no capítulo 2. A metodologia computacional está descrita no capítulo seguinte, e os parâmetros utilizados no capítulo 4. Depois, os resultados são apresentados e discutidos no capítulo 5. 33 2 MODELO MATEMÁTICO 2.1 EQUACIONAMENTO ESTRUTURAL Nesse capítulo, será apresentado o modelo matemático adotado no estudo, utilizando dois graus de liberdade acoplados. A Figura 4 abaixo ilustra isso. Os deslocamentos translacional e angular são representados, respectivamente, por h e θ. Figura 4 – Representação da Seção Aeroelástica fonte: (SANTOS; PEREIRA; MARQUES, 2017) - adaptado pelo autor O comprimento característico é a corda (c), e a localização de alguns pontos importantes são o centro aerodinâmico (xac), eixo elástico (xea) e centro de gravidade (CG). Os carregamentos aerodinâmicos são sustentação (L) e momento de arfagem (Mα), no centro aerodinâmico, e U é a velocidade do escoamento. Tomando como base a formulação de Lagrange proposta em (ERTURK A.; INMAN, 2011), desconsiderando os efeitos dissipativos, o princípio de Hamilton é dado por:∫ t2 t1 (δK − δP + δWnc)dt = 0 (2.1) onde K e P são a energia cinética e potencial, respectivamente. Wnc é o trabalho virtual das compo- nentes não-conservativas. A equação final de Lagrange é então representada da seguinte forma: d dt ( δK δq̇i ) + δP δqi = Qi (2.2) em que i representa cada grau de liberdade. Assim, tem-se: 34 d dt ( δK δθ̇ ) + δP δqθ = Mα d dt ( δK δḣ ) + δP δqh = −L (2.3) onde P = 1 2 khh 2 + 1 2 kθθ 2 (2.4) K = 1 2 mtḣ 2 + 1 2 Iθθ̇ 2 +mxθbḣθ̇ (2.5) (HODGES; PIERCE, 2002) em que mt é a massa translacional e m a massa do aerofólio1, Iθ o momento de inércia ao redor do eixo elástico, kθ a rigidez torcional, kh rigidez linear, b = c 2 e xθ = xea − CG. Após adimensionalização das equações aeroelásticas de movimento (SOUSA, 2016), elas são então descritas da seguinte forma: [ µe xθ xθ r2θ ]{ ¨̄h θ̈ } + [ d1,1 d1,2 d2,1 d2,2 ]{ ˙̄h θ̇ } + [ ω2 h 0 0 r2θω 2 θ ]{ h̄ θ } = ρU2 m { −Cl 2Cmea } (2.6) onde h̄ = h b é um deslocamento adimensional, µe é a razão entre a massa translacional mt e a massa do aerofólio m, rθ o raio de giração ao redor do eixo elástico, ωh e ωθ as frequências naturais de translação e rotação respectivamente, ρ a densidade do ar, Cl e Cmea os coeficientes de sustentação e momento, respectivamente, no eixo elástico, e di,j são os coeficientes de amortecimento viscoso de Rayleigh. O ângulo de ataque α(t), relacionado com o carregamento aerodinâmico pode ser expresso em função de θ(t) (relacionado com deslocamentos estruturais) de acordo com a seguinte relação: α(t) = θ(t) + arctan [ ḣ(t) U ] (2.7) Por sua vez, uma taxa de arfagem é conveniente para a modelagem aerodinâmica, e calculada no eixo elástico é dada por: q(t) = θ̇(t) c U (2.8) 2.2 EQUACIONAMENTO AERODINÂMICO - MODELAGEM DE BEDDOES-LEISHMAN A modelagem proposta por Beddoes-Leishman pode ser convenientemente expressa na represen- tação em espaço de estados. Essa formulação é expressa por oito variáveis de estado considerando aerodinâmica linear, não estacionária, e quatro considerando a contribuição não linear do carregamento aerodinâmico (SANTOS; PEREIRA; MARQUES, 2017), (CHANTHARASENAWONG, 2007). 1 Apesar de um mesmo perfil, as massas são diferentes porque no movimento translacional também é contabilizada a parcela referente à estrutura. 35 2.2.1 Equacionamento Aerodinâmica Linear Os estados aerodinâmicos lineares são dados da seguinte forma (LEISHMAN; CROUSE, 1989): {ẋ} = [A]{x}+ [B] { α q } (2.9) onde {x} é o vetor das variáveis de estados aerodinâmicas lineares (x1, . . . , x8) e [A] uma matriz diagonal dada por: [A] = diag[a11 a22 a33 a44 a55 a66 a77 a88] (2.10) com os coeficientes de acordo com a11 = −b1β 2 ( 2U c ) a22 = −b2β 2 ( 2U c ) a33 = − 1 KαTI a44 = − 1 KqTI ,  a55 = − 1 b3KαMTI a66 = − 1 b4KαMTI a77 = −b5β 2 ( 2U c ) a88 = − 1 KqMTI (2.11) onde β = √ 1−M2, M = U/Vs é o número de Mach, Vs = √ γR T é a velocidade do som, com γ = 1, 4, R = 286 J/(Kg.K), T a temperatura em Kelvin, e TI = c/Vs. Kα = 0, 75/[(1−M) + πβ2M2 (A1b1 + A2b2)] Kq = 0, 75/[(1−M) + 2πβ2M2 (A1b1 + A2b2)] KαM = (A3b4 + A4b3)/[b3b4 (1−M)] KqM = 7/[15 (1−M) + 3πβM2b5] , (2.12) com os seguintes parâmetros constantes de modelagem A1 = 0, 3, A2 = 0, 7, A3 = 1, 5, A4 = −0, 5, b1 = 0, 14, b2 = 0, 53, b3 = 0, 25, b4 = 0, 1, e b5 = 0, 5. A matriz [B] em (2.9) é dada por: [B] = [ 1 1 1 0 1 1 0 0 0.5 0.5 0 1 0 0 1 1 ]T (2.13) A equação a seguir apresenta a contribuição linear do carregamento aerodinâmico total.{ Cp n Cp m } = [C]{x}+ [D] { α q } (2.14) A modelagem é completada com os coeficientes de força normal potencial (Cp n) e momento de arfagem (Cp m), lineares, utilizados no equacionamento aerodinâmico não linear. A matriz [C] é definida da seguinte maneira: 36 [C] = [ c11 c12 c13 c14 0 0 0 0 c21 c22 0 0 c25 c26 c27 c28 ] (2.15) em que os coeficientes são dados por: c11 = Cnαβ 2 ( 2U c ) A1b1 c12 = Cnαβ 2 ( 2U c ) A2b2 c13 = − 4 M ( 1 KαTI ) c14 = − 1 M ( 1 KqTI ) c21 = c11(0, 25− xac) ,  c22 = c12(0, 25− xac) c25 = 1 M ( A3 b3KαMTI ) c26 = 1 M ( A4 b4KαMTI ) c27 = −Cnα 16 b5β 2 ( 2U c ) c28 = 7 12M ( 1 KqMTI ) (2.16) onde Cnα = 2π/β rad−1 é o coeficiente angular da força normal por ângulo de ataque. Em (2.14), a matriz [D] é tal que: [D] =  4 M 1 M − 1 M − 7 12M  (2.17) 2.2.2 Equacionamento Aerodinâmica Não-linear Para inserir as características não lineares devido aos efeitos de estol dinâmico, é preciso incluir mais estados na modelagem (x9, ..., x12). Esses novos estados dependem de constantes adimensionais de tempo, contudo, podem ser convertidos para uma abordagem de tempo dimensional ao multiplicar pelo fator U/b. O primeiro estado está relacionado com o estol no bordo de ataque, sendo uma parcela de atraso de fase em Cn, agora C ′ n, considerando a constante de tempo Tp: ẋ9 = ( U b ) Cp n(t)− x9 Tp e x9 = C ′ n(t), (2.18) Os dois estados seguintes, x10 e x11, estão relacionados com a separação do escoamento, baseados na teoria clássica de Kirchoff para escoamento bidimensional em uma placa plana. O ponto de separação f , normalizado pela corda, em função de um ângulo arbitrário α̂ é tal que: f = { 1− 0.3e(|α̂|−α1)/ s1 se |α̂| ≤ α1 0.04 + 0.66e(α1−|α̂|)/ s2 se |α̂| > α1 (2.19) em que α1 é o ângulo de estol estático, tal que f(α1) = 0, 7 e s1 e s2 são constantes obtidas experimentalmente. Sob condições dinâmicas, o atraso na separação do bordo de fuga é modelado com uma constante de tempo Tf0 , e o estado x10 relacionado com a distribuição de pressão na superfície do aerofólio (C ′ n) 37 sob um ângulo de ataque equivalente: ẋ10 = ( U b ) f(C ′ n(t) Cnα )− x10 Tf e x10 = fa(t) (2.20) onde fa(t) é o ponto de separação do escoamento, considerando o ângulo ( C′ n(t) Cnα ) . Para o estado x11, relacionado com o ponto de separação considerando o ângulo de ataque ins- tantâneo, representa variações no momento de arfagem pelo deslocamento do centro de pressão do aerofólio em decorrência do desprendimento de vórtices. ẋ11 = ( U b ) 2(f(α)− x11) Tf0 e x11 = fb(t) (2.21) sendo fb(t) similar à fa(t). Com o apresentado acima, a contribuição dos coeficientes aerodinâmicos são dadas por: Cf n(t) = Cnα ( 1 + √ fa(t) 2 )2 αE(t) Cf m(t) = {K0 +K1(1− f̂(t)) +K2sen(πf̂(t) 2)}CnααE(t) + Cm0 Cf t (t) = ηCnα √ fa(t)α 2 E(t) (2.22) em que αE = β2V b (A1b1x1 + A2b2x2) , (2.23) e f̂(t) = max[fa(t), fb(t)]. K0, K1 e K2 são constantes empíricas, Cm0 é o coeficiente de momento de arfagem quando a força normal é nula, e η representa uma porcentagem da força tangencial medida experimentalmente em condições de escoamento colado, comparado com valor previsto hipoteticamente pela teoria potencial. O último dos quatro estados não lineares está associado com a evolução do fenômeno de estol dinâmico, em decorrência do desprendimento de vórtices no bordo de ataque, percorrendo a superfície do aerofólio até sair pelo bordo de fuga. Quando |C ′ n| ≥ Cn1 , com Cn1 sendo o valor de Cn na iminência do estol sob condições estáticas, é então que a fase de desprendimento dos vórtices é provocada (sob outra condição, o escoamento é colado). Nessa condição, a parcela do coeficiente de força normal relacionada com o desprendimento dos vórtices é dado pelo estado x12 tal qual apresentado a seguir. ẋ12 = ( U b ) Ċυ − x12 Tυ e x12 = Cυ n(t) (2.24) em que Cυ = { Cnα [1− 0, 25(1 + √ fa(t)) 2αE(t)] se τυ ≤ 2Tυl 0 se τυ > 2Tυl (2.25) onde Tυl é um parâmetro experimental, definido como janela de tempo adimensional de desprendimento dos vórtices entre o bordo de ataque e bordo de fuga. 38 Durante o deslocamento do vórtice sobre a superfície do aerofólio, a posição do centro de pressão é alterada, e uma das parcelas do coeficiente do momento de arfagem é então expressa da seguinte maneira: Cυ m(t) = −0.25 [ 1− cos ( πτυ Tυl )] Cυ n (2.26) Além disso, durante a fase de desprendimento dos vórtices, as constantes de tempo Tf e Tυ e o ângulo α1 (cf. eq. (2.19)) são modificados, como segue abaixo: Tf =  Tf0 se 0 ≤ τυ ≤ Tυl e αα̇ ≥ 0 1 3 Tf0 se Tυl < τυ ≤ 2Tυl e αα̇ ≥ 0 1 2 Tf0 se 0 ≤ τυ ≤ 2Tυl e αα̇ < 0 4Tf0 se 2Tυl < τυ (2.27) Tυ =  Tυ0 se 0 ≤ τυ ≤ Tυl e αα̇ ≥ 0 1 4 Tυ0 se Tυl < τυ ≤ 2Tυl e αα̇ ≥ 0 1 2 Tυ0 se 0 ≤ τυ ≤ 2Tυl e αα̇ < 0 9 10 Tυ0 se 2Tυl < τυ (2.28) α1 = { α10 se αα̇ ≥ 0 α10 − (1− x10) 0.25δα1 se αα̇ < 0 . (2.29) onde τυ é um contador, e δα1 a variação máxima de α1. O subscrito 0 denota o valor de cada parâmetro antes do desprendimento do vórtice. Após o término da fase de desprendimento de vórtice, o escoamento é recolado (|C ′ n| < Cn1) e |α(t)| diminui, de forma que α1 = α10 , Tυ = Tυ0 , e a variação em Tf representada por: Tf = { Tf0 se fa(t) ≥ 0, 7 2Tf0 se fa(t) < 0, 7 (2.30) Diante de todo o apresentado, os coeficientes aerodinâmicos podem então ser expressos por: Cn(t) = Cp n(t) + Cf n(t) + Cυ n(t)− CnααE(t) Cm(t) = Cp m(t) + Cf m(t) + Cυ m(t) Ct(t) = Cf t (t) (2.31) Por fim, decompondo as forças normal e tangencial, o coeficiente de sustentação é dado: Cl(t) = Cn(t)cos[α(t)]− Ct(t)sen[α(t)], (2.32) e o coeficiente do momento de arfagem no eixo elástico é tal que: Cmea = Cn(t)(xea − xac) + Cm(t) (2.33) 39 3 METODOLOGIA 3.1 MÉTODO COMPUTACIONAL Para a solução do problema no tempo, uma vez que a modelagem aerodinâmica proposta por Bedoes-Leishman é dada por equações diferenciais ordinárias de primeira ordem, em espaços de estados, necessita-se de um método computacional capaz de calcular a resposta no estado seguinte tendo como base o estado atual. Para isso, um método numérico adequado com baixo erro é o Runge-Kutta de quarta ordem. A solução no passo de tempo seguinte é dada pelo estado atual somado o produto do passo de tempo e uma inclinação estimada ponderada (CENGEL; Palm III, 2014). Matematicamente: yn+1 = yn + ts(φ1 + 2φ2 + 2φ3 + φ4) 6 (3.1) onde ts é o passo de tempo, e φi são estimativas da inclinação, dadas por: φ1 = f(xn, yn) φ2 = f(xn + ts 2 , yn + ts 2 φ1) φ3 = f(xn + ts 2 , yn + ts 2 φ2) φ4 = f(xn + ts, yn + tsφ3) (3.2) Destaca-se que a formulação acima apresentada também é válida vetorialmente. A resposta estrutural também é calculada com base nesse método, tendo como entradas os dados da saída aerodinâmica da modelagem de BL. Essa resposta, por sua vez, é utilizada para o cálculo do estado aerodinâmico seguinte, repetindo-se o processo até o final do laço de tempo de simulação desejado. Como solução do problema, são então apresentadas as respostas de deslocamento (h), ângulos de incidência (θ) e de ataque (α) do perfil, ao longo do tempo. Os cálculos foram realizados através do software MATLAB, que tem ferramentas apropriadas para isso. 3.2 FLUXOGRAMA DE CÁLCULO A Figura 5 ilustra o fluxograma seguido na simulação para resolução das equações diferenciais ordinárias apresentadas. 40 Figura 5 – Fluxograma Resolução Simulação fonte: Elaborado pelo autor 41 4 PARÂMETROS Os parâmetros estruturais do aerofólio e aparato estão dispostos na tabela 1 a seguir. Tabela 1 – Parâmetros Estruturais do Modelo Parâmetros Estruturais Parâmetro Variável Valor Corda (m) c 0,25 Semicorda (m) b 0,125 Massa do aerofólio (kg) m 1,5 Massa total (kg) mt 3,6733 Frequência natural de translação (rad/s) ωh 21,77 Frequência natural de rotação (rad/s) ωθ 24,85 Distância adimensional entre eixo elástico e CG xθ 0,66 Posição relativa do eixo elástico xea 0,25 Posição relativa do centro aerodinâmico xac 0,25 Raio de giração rθ 0,7303 Elementos da matriz de amortecimento de Rayleigh di,j [ 5, 49 10, 03 10, 03 32, 48 ] fonte: (SANTOS; PEREIRA; MARQUES, 2017) A tabela 2 abaixo apresenta os parâmetros relacionados à aerodinâmica. Tabela 2 – Parâmetros Aerodinâmicos Parâmetros Aerodinâmicos Parâmetro Variável Valor Densidade do ar (kg/m3) ρ 1,225 Velocidade crítica (m/s) Vcr 14,2 Velocidade do som (m/s) Vs 343 fonte: (SANTOS; PEREIRA; MARQUES, 2017) - Adaptado. Por fim, na tabela 3 a seguir são mostradas as constantes empíricas relacionadas à parte não linear do modelo aerodinâmico de BL. 42 Tabela 3 – Constantes Experimentais: Modelagem Aerodinâmica Não-Linear de Beddoes-Leishman Parâmetros Aerodinâmica Não-Linear BL Parâmetro Variável Valor Ângulo onde f(α) = 0, 7 (◦) α10 15,25 Constante (◦) s1 3,0 Constante (◦) s2 2,3 Constante K0 0,0025 Constante K1 -0,135 Constante K2 0,04 Fator de ponderação força tangencial η 0,965 Constante de tempo 1º estado não linear Tp 1,7 Constante de tempo - início Tf Tf0 3,0 Constante de tempo - início Tυ Tυ0 6,0 Janela de tempo viagem do vórtice BA-BF Tυl 7,0 Cn na iminência do estol Cn1 1,45 Máxima variação de α1 ( ◦) δα1 2,1 Coeficiente de momento de arfagem em sustentação zero Cm0 0 fonte: (SANTOS; PEREIRA; MARQUES, 2017) Em relação à simulação numérica, é importante ressaltar que o passo de tempo deve ser escolhido com atenção, pois valores grandes podem não levar à convergência da resposta, enquanto que steps muito pequenos podem demorar a finalizar as iterações. Assim, o passo utilizado nesse trabalho é de 0, 75 · 10−4 s. Para inicio dos cálculos, o deslocamento inicial dado é h = 0, 01 m. 43 5 RESULTADOS 5.1 VERIFICAÇÃO PRELIMINAR As figuras a seguir apresentam as respostas nos graus de liberdade modelados em (2.6), desloca- mento h e ângulo θ, considerando uma entrada impulso de h0 = 10 mm. As velocidades simuladas foram menores, iguais e maiores a Vcr. Na Figura 6, observa-se que para velocidades do escoamento menores que a velocidade crítica do flutter, as oscilações tendem a diminuir. Isso se deve ao fato dos esforços estruturais serem maiores que os aerodinâmicos. Figura 6 – Respostas Estados Estruturais, U < Vcrit (a) Resposta h(t) (b) Resposta θ (t) fonte: Elaborado pelo autor A Figura 7 apresenta a condição em que a velocidade do escoamento se iguala à crítica, o que reflete o equilíbrio dos esforços estruturais e aerodinâmicos. Isso resulta na manutenção da oscilação no decorrer do tempo, com a mesma amplitude inicial. Figura 7 – Respostas Estados Estruturais, U = Vcrit (a) Resposta h(t) (b) Resposta θ (t) fonte: Elaborado pelo autor 44 Por sua vez, com a Figura 8 é possível verificar como o modelo se comparta sob um escoamento com velocidade superior à Vcr. Tal como esperado, as oscilações crescem com o passar do tempo. Ressalta- se a criticidade de tal fenômeno, dadas as amplitudes do ângulo torcional, condição estruturalmente instável e perigosa. Figura 8 – Respostas Estados Estruturais, U > Vcrit (a) Resposta h(t) (b) Resposta θ (t) fonte: Elaborado pelo autor Ainda, na Figura 9 abaixo está ilustrado um dos efeitos capturados pela modelagem aerodinâmica não linear. A amplitude das oscilações, apesar de crescer, têm um limite e se estabilizam. Isso é causado pelo estol dinâmico, contabilizando as perdas nos esforços aerodinâmicos em razão do descolamento do escoamento e desprendimento de vórtices, aproximando assim o resultado de uma resposta realista. Figura 9 – Resposta h - Efeito Estol Dinâmico fonte: Elaborado pelo autor Por fim, a Figura 10 traz uma comparação entre as modelagens de aerodinâmica, linear e não linear, descrito pelos últimos quatros estados propostos por BL (x9, ..., x12). 45 Figura 10 – Resposta θ - Aerodinâmica Linear × Não Linear fonte: Elaborado pelo autor A figura acima ilustra o comportamento esperado com a utilização da abordagem aerodinâmica não linear, em que as oscilações, por mais que sejam crescentes, se estabilizam em um limite máximo. Diante do acima exposto, é possível verificar a funcionalidade do método, modelagens, em consonância com o apresentado na literatura, apesar de uma pequena diferença na velocidade crítica encontrada por Santos et al. (SANTOS; PEREIRA; MARQUES, 2017). 5.2 COMPORTAMENTOS AERODINÂMICOS 5.2.1 Bifurcação α vs U A Figura 11 ilustra o comportamento da amplitude do ângulo de ataque aerodinâmico (α) em função da velocidade do escoamento, para os casos lineares e não lineares abordados na modelagem conforme BL. Figura 11 – Amplitude Ângulo de Ataque Aerodinâmico × Velocidade do Escoamento fonte: Elaborado pelo autor 46 5.2.2 Ângulo de Ataque (α) vs Tempo Na Figura 12, é apresentada a resposta no tempo do ângulo de ataque, para diferentes velocidades de escoamento, em flutter. A frequência de oscilação aumenta com a velocidade do escoamento, o que está de acordo com a literatura. Figura 12 – Ângulo de Ataque × Tempo fonte: Elaborado pelo autor. 5.2.3 Ângulo de Ataque Aerodinâmico (α) e Ângulo de Pitch Estrutural (θ) vs Tempo Da Figura 13 a seguir, é possível verificar a relação entre os ângulos de ataque aerodinâmicos e de pitch estrutural, e a maneira como parcela da velocidade ḣ contribui na composição daquele. Figura 13 – Ângulo de Ataque Aerodinâmico e Ângulo de Pitch Estrutural × Tempo fonte: Elaborado pelo autor. 47 5.2.4 Coeficientes de Força Normal e Momento de Arfagem 5.2.4.1 Resposta no Tempo Figura 14 – Coeficientes de Força Normal e Momento de Arfagem, U = Vcrit (a) Cn (b) Cmea fonte: Elaborado pelo autor Sob condição de flutter na velocidade crítica, observa-se que a componente potencial do coeficiente de força normal, dada pela parcela da aerodinâmica linear, é em sua totalidade equivalente ao coeficiente total com as demais contribuições. Nessa situação ainda é possível verificar o coeficiente relativo à separação do escoamento, levemente defasado, que é prontamente compensado. Também ressalta-se a não influência do desprendimento de vórtices, ainda, devido ao baixo módulo de α. Em relação ao coeficiente de momento, os efeitos não lineares são praticamente desprezíveis, sendo o coeficiente total quase todo correspondente à parcela potencial. O ângulo de ataque aerodinâmico (α) está ilustrado na Figura 15 simplesmente para efeitos comparativos e fornecer uma compreensão física mais clara. 48 Figura 15 – Ângulo de Ataque Aerodinâmico × Tempo, U = 14, 48m/s fonte: Elaborado pelo autor As figuras a seguir trazem a resposta da simulação para velocidades do escoamento superiores à velocidade crítica, em que as oscilações atingem amplitudes maiores. Os efeitos em cada componente são discutidos posteriormente. Figura 16 – Coeficientes de Força Normal e Momento de Arfagem, U = 17m/s (a) Cn (b) Cmea fonte: Elaborado pelo autor 49 A Figura 17 apresenta o comportamento em regime permanente dos mesmos coeficientes, em detalhe, onde é possível identificar que o descolamento do escoamento da superfície do aerofólio ocasiona uma redução no coeficiente de força normal e amplificação no de momento de arfagem no eixo elástico. Figura 17 – Coeficientes de Força Normal e Momento de Arfagem - Regime Permanente, U = 17m/s (a) Cn (b) Cmea fonte: Elaborado pelo autor Na Figura 18, o ângulo de ataque aerodinâmico é novamente apresentado para efeito de com- preensão física. A diferença, entretanto, é que nessa velocidade simulada, o referido ângulo atinge magnitudes mais elevadas, superando o ângulo de estol estático, o que justifica os comportamentos observados em detrimento dos desprendimentos de vórtices. A Figura 19, traz em detalhes a contribuição do desprendimento de vórtices nos coeficientes de força normal e momento de arfagem. Apesar de pouco significativos, não são nulos. 50 Figura 18 – Ângulo de Ataque Aerodinâmico × Tempo, U = 17m/s fonte: Elaborado pelo autor. Figura 19 – Coeficientes de Força Normal e Momento de Arfagem, Induzidos por Desprendimento de Vórtices, U = 17m/s fonte: Elaborado pelo autor. 5.2.4.2 Envelope de α Na Figura 20 e Figura 21 estão descritos os comportamentos dos coeficientes de força normal e momento de arfagem em regime permanente, totais e referente às suas parcelas, respectivamente. Os coeficientes totais estão em plena consonância com a literatura, sendo a resposta de ambos como apresentado na Figura 3. O comportamento do coeficiente de força normal é coerente com a realidade ao passo que cresce com o aumento do ângulo de ataque, até determinado ponto, em que incremento em α reduz o coeficiente. Isso condiz com uma condição de estol severo. Já o coeficiente de momento tem um aspecto predominantemente elíptico, exceto para as regiões próximas ao máximo e mínimo ângulo de ataque, que como é possível observar nas suas componentes constituintes, deve-se às parcelas de separação e desprendimento de vórtices. 51 Figura 20 – Coeficientes Totais de Força Normal e Momento de Arfagem × α, Regime Permanente, U = 17m/s (a) Cn (b) Cm fonte: Elaborado pelo autor. A tabela 4 abaixo resume as principais estruturas do escoamento e descrição das forças e momentos nos pontos identificados na Figura 20. Tabela 4 – Descrição Escoamento e Efeitos Carregamento Aerodinâmico Ponto Estrutura do Escoamento Força e Momento 1 Reversão do escoamento dentro da camada limite, formação do vórtice Ultrapassa sustentação máxima, extrapolação da faixa linear 2 Desprendimento do vórtice e deslocamento sobre o aerofólio Divergência do momento de arfagem, vórtice presente 3 Vórtice passa pelo BF, desenvolvimento completo do estol Sustentação máxima, rápido decaimento, Cmea mínimo e máximo arrasto 4 Recolamento do escoamento Reajuste à faixa linear fonte: Elaborado pelo autor. O ponto 1 no coeficiente de força normal pode ser explicado com o auxílio da Figura 16 (a), Figura 17 (a), Figura 18 e Figura 19. Em aproximadamente 2,2 s ocorre a formação do vórtice, em que o ângulo de ataque é quase 20°. É possível observar o fenômeno também em 3,2 s, entre outros momentos. Em seguida, ocorre o desprendimento dos vórtices, marcados pelo ponto 2 , que provocam a divergência do momento de arfagem. No ponto 3 têm-se como característica a máxima sustentação, seguida de rápido decaimento, comportamento referente ao estol completamente desenvolvido. Além disso, o coeficiente de momento de arfagem é mínimo nesse instante, enquanto que o arrasto é máximo. O arrasto pode ser visualizado como sendo a componente de força tangente multiplicada pelo cosseno do ângulo de ataque. A referida força está apresentada na Figura 25. 52 Por fim, o recolamento do escoamento se dá no ponto 4 . Na Figura 17 (a) isso pode ser verificado em aproximadamente 3,1 s e 9°. Esse movimento pode ser entendido, ainda, como a aproximação do coeficiente Cf n no coeficiente potencial Cp n, proveniente da aerodinâmica linear. Figura 21 – Parcelas dos Coeficientes de Força Normal e Momento de Arfagem × α, Regime Perma- nente, U = 17m/s (a) Cn (b) Cm fonte: Elaborado pelo autor. Em relação às parcelas constituintes do coeficiente de força normal, a potencial segue o proposto pela aerodinâmica linear, uma vez que sua modelagem assim é feita. Por sua vez, a parcela referente à separação tem como causa a redução do coeficiente total. Além disso, o impacto dos vórtices, apesar de baixa magnitude, não se anula devido ao fato de novos desprendimentos ocorrerem antes que os anteriores tenham tido seus efeitos dissipados após o BF, conforme proposto pelo tempo adimensional de 2Tυl. Ambas a Figura 22 e Figura 23 trazem o desenvolvimento dos coeficientes do regime transiente ao regime permanente. 53 Figura 22 – Cn × α: Regime Transiente + Permanente, U = 17m/s fonte: Elaborado pelo autor. Figura 23 – Cm × α: Regime Transiente + Permanente, U = 17m/s fonte: Elaborado pelo autor. 54 5.2.5 Coeficiente de Sustentação 5.2.5.1 Resposta no Tempo Figura 24 – Coeficientes de Sustentação, U = 17m/s fonte: Elaborado pelo autor. 5.2.5.2 Envelope α Figura 25 – Coeficiente de Sustentação × α, U = 17m/s fonte: Elaborado pelo autor. Na figura acima está apresentado o coeficiente de força normal e tangencial, bem como o coeficiente de sustentação, como resultado da combinação desses. Devido ao movimento oscilatório em arfagem, o coeficiente de força normal tem sua orientação de forma alternada, entretanto, isso não se observa para a força tangente. Ainda, esta tem seu pico nos ângulos de ataque máximos, gerando, consequentemente, uma maior força de arrasto nessa condição. 55 5.2.6 Estados Aerodinâmicos Figura 26 – Estados Aerodinâmicos de Beddoes-Leishman, U = 17m/s fonte: Elaborado pelo autor. 5.2.7 Contador τv e Estados Aerodinâmicos Não Lineares O parâmetro de simulação τv é um contador que indica o momento que os efeitos não lineares estão ocorrendo, ou seja, o estol dinâmico. Ele passa a ser contabilizado a partir do momento que o módulo do estado x9, que representa C ′ n ultrapassa o valor de Cn na iminência do estol sob condições estáticas (Cn1). Na Figura 28 tem-se ilustrado o comportamento do referido estado aerodinâmico, que está estri- tamente relacionado com o ângulo de ataque aerodinâmico. Dessa forma, é possível identificar os momentos em que as não linearidades aerodinâmicas estão sendo contabilizadas. 56 Figura 27 – Contador de Tempo Efeitos Não Lineares, U = 21m/s fonte: Elaborado pelo autor. Figura 28 – x9 e α, U = 21m/s fonte: Elaborado pelo autor. Ainda, na Figura 29 abaixo estão apresentados os demais estados aerodinâmicos não lineares, x10, x11, relativos à separação do escoamento, e x12, referente ao desprendimento de vórtices no bordo de ataque. É possível identificar, então, que conforme α cresce, o ponto de separação se desloca em direção ao BA, até que ocorra o desprendimento do vórtice. Esse vórtice, por sua vez, viaja pela superfície do aerofólio até sair pelo BF, e seu efeito nos coeficientes aerodinâmicos é contabilizado por Cυ n e Cυ m, até o instante de 2Tυl. 57 Figura 29 – Estados Aerodinâmicos Não Lineares BL, U = 21m/s fonte: Elaborado pelo autor. 58 5.2.8 Coeficientes de Força Normal e Momento de Arfagem 5.2.8.1 Resposta no Tempo Figura 30 – Coeficientes de Força Normal e Momento de Arfagem, U = 21m/s (a) Cn (b) Cmea fonte: Elaborado pelo autor As análises das figuras 30 e 31, com a velocidade do escoamento a 21 m/s são semelhantes às análises feitas para o escoamento a 17 m/s. 59 5.2.8.2 Envelope de α Figura 31 – Coeficientes Totais de Força Normal e Momento de Arfagem × α, Regime Permanente (a) Cn (b) Cm fonte: Elaborado pelo autor. 5.2.9 Fase Não Linear A Figura 32 a seguir apresenta o comportamento do ângulo de ataque quanto ao início e término da fase não linear em função da velocidade do escoamento. A curva em azul representa o início da fase de desprendimento de vórtices, enquanto que em laranja, o início da fase de recolamento do escoamento. Figura 32 – Início e Fim Fase Não Linear fonte: Elaborado pelo autor É possível identificar que o início da fase de desprendimento de vórtices tem um rápido crescimento inicial, seguido de estabilização e decrescendo de maneira suave posteriormente. Por sua vez, o ângulo no qual a fase de recolamento tem início cai conforme maior é a velocidade do escoamento. Ambos os comportamentos fazem sentido fisicamente. 60 5.2.10 Ponto de separação Para simular o efeito do ponto de separação, dado pela função f (eq. 2.19), foi simulado uma separação forçada do escoamento na metade da corda. Esse tipo de situação fisicamente pode ser obtida com dispositivos sobre o aerofólio. Figura 33 – Estados Aerodinâmicos Não Lineares BL, U = 21m/s, f = 0.5 fonte: Elaborado pelo autor É possível notar que a influência dos vórtices têm um inicio antecipado, como era de se esperar, comparando com a função por partes proposta por BL. Figura 34 – Coeficientes de Força Normal e Momento de Arfagem, U = 21m/s, f = 0.5 (a) Cn (b) Cmea fonte: Elaborado pelo autor Em relação ao coeficiente de força normal, não há uma alteração significativa em seu módulo, porém observa-se que a parcela referente ao desprendimento de vórtices está mais suave. A componente devido a separação não tem o comportamento característico de recolamento, como discutido anteriormente na Tabela 4, também sendo visível na Figura 30 (a), conforme esperado. 61 O coeficiente de momento de arfagem, por sua vez, apresenta módulo consideravelmente inferior com a separação na metade da corda, com o comportamento da parcela referente à separação em oscilações mais regulares. A Figura 35 apresenta a evolução do ângulo de ataque aerodinâmico em função do ponto de separação. Com ela é possível identificar que a separação na metade da corda causa um efeito de redução na velocidade com que as oscilações crescem. O pico atingido, entretanto, não sofre alteração. Figura 35 – Ângulo de Ataque Aerodinâmico × ponto de separação, U = 21m/s fonte: Elaborado pelo autor Figura 36 – Coeficiente de Sustentação × ponto de separação, U = 21m/s fonte: Elaborado pelo autor Sobre o coeficiente de sustentação apresentado pela Figura 36, o comportamento é semelhante ao analisado para o ângulo de ataque. Contudo, sua amplitude tem uma redução, o que condiz com 62 o esperado, uma vez que a distribuição de pressão sobre o aerofólio está alterada com a separação ocorrendo na metade da corda. 63 6 CONCLUSÃO Diante de todo o apresentado, o presente trabalho buscou analisar um perfil bidimensional em condição de flutter e post-flutter, usando a modelagem desenvolvida por Beddoes-Leishman. Modela- gem esta com diversas características favoráveis para uma análise aeroelástica, tais como utilização de aerodinâmica não estacionária e não linear, particularmente na questão do estol dinâmico. Além disso, a formulação em espaços de estados é atraente para resolver simulações no domínio do tempo uma vez que são facilmente integradas. Vale destacar a acurácia dessa modelagem, especialmente em grandes deslocamentos e elevados ângulos de ataque. Especificamente em relação aos coeficientes de força normal e momento de arfagem, foi possível esclarecer a contribuição individual e o comportamento de cada parcela na composição dos respetivos coeficientes totais, sobretudo das componentes não lineares, particularidade da modelagem de BL. Além disso, validou-se a estrutura do escoamento em regime permanente sob flutter com o previsto na literatura, identificando a influência dos coeficientes de força normal e momento de arfagem nela. Por fim, ainda foi feita uma análise dos estados aerodinâmicos não lineares e da fase em que as não linearidades estavam presentes, além de uma breve discussão sobre o efeito do ponto de separação do escoamento. Ressalva-se, apenas, que a velocidade crítica de flutter encontrada pela simulação diverge ligei- ramente da obtida tanto experimentalmente quanto da simulação na literatura de referência. Isso deve-se provavelmente a algum parâmetro de entrada divergente, uma vez que a velocidade crítica é característica do sistema linear. Contudo, os resultados obtidos não sofreram impactos significativos em decorrência disso. 64 b 65 7 CONSIDERAÇÕES FINAIS Sem mais, o trabalho em questão atingiu os objetivos propostos inicialmente, demonstrando a efetividade da modelagem. Destaca-se que os parâmetros estruturais, bem como os relacionados à aerodinâmica não linear, são características particulares e importantes a serem levados em conta. Como sugestão para trabalhos futuros, sugere-se repetir os estudos aqui propostos para outro perfil, e comparar os resultados obtidos. Também pode-se identificar a influência dos parâmetros aerodinâmicos não lineares, que são dependentes do número de Mach do escoamento. Ainda, obter novos dados experimentais no túnel de vento da Faculdade de Engenharia de São João da Boa Vista e realizar a simulação para validação é uma pesquisa relevante a se considerar. 66 b 67 REFERÊNCIAS ABBOTT, I. H.; DOENHOFF, A. E. von. Theory of Wing Sections, Including a Summary of Airfoil Data. New York: Dover Publications, 1959. BISPLINGHOF, R. L.; ASHLEY, H.; HALFMAN, R. L. Aeroelasticity. New York: Dover Publications, 1955. CHANTHARASENAWONG, C. Nonlinear Aeroelastic Behaviour of Aerofoils Under Dynamic Stall. Tese (Doutorado) — University of London, South Kensington Campus, 2007. Department of Aeronautics - Imperial College London. COLLAR, A. R. The First Fifty Years of Aeroelasticity. 1978. CENGEL, Y. A.; Palm III, W. J. Equações Diferenciais. Porto Alegre: AMGH Editora LTDA, 2014. DOWELL, E. H. A Modern Course in Aeroelasticity. Fifth Revised and Enlarged Edition. New York: Springer, 2015. ERTURK A.; INMAN, D. J. Piezoelectric energy harvesting. United Kingdom: John Wiley Sons, Ltd., 2011. HODGES, D. H.; PIERCE, G. A. Introduction to Structural Dynamics and Aeroelasticity. Cambridge, UK: Cambridge University Press, 2002. LEISHMAN, J.; CROUSE, G. A state-space model of unsteady aerodynamics in a compressible flow for flutter analyses. American Institute of Aeronautics and Astronautics, 1989. LEISHMAN, J. G.; BEDDOES, T. S. A semi-empirical model for dynamic stall. Journal of the American Helicopter Society, v. 34, p. 3–17, 1989. MULLENERS K., R. M. Dynamic stall development. Exp Fluids, v. 54, p. 1469–1477, 2013. NELSON, R. C. Flight Stability and Automatic Control. Boston, MA - USA: McGraw Hill, 2ed., 1997. SANTOS, C. R. dos. Modeling the nonlinear electro-aeroelastic behavior of energy harvesters from stall-induced oscillations. Tese (Doutorado) — School of Engineering of the University of São Paulo, São Carlos, 2021. SANTOS, C. R. dos; PEREIRA, D. A.; MARQUES, F. D. On limit cycle oscillations of typical aeroelastic section with different preset angles of incidence at low airspeeds. Journal of Fluids and Structures, v. 74, p. 19–34, 2017. SOUSA, V. C. de. Effects of superelastic shape memory springs on the aeroelastic behavior of a typical airfoil section: passive vibration attenuation and energy harvesting applications. Tese (Doutorado) — School of Engineering of the University of São Paulo, São Carlos, 2016. 68 b 69 APÊNDICE A – RELAÇÕES UTILIZADAS NA ADIMENSIONALIZAÇÃO DA MODELAGEM ESTRUTURAL rθ = √ Iθ m b2 ωh = √ kh m ωθ = √ kθ Iθ 70 b 71 ANEXO A – DADOS EXPERIMENTAIS NA LITERATURA NACA 0012 Figura 37 – Coeficiente de Sustentação e Momento NACA 0012 fonte: (ABBOTT; DOENHOFF, 1959) a7cb6b9607ea2d6babeffc233afc0f9ceae706a7e5129f3c6068d36f2e7b0169.pdf Folha de rosto 1503bab167e728fc34a13be4d92703592b22d61b38bf9d20bc9748a70882c46e.pdf ACFrOgAoGhXuz_vpXlZ0__hAFRquLkwSSxeOZrM...GidBMGIbPK5hKZtuL6Yk0La-1NUVhnCToYcTyE a7cb6b9607ea2d6babeffc233afc0f9ceae706a7e5129f3c6068d36f2e7b0169.pdf Dedicatória Agradecimentos Epígrafe Resumo Abstract Lista de abreviaturas e siglas Lista de símbolos Introdução OBJETIVOS ESTRUTURA DO TRABALHO Modelo Matemático EQUACIONAMENTO ESTRUTURAL EQUACIONAMENTO AERODINÂMICO - MODELAGEM DE BEDDOES-LEISHMAN Equacionamento Aerodinâmica Linear Equacionamento Aerodinâmica Não-linear METODOLOGIA MÉTODO COMPUTACIONAL FLUXOGRAMA DE CÁLCULO PARÂMETROS RESULTADOS VERIFICAÇÃO PRELIMINAR COMPORTAMENTOS AERODINÂMICOS Bifurcação vs U Ângulo de Ataque () vs Tempo Ângulo de Ataque Aerodinâmico () e Ângulo de Pitch Estrutural () vs Tempo Coeficientes de Força Normal e Momento de Arfagem Resposta no Tempo Envelope de Coeficiente de Sustentação Resposta no Tempo Envelope Estados Aerodinâmicos Contador v e Estados Aerodinâmicos Não Lineares Coeficientes de Força Normal e Momento de Arfagem Resposta no Tempo Envelope de Fase Não Linear Ponto de separação CONCLUSÃO CONSIDERAÇÕES FINAIS Referências Relações Utilizadas na Adimensionalização da Modelagem Estrutural Dados Experimentais na Literatura NACA 0012