Método de Painéis de Vórtice 2D - Aerodinâmica - Cálculo Numérico MATLAB

Trabalho N1 PME 2557 Aerodinmica Professor Julio Meneghini 1 Semestre / 2013 Mtodo de Painis de Vrtice

Programa Computacional de Clculo Para Qualquer Perfil

Trabalho realizado por: Andr Valera N 8197662 So Paulo, Abril de 2013

Mtodo de Painis de Vrtice 2

ndice 1. Introduo .............................................................................................................. 3 2. Mtodo de Painis de Vrtice ........................................................................ 4 Formulao Matemtica .............................................................................................................. 4 3. Perfil NACA 4-series Gerao das Coordenadas do Perfil ............... 8 4. Para o Utilizador Programa em MATLAB ............................................ 10 5. Avaliao do Mtodo Anlise de Perfis ................................................ 11 Perfil de Cilindro Infinito ......................................................................................................... 11 Perfil Triangular .......................................................................................................................... 14 Perfil NACA simtrico NACA0012 ....................................................................................... 17 Perfil NACA no simtrico NACA4212 ............................................................................... 21 Comparao entre NACA0012 e NACA4212 ................................................................... 24 6. Concluso .............................................................................................................. 25 7. Bibliografia ........................................................................................................... 25 Livros ................................................................................................................................................ 25 Stios da Internet ......................................................................................................................... 25 8. Anexo Cdigo do Programa em MATLAB ............................................ 26


1. Introduo Este trabalho tem como objetivo a compreenso e aplicao do mtodo numrico de painis de vrtice a duas dimenses (2D). Este mtodo permite calcular numericamente o escoamento em torno de qualquer geometria 2D, sendo naturalmente interessante a sua aplicao a sees de asa, para estudo de propriedades como a sustentao de uma aeronave. Este mtodo, sendo um mtodo de elementos de contorno, discretiza a superfcie do corpo em avaliao, ao invs do fludo que escoa em seu redor, como seria feito no caso do mtodo de elementos finitos. Neste estudo numrico, a compressibilidade e a viscosidade so desprezadas. Uma descrio mais pormenorizada do mtodo ser apresentada no Captulo 2. Para verificar a aplicabilidade do mtodo, o principal elemento do trabalho a elaborao de um programa que o aplique e que apresente resultados pertinentes e teis, principalmente no que diz respeito a propriedades aerodinmicas das sees avaliadas. Assim, utilizou-se o MATLAB como ambiente de programao, sendo o cdigo do programa apresentado em anexo e em formato digital no CD para possvel utilizao. Este programa permite ao utilizador avaliar a seo de um cilindro, a seo de um perfil triangular com qualquer declive e qualquer seo simtrica ou no simtrica da srie de perfis aeronuticos NACA #### (4-series). Com essa avaliao e resultados numricos pretende-se validar o mtodo, comparando os dados obtidos com dados experimentais e tericos.


2. Mtodo de Painis de Vrtice O mtodo de painis de vrtice um mtodo numrico que permite avaliar o comportamento de um qualquer corpo bidimensional quando exposto a um dado escoamento. Sendo baseado na discretizao da seo a avaliar , resulta na substituio da superfcie do corpo por m painis linha que interligam m+1 pontos pertencentes superfcie do corpo formando um sistema de painis de singularidade, tais como fontes ou sorvedouros, dipolos e vrtices. Qualquer singularidade pode ser utilizada, tendo sido escolhida a singularidade de vrtice. O programa, na sua sntese, calcula a intensidade (constante) do vrtice a partir das condies de contorno do problema. O tipo de painel utilizado no programa o mais simples, o painel linear, para mais fcil implementao e programao. Outro painis, como os quadrticos, so verses otimizadas do mtodo, mais complexas, que procuram uma melhor aproximao do sistema de painis geometria da seo. Sobrepondo um escoamento uniforme aos vrtices distribudos pelos painis, tem-se uma soluo da equao de Laplace. Tendo sido admitido um escoamento de fludo no-viscoso, a soluo deve obedecer condio de impenetrabilidade da fronteira. Para isso aplicou-se a condio de contorno de Neumann. Finalmente, necessrio satisfazer a condio de Kutta, que a imposio da unicidade da circulao, resultando na igualdade das velocidades superior e inferior no bordo de fuga da seo, que equivalente a igualar as presses superior e inferior. A soluo de tudo o descrito acima, ou seja, encontrando a distribuio de vrtices que satisfaz as condies de contorno e a condio de Kutta, vai permitir o clculo imediato de importantes coeficientes aerodinmicos caractersticos da seo avaliada coeficiente de presso, que por sua vez permite obter os coeficientes de arrasto e de sustentao. Este ltimo o alvo de todo o mtodo e de todo o programa, pois avalia a capacidade de sustentao oferecida pelo aeroflio, fulcral no projeto de aeronaves. Segue-se a formulao matemtica foi implementada em cdigo no programa. Formulao Matemtica Considerando um aeroflio com m painis lineares de vrtice com intensidade constante j, o potencial induzido em qualquer ponto (x,y) no campo de escoamento pelo painel j :

! , = !2 tan!! ! ! ! ! !!!!


Onde o painel tipo j de comprimento sj est representado abaixo:

Incluindo o escoamento uniforme de velocidade U e somando as contribuies de todos os painis, o potencial total no ponto (x,y) : ! , = cos + sin !2 tan!! ! ! ! ! !!!!!!!! Tendo sido adotados painis lineares: ! ! = !!! !! + ! = cos!! + ! ! ! = !!! !! + ! = sin!! +! Ficando o potencial:

! , = cos + sin !2 tan!! sin!! ! cos!! ! !!!!!!!! Aplicando a condio de contorno de Neumann: = 0 Esta condio aplicada a cada ponto de controlo de modo que: ! ! ,!! = ! cos + ! sin !2 ! tan!! sin!! ! cos!! ! !!!! !!!!!!!! Resultando: !2 ! tan!! sin!! ! cos!! ! !!!! ! =!!!!!!! sin(! )


Desenvolvendo a derivada no integral, conclui-se que a condio fronteira se resume a: !2 !"!"2 + !"!" =!!!! sin(! ) Onde: !" = (! !) cos! (! !) sin! !" = (! !)! + (! !)! !" = sin ! ! , !" = cos ! ! , !" = ln 1+ !! + 2!"!!" !" = (! !) sin! (! !) cos! !" = tan!! !"!!"! + !" , !!!!2 + !!!! = !" = ! ! sin ! 2! + (! !) cos ! 2! !" = ! ! cos ! 2! + (! !) sin ! 2! Assim, em notao matricial, para obter : ! = sin(! ) E a velocidade tangencial ao painel i: !" = cos ! + ! Com: ! = 2 + 2! ( + )! ! = 2 + 2! + ( )! Impondo ainda que: !"" = 0, !"! = !!!!2 !!!! !! + !!!!2! !!(!!!! + !!!!)! !"(!!!) = !" + !"!"2! !"(!"!" + !"!")!


!"" = , !"! = !!!!2 !!!! !! + !!!!2! !!(!!!! + !!!!)! !"(!!!) = !" + !"!"2! !"(!"!" + !"!")! E, pela condio de Kutta: ! = ! 1 ! = !(1) Resultando: ! !!! ! = 1, ! !!! (!!!) = 1, sin(!!! ) = 0 Todos os clculos acima resolvem a matriz do sistema montado pela equao de Laplace, sendo agora necessrio aplicar as suas solues s frmulas que permitem obter os coeficientes aerodinmicos pretendidos: !" = 1 !" !, ! = !"! sin! , ! = !"! cos! !"#$ = ! cos ! sin !"#$ = ! sin + ! cos De referir que na formulao acima consta a velocidade do escoamento uniforme U, sendo que no programa se trabalha com valores adimensionais, ou seja, a velocidade tangencial calculada normalizada pela velocidade U. Est assim concluda a formulao para obter os coeficientes aerodinmicos de sustentao e arrasto, que foi aplicada no cdigo do programa, com pequenas manipulaes matriciais. No que diz respeito ao nmero de painis e, portanto, nmero de pontos de controle a avaliar nos perfis, foi escolhido um nmero adequado a cada tipo de perfil suficiente para aproximar a geometria deste. Em todos os perfis utilizado um nmero maior de painis nos bordos de ataque e fuga, de modo a ter uma maior densidade de controle (nmero de pontos de controle/percentagem da corda) nessas zonas sensveis de anlise aerodinmica, resultando numa anlise mais fina. Particularmente nos perfis NACA, utilizou-se 200 painis, ou seja, densidade de controle de 200/1 = 200. A densidade de pontos de controle no bordo de ataque de 40/0,03 = 1333,33 e no bordo de fuga 60/0,1 = 600.


3. Perfil NACA 4-series Gerao das Coordenadas do Perfil O perfil NACA 4-series vai ser o principal perfil estudado pelo programa, sendo por isso importante apresentar a formulao matemtica de obteno da geometria e coordenadas deste tipo de perfil. Essa formulao faz parte do cdigo do programa para definir os painis e pontos de controle dos perfis NACA. Segue-se um resumo e apresentao dessa formulao. Esta metodologia comeou nos anos 30 com a publicao de um relatrio da NACA, em que os seus autores relacionaram os aeroflios de melhor desempenho com as duas principais variveis que afetava esses perfis, o declive da linha media do perfil e a distribuio de espessura acima e abaixo desta linha. Apresentaram um conjunto de equaes que incorporavam estas duas variveis, que podiam ser usadas para gerar uma famlia de perfis. A famlia de perfis que se ir estudar a famlia 4-series. O esquema geomtrico de um perfil NACA o seguinte.

Fig. 1 Construo geomtrica de um perfil NACA A famlia 4-series um perfil NACA composto por quatro nmeros, NACA####. O primeiro dgito especifica o arqueamento mximo (m) em percentagem da corda, o segundo indica a posio do arqueamento mximo (p) em dezenas da corda e os dois ltimos so relativos mxima espessura (t) em percentagem da corda. Utilizando os valores de m, p e t possvel obter as coordenadas usando as frmulas seguintes. Para valores de x entre 0 e c, a equao da linha mdia de arqueamento (yc) : ! = ! 2 ! , 0 ! = 1 ! 1 2 + 2 ! ,


A distribuio de espessura, acima (+) e abaixo () da linha mdia de arqueamento, pode ser calculada introduzindo o valor de t para cada x: ! = 0,2 0,2969 0,1260 0,3516! + 0,2843! 0,1036! As coordenadas finais podem ser calculadas para a superfcie superior do perfil (xs, ys) e superfcie inferior (xi, yi) usando as seguintes relaes: ! = ! sin ! = + ! cos ! = + ! sin ! = ! sin Onde: = tan!! ! Assim, j podem ser geradas as coordenadas de todo o perfil de qualquer NACA da famlia 4-series. Esta formulao utilizada no programa para gerar o perfil e tambm os painis.


4. Para o Utilizador Programa em MATLAB O programa elaborado no ambiente MATLAB permite, seguindo o mtodo atrs descrito, calcular os coeficientes aerodinmicos para qualquer perfil, tendo sido escolhidas 3 geometrias diferentes para o utilizador analisar. O programa est organizado num ficheiro principal, mainNACA.m, de chamada das funes .m necessrias para desenvolver o clculo. O utilizador pode escolher entre avaliar um perfil de cilindro (1), um perfil triangular (2) ou um perfil de aeroflio NACA #### (4series).

Depois de selecionado o perfil, escolhido o ngulo de ataque e, no caso do perfil triangular, o declive, no caso do perfil NACA, o nmero NACA pretendido.

De seguida , o programa faz todo o clculo e devolve quatro grficos com as caractersticas do perfil: perfil e distribuio dos painis, pontos de controle e coeficiente de presso nos pontos de controle, variao do coeficiente de presso ao longo da corda e/ou variao do quadrado da velocidade adimensional com a corda.


Adicionalmente, exportado um ficheiro texto com os resultados relevantes da anlise. Este ficheiro guardado na pasta output_perfis, constando no nome do ficheiro o tipo de perfil e o ngulo de ataque. No ficheiro texto, o utilizador encontra os valores dos coeficientes de sustentao e arrasto e as coordenadas dos pontos de controle. Todo o cdigo do programa apresentado em anexo. 5. Avaliao do Mtodo Anlise de Perfis Neste captulo apresentam-se os resultados obtidos para os vrios perfis considerados. Ser interessante comparar os valores resultantes da anlise numrica com os valores tericos e experimentais. Os perfis analisados sero o perfil de cilindro infinito, perfil triangular e dois perfis NACA, um simtrico e um no simtrico. Para cada perfil mostra-se a sua forma geomtrica e os painis delimitados. Depois, para diversos ngulos de ataque, apresentam-se os pontos de controle, as intensidades do coeficiente de presso (Cp) em cada ponto de controle e consideraes baseadas nessas intensidades e nos grficos da distribuio deste coeficiente e/ou do quadrado da velocidade tangencial adimensionalizada. Perfil de Cilindro Infinito Este perfil corresponde na sua geometria a um crculo de dimetro unitrio, ou seja, a corda do perfil a unidade. Para anlise numrica consideraram-se 146 painis, com densidade de controle nos bordos de ataque e fuga igual a 42/0,02 = 2100.

Fig. 2 Perfil de cilindro infinito e painis delimitados painis a verde nos bordos de ataque e fuga


ngulo de ataque = 0

Fig. 3 Pontos de controle a amarelo e intensidade do Cp em cada ponto Como se pode ver no grfico da fig. 2, o Cp negativo na superfcie superior o vetor com a sua intensidade est direcionado para dentro do perfil e positivo na superfcie inferior vetor direcionada para fora do perfil. Junto dos bordos ele so aproximadamente simtricos. Desta anlise, imediatamente se comprova que o cilindro infinito, quando num escoamento uniforme com ngulo de ataque nulo, tem sustentao nula, ou seja, pois o somatrio das foras integrao do Cp ao longo de toda a superfcie igual a zero, valor que o programa devolve para o coeficiente de sustentao.

Fig. 4 Cp vs j (inclinao do painel j) com circulao comeando no bordo de fuga O grfico da fig. 3 mostra o Cp ao longo da corda em funo da inclinao de cada painel, onde as distribuies nas superfcies superior e inferior esto sobrepostas, uma vez que so iguais.


O grfico obtido, e portanto os resultados, so exatamente iguais aos resultados tericos, como se observa no grfico seguinte, de J.D. Anderson (2001).

Fig. 5 Cp vs j com circulao comeando no bordo de ataque at bordo de fuga O Cp varia do mesmo modo e na mesma faixa de valores, entre 1 e -3 (soluo terica), comprovando a boa aplicao, computao e aproximao do mtodo numrico. ngulo de ataque = 10

Fig. 6 Pontos de controle a amarelo e intensidade do Cp em cada ponto, para = 10 Nesta situao, observa-se que a intensidade do Cp maior e negativa na superfcie superior, enquanto na superfcie inferior positiva, mas de intensidade inferior. Este desequilbrio entre as superfcies comprova a existncia de um gradiente negativo de presso na face superior, fornecendo sustentao ao perfil. Nesta situao, tem-se Clift = 2,1439.


Fig. 7 Cp vs j (inclinao do painel j), para = 10 No grfico da fig. 6 verifica-se, relativamente ao caso anterior de ngulo de ataque nulo, uma assimetria da distribuio do Cp. A sua intensidade mxima negativa deslocou-se para o bordo de ataque, sendo que na superfcie inferior a variao muito tnue, como j se podia ter observado na fig. 5. Est mais uma vez demonstrado que o resultado da integrao do Cp ir resultar em sustentao positiva. Perfil Triangular Este perfil apresentado pela sua geometria particular, que apresenta um bordo de fuga plano, do tamanho da espessura mxima do perfil e perpendicular linha mdia. No programa, possvel escolher o declive dos lados iguais que constituem o tringulo issceles de corda unitria. Segue-se o perfil com declive de 0,25.

Fig. 8 Perfil triangular e painis delimitados painis a verde nos bordos de ataque e fuga


ngulo de ataque = 0

Fig. 9 Pontos de controle a amarelo e intensidade do Cp em cada ponto Aqui, pode ser feita a mesma anlise que no caso do cilindro infinito com ngulo de ataque nulo. Mais uma vez, calculado zero para o Cp, que vem comprovar a sustentao nula para perfis simtricos quando o ngulo de ataque igual a zero. A caracterstica deste perfil de ter uma descontinuidade geomtrica no bordo de fuga, resulta numa perturbao do escoamento diferente em relao aos outros perfis, da o coeficiente de presso ter maior intensidade nesse ponto quando comparado com a restante superfcie do perfil. No grfico da fig. 9, a distribuio do Cp naturalmente igual em ambas as superfcies, da se a sustentao ser nula.

Fig. 10 Cp vs corda, para = 0


ngulo de ataque = 10

Fig. 11 Pontos de controle a amarelo e intensidade do Cp em cada ponto Na situao de ngulo de ataque positivo em 10 graus, expectvel sustentao, no sendo totalmente perceptvel essa situao no grfico da fig. 10. No entanto, observando o grfico da fig. 11, possvel constatar a diferena de distribuio de Cp entre a superfcie superior e inferior, concluindo a existncia de sustentao positiva, calculada no valor de Clift = 1,5183.

Fig. 12 Cp vs corda, para = 10 Outras geometrias poderiam ser analisadas, mas avaliaram-se os perfis de cilindro infinito e triangular para validar o mtodo numrico, e assim considerar os resultados para perfis NACA, que j so perfis com utilidade aeronutica.


Perfil NACA simtrico NACA0012 O perfil NACA0012 um perfil clssico de asa, simtrico e j amplamente estudado, sendo por isso o perfil ideal a ser avaliado pelo programa, testando a os resultados obtidos por este meio com os resultados tericos e experimentais, dada a ampla informao disponvel sobre ele. Para obter a geometria deste perfil e de qualquer outro perfil NACA, foi seguida a formulao matemtica apresentada no captulo 3, utilizando-se uma corda unitria, para facilitar clculos e visualizao grfica na relao x/c, c = 1. Verifica-se, na fig. 13, a simetria do NACA 0012 e a distribuio dos painis com uma maior densidade nos bordos de ataque e de sada, para melhor preciso de clculo, como j referido anteriormente.

Fig. 13 Perfil NACA0012 e painis delimitados painis a verde no bordos de ataque e a vermelho no bordo de fuga Segue-se agora a apresentao dos resultados para ngulo de ataque nulo e a sua comparao com dados tericos, sendo depois avaliada a evoluo do Cp com o ngulo de ataque e a sua comparao com dados experimentais.


ngulo de ataque = 0

Fig. 14 Pontos de controle a amarelo e intensidade do Cp em cada ponto



Naturalmente, observa-se coeficiente de presso negativo na superfcie superior e positivo na superfcie inferior. Mas, sendo o perfil simtrico e o ngulo de ataque nulo, verifica-se a sobreposio da distribuio do Cp, como a fig. 15 mostra. Verificando ainda o grfico de velocidade tangencial adimensionalizada mesma distribuio do Cp, lembrando que Cp = 1 (vt/U)2 e sobrepondo os valores experimentais retirados das tabelas em anexo de I.H. Abbott (1959), verifica-se uma correspondncia de resultados muito boa, como se pode ver na fig. 16 e 17.

Fig. 16 (vt/U)2 vs corda, para = 0

Fig. 17 (vt/U)2 vs corda, para = 0, sobreposio do terico e mtodo numrico


Esta correspondncia de resultados comprova mais uma vez a validao do programa. Assim, admitindo que os resultados obtidos pelo programa tem uma muito boa aproximao dos tericos, interessa agora comparar o comportamento dos perfis NACA em diversas condies, seja variando o perfil, seja variando o ngulo de ataque. Fazendo uma anlise da variao do Cp com o ngulo de ataque, para o perfil NACA0012 tem-se:

Fig. 18 Clift vs = 0, sobreposio do experimental e mtodo numrico Na figura 18, tem-se a vermelho o Clift obtido pelos dados do programa, sobrepondo os dados experimentais retirados de I.H. Abbott (1959). Os dados do programa foram retirados depois de testar o perfil para 16 ngulos de ataque diferentes na faixa [-10;10]. Esses valores foram inseridos num grfico e sobrepostos ao grfico apresentado.


Conclui-se um comportamento linear, idntico entre o mtodo numrico e os valores experimentais na faixa [-10;10] graus de ngulo de ataque. Experimentalmente, para valores de mais elevados, acontece o stall, fenmeno de queda de sustentao quando o valor crtico de atingido por voltas dos 15 graus. O mtodo numrico no consegue seguir este fenmeno, resultado da integrao a continuao do comportamento linear da relao Cp vs . Para todos os perfis analisados a relao ser linear, sendo agora importante comparar o comportamento de perfis diferentes. Perfil NACA no simtrico NACA4212 Aqui, sero apresentados os resultados para este perfil no simtrico, com uma maior sustentao oferecida pela geometria. O objetivo passa por obter igualmente a evoluo do Cp deste perfil com o ngulo de ataque e verificar qual perfil NACA tem maior Cp.

Fig. 19 Perfil NACA4212 e painis delimitados painis a verde no bordos de ataque e a vermelho no bordo de fuga


ngulo de ataque = 0

Fig. 20 Pontos de controle a amarelo e intensidade do Cp em cada ponto



Da fig. 20 e 21 conclui-se que para um ngulo de ataque nulo, este perfil tem fora de sustentao positiva, pois o coeficiente de presso na superfcie superior de intensidade negativa, que em mdulo muito maior que a intensidade na superfcie inferior. visvel que pouco depois do bordo de ataque, quando a corda cerca de 0,1, o Cp atinge intensidades negativas muito acentuadas no perfil superior. Na integrao do Cp resulta um Clift = 0,4518. Assim, prev-se desde j que o Cp variando com o ngulo de ataque ser superior neste perfil quando comparado com o perfila NACA simtrico analisado anteriormente, pois em = 0 ir cruzar o eixo do Cp num valor j positivo, o que se verificar mais frente, quando se compararem os dois perfis.

ngulo de ataque = 10

Fig. 23 Pontos de controle a amarelo e intensidade do Cp em cada ponto Nesta situao de ngulo de ataque negativo, a sustentao ser naturalmente negativa, registando-se um Clift = 0,7483. O utilizador do programa pode efetuar esta avaliao feita ao longo do relatrio, observando todos os grficos para vrios perfis e os resultados guardados nos ficheiros de texto.


Comparao entre NACA0012 e NACA4212

Fig. 24 Clift e Cdrag dos perfis NACA4212 e NACA0012 Na fig. 24 verifica-se que a relao linear Clift com do NACA4212 est deslocada cerca de 50% acima da do NACA0012. portanto um perfil mais eficaz aerodinamicamente. importante fazer tambm uma referencia ao coeficiente de arrasto. Este coeficiente aerodinmica no explorado neste trabalho, uma vez que este mtodo numrico no o mais apropriado para o seu clculo, sendo at que teoricamente tem que se nulo nesta anlise. No entanto, o programa devolve o seu valor, pois possvel integrar o Cp de modo a calcular o Cdrag. O clculo aproxima muito bem a teoria, pois devolve valores nulos nas duas primeiras casas decimais, podendo observar-se que praticamente nulo e totalmente desprezvel no grfico da fig. 24.


6. Concluso A principal observao que se pode fazer sobre este trabalho a validao do mtodo numrico como uma excelente aproximao da integrao terica para aeroflios. Os resultados obtidos aproximam muito bem os resultados experimentais, sendo por isso um mtodo prtico e fcil para a modulao de aeroflios, nomeadamente asas de avio. O objetivo deste trabalho passava por entender o mtodo dos painis, aplicando-o a diversas formas geomtricas, o que foi cumprido tendo o programa respondido muito bem na sua aplicao do mtodo. Foi possvel observar as diferenas aerodinmicas inerentes variabilidade geomtrica dos vrios perfis considerados e confirmar o seus resultados com dados experimentais e tericos. O programa de muito fcil utilizao e com ele imediata a obteno de resultados para qualquer NACA 4-series pretendido, concluindo-se que vivel e confivel.

7. Bibliografia

Livros

Anderson Jr.,J.D., Fundamentals Of Aerodynamics, McGraw-Hill, 3ed, New York, USA, 2001 Abbott,I.H., von Doenhoff,A.E., Theory Of Wing Sections, Dover Publications, 2ed, New York, USA, 1959 Kuethe,A.M, Chow,C., Foundations OF Aerodynamics: Bases Of Aerodynamic

Design, John Wiley & Sons Inc, 5ed, Exeter, UK, 1997 Stios da Internet

AEROSPACEWEB (http://www.aerospaceweb.org/question/airfoils/q0041.shtml)


8. Anexo Cdigo do Programa em MATLAB

Documents

Método de Painéis de Vórtice 2D - Aerodinâmica - Cálculo Numérico MATLAB