METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII...

37
Universitatea din Bucure¸ sti METODE ASIMPTOTICE S ¸I NUMERICE ˆ IN STUDIUL ARIPII PORTANTE ˆ IN CURENT SUBSONIC Tezˇ a de doctorat Coordonator ¸ stiint ¸ific: Adrian STOICA Prof. dr. Sanda T ¸ IGOIU septembrie 2015

Transcript of METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII...

Page 1: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Universitatea din Bucuresti

METODE ASIMPTOTICE SI NUMERICE IN

STUDIUL ARIPII PORTANTE IN CURENT

SUBSONIC

Teza de doctorat

Coordonator stiintific:

Adrian STOICA Prof. dr. Sanda TIGOIU

septembrie 2015

Page 2: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Introducere

Scopul acestei lucrari este de a rafina metodele matematice de studiu ale miscarii in-

compresibile ın jurul aripilor planare si extinderea lor la cazul aripilor nonplanare. Daca

prin aripa portanta planara ıntelegem o suprafata ”putin” diferita de o suprafata plana,

care poate genera o forta de sustentatie, prin aripa nonplanara vom ıntelege o suprafata

”putin” diferita de o suprafata cilindrica avand generatoarea paralela cu directia curentului

neperturbat. In categoria aripilor nonplanare se ıncadreaza aripioarele(winglets), biplanele,

parasutele, aripile inelare. Unul din motivele importante pentru care sunt studiate aripile

nonplanare rezida ın faptul ca acestea au caracteristici aerodinamice, ın termeni de raport

ıntre coeficientul de rezistenta la ınaintare si patratul coeficientul de portanta, mai bune

decat aripile planare, cu alte cuvinte la anverguri si portante egale, exista aripi nonplanare

a caror rezistenta la ınaintare este mai mica decat cea a aripii planare de rezistenta minima.

S-a constatat experimental ca montarea cate unei aripioare avand forma unei spire(spiroid

winglet) la bordurile laterale ale aripilor unui avion Falcon 50, a condus la o micsorare a

coeficientului de rezistenta la ınaintare cu peste 11% cu consecinte proportionale ın ceea ce

priveste consumul de carburant.

Rezultatele originale din aceasta teza s-au concretizat prin publicarea sau trimiterea spre

publicare a 4 articole stiintifice ın reviste de specialitate:

1)A.Stoica, A Numerical Solution for the Lifting Surface Integral Equation, INCAS BUL-

LETIN, Volume 6, Special Issue 1/ 2014

2)A.Stoica, A rigorous derivation of certain equations arising in the lifting wing theory,

Bulletin Mathematique de la Societe des Sciences Mathematiques de Roumanie (trimis spre

publicare)

3)A.Stoica, An exact solution for the minimum induced drag of an circular arched wing

and the validation of a conjecture, Journal of Engineering Mathematics (trimis spre publi-

care)

1

Page 3: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

4)A.Stoica, A singular integral equation for curved wings of high aspect ratio, Proc. Ro.

Acad. Series A (trimis spre publicare)

2

Page 4: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Capitolul 1

Miscarea ın jurul aripii portante

In primul capitol extindem la cazul aripii nonplanare metoda solutiilor fundamentale

ale sistemului aerodinamicii liniarizate, dezvoltata de Dragos si Homentcovschi ın [6], [7],

[9] pentru aripa planara. Prin aceasta metoda, care utilizeaza aparatul matematic al teoriei

distributiilor, se arata pe de o parte ca perturbatia este potentiala ın spatele aripii, exceptand

o suprafata plana (cilindrica ın cazul aripii nonplanare) care se ıntinde spre infinit aval pe

directia curentului neperturbat, iar pe de alta parte se determina perturbatia presiunii prin

intermediul unui produs de convolutie care reprezinta de fapt o superpozitie continua de

dublete avand ın fiecare punct directia data de normala la suprafata cilindrica de referinta

si intensitatea ın fiecare punct ca fiind diferenta f de presiune de pe extrados si intrados.

Lucrurile acestea doua se postuleaza ın cadrul teoriei clasice Multhopp [12], Ashley si Landahl

[1]. Am aflat de asemenea reprezentarea componentelor campului de viteze de perturbatie

si am calculat coeficientii aerodinamici cu ajutorul carora se determina actiunea fluidului

asupra aripii. Determinarea saltului de presiune f, atrage dupa sine cunoasterea miscarii

perturbate precum si caracteristicile aerodinamice ale aripii.

3

Page 5: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Se considera miscarea stationara a unui fluid incompresibil nevascos, in jurul unei aripi

fara grosime, de anvergura finita. La mari distante, ın amonte, miscarea este uniforma cu

viteza constanta V∞. Alegem un sistem cartezian ortogonal OXY Z si {O, i, j,k} reperul

asociat astfel ıncat axele OX si OZ coincid cu directia curentului neperturbat respectiv

directia portantei. Notam cu L0, V∞, p∞, ρ respectiv lungimea maxima a semicorzii aripii,

modulul vitezei, presiunea si densitatea fluidului la infinit amonte. De asemenea, vom nota

cu x, y, z coordonatele adimensionale, legate de cele cu dimensiune x1, y1, z1, prin relatiile

(x1, y1, z1) = L0(x, y, z). Daca nu se specifica altfel, peste tot ın acest capitol prezenta in-

dicelui inferior ”1” arata fie ca obiectele matematice sunt descrise ın raport cu variabilele cu

dimensiune sau ca marimile respective au dimensiune. Introducem perturbatiile adimension-

ale p si v ale presiunii si respectiv vitezei, legate de presiunea p1 si viteza v1 din miscarea

perturbata (ın raport cu miscarea uniforma), prin relatiile:

p1 = p∞ + ρV 2∞p, v1 = V∞(i + v). (1.1)

Facem urmatoarele ipoteze:

1.Geometria aripii(Figura 1.1)

Suprafata suport W a aripii este putin deviata fata de o suprafata cilindrica Σ, directia

generatoarei fiind cea a curentului neperturbat, iar curba directoare C, presupus a fi simpla,

neteda de clasa C2, regulata, situata ıntr-un plan perpendicular pe OX, este descrisa de

ecuatiile:

y = ys (η) , z = zs (η) , −θ ≤ η ≤ θ. (1.2)

Daca zs ≡ 0, suntem ın cazul clasic al aripii planare, ın caz contrar fiind vorba despre

aripa nonplanara. Aceasta ultima situatie corespunde aripii ınchise, daca C este o curba

inchisa. Notam cu S proiectia ortogonala a suprafetei W a aripii pe suprafata cilindrica Σ.

Pentru fiecare η ∈ [−θ, θ], desemnam prin x−(η) si x+(η) componenta dupa OX a punc-

tului situat pe bordul de atac, respectiv bordul de fuga. Ecuatiile care descriu suprafetele

W si S sunt:

W :

x = ξ

y = ys (η) + yε (ξ, η) , (ξ, η) ∈ D,z = zs (η) + zε (ξ, η)

(1.3)

4

Page 6: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Curba CProiectia

suprafetei W

Y

Z

(a) Proiectia aripii si a suprafetei de referinta Σpe planul Y OZ

X

Y

Bordul de atac

𝒙 = 𝒙−

Bordul de fugă

𝒙 = 𝒙+

(b) Proiectia aripii pe planul XOY

Figura 1.1: Geometria aripii

S :

x = ξ

y = ys (η) , (ξ, η) ∈ D,z = zs (η)

(1.4)

unde

D = {(ξ, η) | x− (η) ≤ ξ ≤ x+ (η) ,−θ ≤ η ≤ θ}, (1.5)

5

Page 7: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

iar ε � 1 este un parametru cu ajutorul caruia descriem deviatia suprafetei W a aripii de

la suprafata cilindrica Σ.

Admitem ca x−, x+, yε si zε sunt functii netede astfel ıncat |yε| , |zε|,∣∣∣∣∂yε∂ξ

∣∣∣∣ , ∣∣∣∣∂yε∂η∣∣∣∣ , ∣∣∣∣∂zε∂ξ

∣∣∣∣ , ∣∣∣∣∂zε∂η∣∣∣∣

sunt de ordinul O(ε) pe D.

2. Ipoteza micilor perturbatii

Datele problemei ne permit sa admitem ipoteza micilor perturbatii |p| � 1, |v| � 1 ın tot

domeniul ocupat de fluid, mai mult presupunem ca perturbatia p a presiunii si componentele

u, v, w ale perturbatiei vitezei sunt de ordinul O(ε) si ca perturbatia dispare la infinit

lim v = 0, lim p = 0, cand ‖(x, y, z)‖ → ∞. (1.6)

3. Aripa portanta

Exista o diferenta ıntre presiunile de pe extradosul si intradosul aripii, astfel ca daca

notam cu p+ si p− valorile limita ale presiunii pe fata exterioara respectiv interioara a

suprafetei W, atunci

[[p]] = p+(x, y, z)− p−(x, y, z) = f(x, y, z), (∀)M(x, y, z) ∈ W. (1.7)

Ecuatiile de miscare liniarizate

In variabile adimensionale, dupa liniarizarea ın jurul miscarii neperturbate, ecuatia de

continuitate si ecuatia Euler pentru miscarea stationara formeaza sistemul de ecuatii al

aerodinamicii liniarizate stationare (Homentcovschi [9], Dragos [6]):

divv = 0,∂v

∂x+ gradp = 0.

(1.8)

Conditia de alunecare liniarizata:

y′

sw − z′

sv =∂zε∂ξ

y′

s −∂yε∂ξ

z′

s pe S. (1.9)

6

Page 8: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Ecuatia de salt liniarizata a perturbatiei presiunii se scrie:

[[p]] = p+(x, y, z)− p−(x, y, z) = f(x, y, z), (∀)M(x, y, z) ∈ S. (1.10)

Homentcovschi [9] a aratat ca daca S este o suprafata (solida) de discontinuitate a

campului de presiune, atunci sistemul ın distributii asociat sistemului aerodinamicii liniarizate

(1.8) si conditiei de salt (1.10) este:divv = 0∂v

∂x+ gradp = fδSn,

(1.11)

unde n = (nx, ny, nz) este versorul normalei exterioare la suprafata S, iar fδSn reprezinta

distributia de simplu strat pe suprafata S cu densitatea f . Aici nx=0.

Definitia 1. Solutiile ın distributii ale sistemuluidivv = 0∂v

∂x+ gradp = f

(1.12)

se numesc solutii fundamentale (Homentcovschi [9] si Dragos [6]) ale sistemului aerodina-

micii liniarizate.

Folosind solutia fundamentala a sistemului (1.12) (Dragos [7], Cap 2) pentru cazul ın

care f = fδSn obtinem ca

v = H(x) ∗ fδSn + gradϕ, (1.13)

unde H este distributia Heaviside, ϕ este potentialul perturbatiei vitezei,

ϕ (x, y, z) = − 1

∫∫S

f (x1, y1, z1)(r− r1) · n1

(y − y1)2 + (z − z1)2

(1 +

x− x1|r− r1|

)dσ1, (1.14)

iar r = xi + yj + zk si r1 = x1i + y1j + z1k. Se observa din ecuatia (1.13), ca miscarea este

potentiala ın afara suprafetei Σ′ ⊂ Σ care pleaca de la bordul de atac si continua spre infinit

aval.

Teorema 1. Perturbatia miscarii uniforme, datorata prezentei aripii portante, este descrisa,

ın orice punct M(x, y, z) ∈ R3 \Σ′, prin campul de viteze v, de componente u = −p, iar v si

7

Page 9: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

w definite prin (1.15) si respectiv (1.16). Distributia de presiune este definita prin formula

(1.17).

v (x, y, z) =− 1

∫∫S

f (x1, y1, z1)n1y

(y − y1)2 + (z − z1)2

(1 +

x− x1|r− r1|

)dσ1

+1

∫∫S

f (x1, y1, z1) (y − y1)2 n1y[(y − y1)2 + (z − z1)2

]2 (1 +x− x1|r− r1|

)dσ1

+1

∫∫S

f (x1, y1, z1) (y − y1) (z − z1)n1z[(y − y1)2 + (z − z1)2

]2 (1 +

x− x1|r− r1|

)dσ1

+1

∫∫S

f (x1, y1, z1) [(y − y1)n1y + (z − z1)n1z]

(y − y1)2 + (z − z1)2(x− x1) (y − y1)|r− r1|3

dσ1

(1.15)

w (x, y, z) =− 1

∫∫S

f (x1, y1, z1)n1z

(y − y1)2 + (z − z1)2

(1 +

x− x1|r− r1|

)dσ1

+1

∫∫S

f (x1, y1, z1) (z − z1)2 n1z[(y − y1)2 + (z − z1)2

]2 (1 +x− x1|r− r1|

)dσ1

+1

∫∫S

f (x1, y1, z1) (y − y1) (z − z1)n1y[(y − y1)2 + (z − z1)2

]2 (1 +

x− x1|r− r1|

)dσ1

+1

∫∫S

f (x1, y1, z1) [(y − y1)n1y + (z − z1)n1z]

(y − y1)2 + (z − z1)2(x− x1) (z − z1)|r− r1|3

dσ1

(1.16)

p(x, y, z) =1

∫∫S

f(x1, y1, z1)∂

∂n1

(1

|r− r1|

)dσ1. (1.17)

8

Page 10: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Reprezentarea data prin formulele (1.17),(1.15) si (1.16) au obtinut-o Ashley si Landahl

[1], prin asimilarea aripii (de fapt a proiectiei S) cu o distributie continua de dublete cu axele

normale la suprafata S. In reprezentarea obtinuta de noi, am plecat de la ipoteza existentei

saltului de presiune pe aripa care creeaza forta de sustentatie.

Actiunea fluidului asupra aripii. Coeficientii aerodinamici

Pentru a calcula actiunea fluidului asupra aripii, integram pe suprafata S1 saltul campului

de presiune din miscarea perturbata. Astfel, forta rezultanta esta data de

R1 = −∫∫

S1

[[p1]]ndσ1 = −ρV 2∞L

20

∫∫S

fndσ,

iar momentul rezultant datorat interactiunii fluidului cu aripa se defineste prin

M1 = −∫∫

S1

r1 × [[p1]]ndσ1 = −ρV 2∞L

30

∫∫S

r× fndσ,

unde n este versorul normalei exterioare la suprafata Σ.

Atunci, coeficientii aerodinamici adimensionali sunt:

Coeficientul de portanta

CL =R1z

1/2ρV 2∞Aria(S1)

= − 2

A

∫∫D

y′s (η) f (ξ, η) dξdη. (1.18)

Coeficientul de rezistenta la ınaintare

CD =R1x

1/2ρV 2∞Aria(S1)

= − 2

A

∫∫D

(∂yε∂ξ

z′

s (η)− ∂zε∂ξ

y′

s (η)

)f (ξ, η) dξdη. (1.19)

Coeficientul fortelor laterale

CY =R1y

1/2ρV 2∞Aria(S1)

=2

A

∫∫D

z′s (η) f (ξ, η) dξdη. (1.20)

unde A = Aria(S) =

∫∫D

√y′2s (η) + z′2s (η)dξdη.

9

Page 11: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Momentul rezultant datorat interactiunii fluidului cu aripa se defineste prin

M1 = −∫∫

S1

r1 × [[p1]]ndσ1 = −ρV 2∞L

30

∫∫S

r× fndσ,

de unde obtinem coeficientii de moment:

Momentul de ruliu

cx =M1x

1/2ρV 2∞Aria(S1)a1

= − 2

Aa0

∫∫D

(zs(η)z

s(η) + ys(η)y′

s(η))

f (ξ, η) dξdη. (1.21)

Momentul de tangaj

cy =M1y

1/2ρV 2∞Aria(S1)a1

=2

Aa0

∫∫D

ξy′

s(η)f (ξ, η) dξdη. (1.22)

Momentul de giratie

cz =M1z

1/2ρV 2∞Aria(S1)a1

=2

Aa0

∫∫D

ξz′

s(η)f (ξ, η) dξdη, (1.23)

unde a1 si respectiv a0 sunt lungimea medie cu dimensiune, respectiv fara dimensiune, a

corzii aripii (a1 = L0a0).

Daca ın formulele 1.18 si 1.19, loc de A, aria suprafetei aripii S, punem A′, aria supafetei

proiectiei aripii pe planul XOY, obtinem coeficientii C ′L si respectiv C ′D cu ajutorul carora

definim coeficientul de eficienta aerodinamica ın raport cu anvergura (span effi-

ciency factor) prin relatia

C ′D =C′2L

πCefAR, (1.24)

unde raportul de aspect AR este definit ca raportul dintre patratul anvergurii si aria

proiectiei aripii pe planul XOY.

Concluzie: Formulele (1.15)-(1.23) arata ca pentru determinarea perturbatiei si a actiunii

fluidului asupra aripii este suficient sa determinam saltul de presiune f .

10

Page 12: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Capitolul 2

Aripa portanta planara de anvergura

finita

In cazul aripii planare, pentru a determina valoarea limita a componentei w a perturbatiei

vitezei, am propus, la ınceputul capitolului al doilea, o metoda riguroasa care evita operatii

formale de permutare a limitei cu integrala sau cu derivata. Prin impunerea conditiei de

alunecare a fluidului la suprafata aripii, se obtine o ecuatie integrala pentru functia necunos-

cuta f . Aceasta ecuatie prezinta o singularitate tare, de ordinul al doilea, de tip Hadamard.

In literatura am ıntalnit o singura metoda numerica de rezolvare a ecuatiei suprafetei

portante a carei convergenta este demonstrata [2]. Schema numerica respectiva este ınsa

lent convergenta, avand ordinul O( 1nα

), 0 < α < 1. In cazul aripii dreptunghiulare, prop-

unem o metoda numerica de rezolvare bazata pe formule de cuadratura de tip Gauss pentru

integrale cu singularitate Cauchy. Simularile numerice sugereaza un ordin de convergenta

de ordinul O(1/n3.5). Metoda se dovedeste a fi eficienta si ın cazul miscarii ın jurul aripii

dreptunghiulare ın prezenta solului.

Pentru linia planara curba, pornind de la dezvoltarea asimptotica a componentei w a

perturbatiei vitezei(data de Guermond [8]), am plasat linia portanta astfel ıncat ın orice

sectiune transversala, momentul act, iunii fluidului fat, a de un un punct situat pe linia por-

tanta sa fie zero. In cazul aripii de incident, a constanta, aripa este plasata la 1/4 din lungimea

corzii fat, a de bordul de atac. Ecuatia obtinuta este o generalizare a ecuatiei liniei portante

drepte a lui Prandtl. Ea surprinde atat curbura cat si ” sageata ” aripii. Am propus de

asemenea o metoda numerica de rezolvare a acestei ecuatii.

11

Page 13: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Pentru aripa planara, suprafata de referinta Σ este planul XOY. Notam cu D proiectia

suprafetei aripii pe planul XOY, numita planforma aripii(figura 2.1). In acest caz relatiile

(1.3) si (1.4) care descriu geometria aripii sunt:

W : z = h(ξ, η), (ξ, η) ∈ D (2.1)

este ecuatia carteziana a suprafetei aripii, iar suprafata S este caracterizata de:

ys(η) = η, zs(η) = 0, (ξ, η) ∈ D,

D = {(ξ, η)|x−(η) ≤ ξ ≤ x+(η),−b ≤ η ≤ b}.(2.2)

Admitem cu a functia h este neteda pe D, iar h,∂h

∂xsi∂h

∂yau ordinul O(ε) pe D, unde

ε� 1 este parametrul care modeleaza deviatia suprafetei de la planul XOY.

In ecuatiile de mai sus b este semianvergura aripii. Am admis ca bordurile laterale ale

planformei fie sunt segmente de dreapta paralele cu directia curentului neperturbat, fie pot

fi reduse la un punct.

𝐕∞

X

Y

an

ver

gu

ra

Bord de atac

Bord de fugă

coarda

-b

b

Figura 2.1: Planforma aripii

12

Page 14: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Campurile de viteza si presiune.

In relatiile (1.17),(1.15) si (1.16) punem

r = xi + yj + zk, r1 = ξi + ηj. (2.3)

Introducem notatiile

x0 = x− ξ, y0 = y − η, R1 = |r − r1| =√x20 + y20 + z2. (2.4)

Deoarece in acest caz n1y = 0, n1z = 1, din formulele generale (1.15),(1.16) si(1.17) regasim

formulele date de Dragos [7] pentru componentele perturbatiei vitezei ın fluid.

Propozitia 1. Perturbatia miscarii uniforme, datorata prezentei aripii portante planare,

este descrisa prin campul de viteze v, de componente u, v, w si distributia de presiune p

definite prin

u(x, y, z) = − z

∫∫D

f(ξ, η)

R31

dξdη, (2.5)

v(x, y, z) =z

∫∫D

f(ξ, η)y0(y20 + z2)2

(1 +

x0R1

)dξdη +

z

∫∫D

f(ξ, η)x0y0(y20 + z2)R3

1

dξdη, (2.6)

w(x, y, z) = − 1

∫∫D

f(ξ, η)y20 − z2

(y20 + z2)2

(1 +

x0R1

)dξdη+

z2

∫∫D

f(ξ, η)x0(y20 + z2)R3

1

dξdη, (2.7)

p(x, y, z) = −u(x, y, z). (2.8)

Tinand cont de (2.1), conditia de alunecare (1.9) devine

w(x, y, 0) =∂h

∂x(x, y), (∀) (x, y) ∈ D. (2.9)

13

Page 15: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Deducerea ecuatiei suprafetei portante

Pentru a afla valorile perturbatiei componentei w pe aripa, trebuie sa trecem la limita

cand z → 0. Prin simpla ınlocuire z = 0 ın formula (2.7), integralele devin singulare.

Definitia 2. Pentru orice functie test ϕ ∈ D (R2), supp ϕ ⊂ [−A,A]× [−A,A], definim⟨FP 1

y2

(1 +

x√x2 + y2

), ϕ

⟩= lim

ε→0+

[∫∫Dε

ϕ (x, y)

y2

(1 +

x√x2 + y2

)dxdy−

−2

ε

∫ A

−Aϕ(x, 0)(1 + sgnx)dx

],

(2.10)

unde Dε = [−A,A]× {[−A,A] \ [−ε, ε]} .

Propozitia 2. Pentru orice functie test ϕ ∈ D (R2), supp ϕ ⊂ [−A,A] × [−A,A], are loc

egalitatea⟨FP 1

y2

(1 +

x√x2 + y2

), ϕ

⟩=

=

A∫−A

A∫−A

ϕ(x, y)− ϕ(x, 0)− y∂ϕ∂y

(x, 0)

y2

(1 +

x√x2 + y2

)dxdy−

− 2

A

A∫−A

ϕ(x, 0)− ϕ(0, 0)

x

(x+√A2 + x2

)dx− 4ϕ(0, 0).

(2.11)

Propozitia 3. In spatiul distributiilor D′(R2) au loc urmatoarele limite:

a) limz→0

y2 − z2

(y2 + z2)2

(1 +

x√x2 + y2 + z2

)= FP 1

y2

(1 +

x√x2 + y2

).

b) limz→0

xz2

(y2 + z2)(x2 + y2 + z2)3/2= 0.

(2.12)

Teorema 2. Daca f este o functie local integrabila cu suport compact, atunci

limz→0

w(x, y, z) = − 1

∗∫∫D

f(ξ, η)

(y − η)2

[1 +

x− ξ√(x− ξ)2 + (y − η)2

]dξdη. (2.13)

Pe baza conditiei de alunecare (2.9) obtinem ecuatia suprafetei portante, pentru functia

14

Page 16: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

necunoscuta f, data prima data de Multhopp (1955):

− 1

∫∫ ∗D

f(ξ, η)

(y − η)2

(1 +

x0R

)dξdη =

∂h

∂x(x, y), (x, y) ∈ D. (2.14)

Definitia 3. Pentru orice functie test ϕ ∈ D (R2), supp ϕ ⊂ [−A,A]× [−A,A], definim⟨PV

[1

y

(1 +

√x2 + y2

x

)], ϕ

⟩= lim

ε→0+

∫∫D1ε

ϕ (x, y)

y

(1 +

√x2 + y2

x

)dxdy (2.15)

unde D1ε = ([−A,−ε] ∪ [ε, A])× ([−A,−ε] ∪ [ε, A]).

Propozitia 4.

∂yPV

[1

y

(1 +

√x2 + y2

x

)]= −FP

[1

y2

(1 +

x√x2 + y2

)], ın D′(R2). (2.16)

Obtinem astfel o alta forma ecuatiei suprafetei portante, care contine doar singularitati

de tip Cauchy:

1

∂y

∫∫ ∗D

f(ξ, η)

y − η

(1 +

R

x0

)dξdη =

∂h

∂x(x, y), (x, y) ∈ D. (2.17)

O metoda numerica ın cazul aripii dreptunghiulare

Vom da o alta forma a ecuatiei suprafetei portante. Pornim de la ecuatia (2.17) si facem

urmatoarele ipoteze:

i)Axa OX este axa de simetrie pentru domeniul D;

ii)D este multime convexa.

Propozitia 5. In ipotezele i) si ii) de mai sus, o alta forma a ecuatiei suprafetei portante

este:

1

∫∫ ∗D

f (ξ, η)

y − η

1 +

√(x− ξ)2 + (y − η)2

x− ξ

dξdη =

∫ y

0

∂h

∂x(x, t)dt, (x, y) ∈

◦D. (2.18)

Daca D=[−1, 1]× [−b, b] atunci ipotezele i) si ii) sunt ındeplinite. Efectuam schimbarea

de variabile

x = x′, y = by′, ξ = ξ′, η = bη′, (2.19)

15

Page 17: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

dar pentru usurinta scrierii vom renunta la semnul ”prim” si vom nota noile variabile la fel

ca pe cele vechi. Tinand cont de conditiile la frontiera pe care trebuie sa le ındeplineasca

saltul de presiune f , cautam solutii(ın vechile variabile) de forma

f(x, y) =

√x+(y)− xx− x−(y)

√b2 − y2F (x, y), (x, y) ∈ D. (2.20)

In noile variabile, ecuatia (2.18) se scrie

1

∫ ∗1−1

√1− η2y − η

(∫ ∗1−1

√1− ξ1 + ξ

F (ξ, η)

[1 +

√(x− ξ)2 + b2(y − η)2

x− ξ

]dξ

)dη =

=

∫ by

0

∂h

∂x(x, t)dt.

(2.21)

Impunem ca ecuatia (2.21) sa fie verificata ın punctele de control (xi, yj), unde

xi = cos

(2i− 1

2m+ 1π

), i = 1,m

yj = cos

(2j − 1

n+ 1

π

2

), j = 1, n+ 1.

(2.22)

Folosim formulele de cuadratura, vezi Dragos [7]:

∫ ∗1−1

√1− ξ1 + ξ

f(ξ)

ξ − xidξ ' 2π

2m+ 1

m∑k=1

1− ξkξk − xi

f(ξk),

∫ ∗1−1

√1− η2 f(η)

η − yjdη ' π

n+ 1

n∑l=1

1− η2lηl − yj

f(ηl),

unde

ξk = cos

(2kπ

2m+ 1

), k = 1,m,

ηl = cos

(2l − 1

n+ 1

π

2

), l = 1, n+ 1

(2.23)

reprezinta punctele de colocatie.

Dupa aplicarea formulelor de cuadratura de mai sus ın ecuatia 2.21, obtinem un sistem

16

Page 18: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

liniar de forma

m∑k=1

n∑l=1

AijklFkl = Hij, i = 1,m, j = 1, n+ 1, (2.24)

unde Fkl = F (ξk, ηl) sunt necunoscute, Hij =

∫ byj

0

∂h

∂x(xi, t)dt sunt cunoscute, iar coeficientii

cunoscuti ai matricei sistemului sunt

Aijkl =π

2(2m+ 1)(n+ 1)

(1− η2l )(1− ξk)yj − ηl

[1 +

√(xi − ξk)2 + b2(yj − ηl)2

xi − ξk

]. (2.25)

Sistemul de ecuatii liniare (2.24) este supradeterminat, avand m(n + 1) ecuatii si mn

necunoscute. Pentru a echilibra numarul de ecuatii cu numarul de necunoscute, introducem

m necunoscute fictive {Ci}i=1,m. astfel ca obtinem sistemul patratic

m∑k=1

n∑l=1

AijklFkl + Ci = Hij, i = 1,m, j = 1, n+ 1. (2.26)

Prin rezolvarea sistemului (2.26) obtinem Ci = 0, i = 1,m, ceea ce arata ca sistemul (2.24)

este compatibil.

Rezultate numerice

Placa plana de incidenta constanta α si raport de aspect AR = 1

Am ales ales acelasi numar n de puncte de colocatie atat pe coarda cat si pe anvergura.

O fitare a datelor sugereaza un comportament de forma

CL/α = 1.4602265 + 0.1169319/n3.5. (2.27)

Valoarea CL/α = 1.4602265 coincide cu cea obtinuta de Tuck [19] pana la a saptea cifra

semnificativa.

Placa plana de anvergura mica (low aspect ratio) si incidenta constanta α

In Dragos [7] este data solutia exacta ın acest caz. Astfel CL/α = πb/2. Pentru b = 0.1,

solutia exacta este CL/α = 0.15708 ın timp ce pentru metoda propusa de noi

CL/α = 0.15702.

17

Page 19: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Linia portanta curba pentru aripi planare

C

Y

X

Bord de fugă

Bord de atac

Linie portantă

curbă

Figura 2.2

Alegem ca lungime caracteristica semianvergura aripii (jumatatea distantei ıntre bor-

durile laterale ale aripii). Atunci domeniul D poate fi descris de:

D = {(x, y)| x−(y) ≤ x ≤ x+(y), y ∈ [−1, 1]} . (2.28)

Facem ipoteza ca anvergura aripii este mult mai mare decat maximul corzii profilului, ceea

ce ınseamna ca

ε = supy∈[−1,1]

|x+(y)− x−(y)| � 1. (2.29)

Aceasta ipoteza permite, ın contextul liniei portante, asimilarea aripii cu un segment de linie

18

Page 20: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

dreapta sau curba. Convenim sa reprezentam bordurile de atac si de fuga prin urmatoarele

ecuatii

x−(y) = x0(y) + εxa(y), x+(y) = x0(y) + εxf (y), y ∈ [−1, 1] , (2.30)

unde functiile y 7→ xa(y) si y 7→ xf (y) sunt de ordinul O(1) si de clasa C1(−1, 1).

x = x0(y) reprezinta ecuatia carteziana a unei linii curbe L care va juca rolul segmentului

de dreapta cu care Prandtl a asimilat aripa. Locul ın care este situata linia L ıl vom preciza

mai tarziu. Admitem ca functia y 7→ x0(y) ımpreuna cu prima si a doua derivata sunt de

ordinul O(1) pe intervalul (-1,1).

Introducem functia necunoscuta Γ, circulatia vitezei de-a lungul profilului obtinut prin

sectionarea aripii cu un plan Y = y. Dupa liniarizare, urmeaza ca

Γ(y) =

∫ x+(y)

x−(y)

(p+ − p−)dx =

∫ x+(y)

x−(y)

f(x, y)dx, (2.31)

Daca bordurile laterale se reduc la un punct, atunci x−(±1) = x+(±1), iar daca bordurile

laterale sunt segmente paralele cu directia curentului neperturbat, atunci din conditia Kutta-

Jukovski f(x,±1) = 0. Prin urmare, pentru ambele tipuri de aripa pe care le luam ın

considerare

Γ(1) = Γ(−1) = 0. (2.32)

Admitem ca functia Γ este derivabila pe intervalul (-1,1), iar derivata sa este functie Holder

pe intervalul [-1,1].

Notam

tgΛ(y) = x0′(y),

c(y) =x+(y)− x−(y)

2,

P (y) =c(y)

2r(y)

[2k(y) + tg2Λ(y) + ln

4

e2c(y) cos2 Λ(y)

],

Q(y) = c(y)

[2k(y) sin Λ(y) + sin Λ(y)ln

4

ec(y) cos2 Λ(y)− ln

1 + sin Λ(y)

cos Λ(y)

],

j(y) = 2

∫ x+(y)

x−(y)

√x− x−(y)

x+(y)− x∂h

∂x(x, y)dx.

(2.33)

19

Page 21: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Teorema 3. Ecuatia integro-diferentiala a liniei portante curbe este

Γ(y)

cos Λ(y)− c(y)

2

∗1∫−1

Γ (η)

(y − η)2

1 +x0 (y)− x0 (η)√

(x0 (y)− x0 (η))2 + (y − η)2

dη+

+P (y)Γ(y) +Q(y)Γ′(y) = j(y), −1 < y < 1,

(2.34)

unde functia necunoscuta Γ reprezinta circulatia vitezei ın jurul profilului obtinut prin intersectia

aripii cu planul η = η0, iar functiile Λ, c, P, Q si j sunt definite de formulele (2.33) de mai

sus.

O solutie numerica a ecuatiei liniei portante curbe

Vom da o metoda numerica de rezolvare a unei ecuatii integro-diferentiale singulare de

forma:

a(y)Γ(y) + b(y)Γ′(y) + c(y)

∫ ∗1−1

Γ(η)

(y − η)2[1 + g(y, η)sign(y − η)] dη = j(y),

pentru − 1 < y < 1, cu conditiile Γ(−1) = Γ(1) = 0.

(2.35)

Cautam solutii de forma Γ(y) =√

1− y2 C(y). Deoarece ın evaluarea numerica a integralei

intervin valorile functiei Γ ın nodurile Cebasev de speta a doua,

yk = cos tk, tk =kπ

n+ 1, k = 1, n,

aproximam functia C cu polinomul de interpolare Lagrange care trece prin aceste noduri.

Asadar vom cauta o solutie aproximativa de forma

Γ(y) =√

1− y2n∑k=1

CkUn(y)

(y − yk)U ′n(yk), Ck = C(yk). (2.36)

Notam

Γk = Γ(yk) =√

1− y2k Ck, ak = a(yk), bk = b(yk),

ck = c(yk), jk = j(yk), gkl = g(yk, yl), k, l = 1, n.(2.37)

Se cunoaste din [14] ca

Γ′(y) = − 2

n+ 1

n∑i=1

Cisin tisin t

n∑l=1

l sin lti cos lt, y = cos t, 0 < t < π. Integrala care apare

20

Page 22: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

ın ecuatia (2.35) contine o discontinuitate de speta ıntai ın η = y. Pentru a o evalua numeric,

mai ıntai o regularizam prin ”extragerea” singularitatii si apoi aplicam formula de integrare

numerica∫ ∗1−1

√1− η2 C(η)

(yk − η)2dη ' π

n+ 1

n∑j 6=k

[1− (−1)j+k]1− y2j

(yk − yj)2Cj −

π(n+ 1)

2Ck.

Astfel ∫ ∗1−1

Γ(η)

(yk − η)2g(yk, η)sign(yk − η)dη =

=

∫ ∗1−1

√1− η2

(yk − η)2[C(η)g(yk, η)− Ckgkk]sign(yk − η)dη+

+ Ckgkk

∫ ∗1−1

√1− η2

(yk − η)2sign(yk − η)dη =

n+ 1

n∑j 6=k

[1− (−1)j+k]1− y2j

(yk − yj)2[Cjgkj − Ckgkk]sign(yk − yj)+

+ Ckgkk

(2yk√1− y2k

− 2 arcsin yk +2yk log(2− 2y2k)√

1− y2k

), k = 1, n.

(2.38)

Impunand ca ecuatia (2.35) sa fie satisfacuta ın punctele yk, k = 1, n, obtinem un sistem

liniar de n ecuatii si n necunoscute (Ck)k=1,n.

21

Page 23: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Capitolul 3

Aripa portanta nonplanara

Capitolul al treilea este dedicat studiului aripii nonplanare. Pornind de la ecuatia

suprafetei portante nonplanare si folosind metoda introdusa de Kida si Miyai [11], am de-

terminat ecuatia liniei portante si am dat interpretarea ın contextul modelului lui Prandtl.

Pentru aceasta ecuatie am determinat doua solutii exacte, pentru aripa avand o distributie

pseudoeliptica a lungimii corzii si pentru aripa inelara circulara. Pentru fiecare din aceste

situatii, solutia asimptotica este aproximarea liniara a solutiei exacte.

Am propus doua metode numerice de rezolvare a ecuatiei liniei portante nonplanare pen-

tru cazul ın care aripa este deschisa si o metoda pentru aripa inelara eliptica. In una din

aceste metode am adus ecuatia la o forma pentru care exista deja ın literatura o teorema de

convergenta a metodei numerice. Testarea acestor metode s-a facut pe solutiile exacte ale

problemei. Simularile numerice dovedesc ca toate cele trei metode au un ordin de convergenta

ridicat. Pentru aripa arcuita circulara avand lungimea corzii constante, am aratat ca solutia

asimptotica diverge la bordurile laterale ale aripii.

Am determinat analitic geometria aripii arcuite circulare avand rezistenta la ınaintare

minima. Cu ajutorul acestei solutii exacte, am validat analitic, pentru clasa aripilor ın forma

de arc de cerc, o conjectura legata de aripa de rezistenta minima.

22

Page 24: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Teorema 4. Ecuatia integrala singulara a suprafetei portante nonplanare este

y′

s(η0)w(ξ0, η0)− z′

s(η0)v(ξ0, η0) =∂zε∂ξ

y′

s −∂yε∂ξ

z′

s, (ξ0, η0) ∈ S, (3.1)

unde componentele v si w ale perturbatiei vitezei pe aripa sunt date de formulele

v (ξ0, η0) =1

∫∫ ∗D

f (ξ, η) z′s (η)

R2s

(1 +

ξ0 − ξR1s

)dξdη−

− 1

∫∫ ∗D

f (ξ, η) z′s (η) [ys (η0)− ys (η)]2

R4s

(1 +

ξ0 − ξR1s

)dξdη+

+1

∫∫ ∗D

f (ξ, η) y′s (η) [ys (η0)− ys (η)] [zs (η0)− zs (η)]

R4s

(1 +

ξ0 − ξR1s

)dξdη−

− 1

∫∫ ∗D

f (ξ, η) {z′s (η) [ys (η0)− ys (η)]− y′s (η) [zs (η0)− zs (η)]}R2s

[ξ0 − ξ] [ys (η0)− ys (η)]

R31s

dξdη,

(3.2)

w (ξ0, η0) = − 1

∫∫ ∗D

f (ξ, η) y′s (η)

R2s

(1 +

ξ0 − ξR1s

)dξdη+

+1

∫∫ ∗D

f (ξ, η) y′s (η) [zs (η0)− zs (η)]2

R4s

(1 +

ξ0 − ξR1s

)dξdη−

− 1

∫∫ ∗D

f (ξ, η) z′s (η) [ys (η0)− ys (η)] [zs (η0)− zs (η)]

R4s

(1 +

ξ0 − ξR1s

)dξdη−

− 1

∫∫ ∗D

f (ξ, η) {z′s (η) [ys (η0)− ys (η)]− y′s (η) [zs (η0)− zs (η)]}R2s

[ξ0 − ξ] [zs (η0)− zs (η)]

R31s

dξdη,

(3.3)

iar

R2s = [ys (η0)− ys (η)]2 + [zs (η0)− zs (η)]2 , (3.4)

R1s =

√(ξ0 − ξ)2 + [ys (η0)− ys (η)]2 + [zs (η0)− zs (η)]2. (3.5)

23

Page 25: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Linia portanta pentru aripa nonplanara

Figura 3.1

Alegem ca lungime caracteristica semianvergura aripii. Notam cu ε maximul lungimii

semicorzii profilului si facem urmatoarele ipoteze:

i) ε� 1 (aripa de anvergura mare);

ii) |ξ| ≤ ε pentru orice punct M(ξ, η) ∈ D (aripa asimilata cu arcul de curba C).(3.6)

Descriem domeniul D si ecuatiile bordului de atac respectiv fuga prin

D = {(ξ, η)| x−(η) ≤ ξ ≤ x+(η), η ∈ [−θ, θ]} ,x−(η) = εxa(η), x+(η) = εxf (η), η ∈ [−θ, θ] ,

(3.7)

unde functiile η 7→ xa(η) si η 7→ xf (η) sunt de ordinul O(1), continue pe [−θ, θ] si de clasa

C1(−θ, θ). Introducem circulatia vitezei ın jurul profilului obtinut prin intersectia aripii cu

planul η = η0,

Γ(η0) =

∫ x+(η0)

x−(η0)

f(ξ, η)dξ. (3.8)

24

Page 26: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Daca C este un arc deschis de curba, la capetele aripii diferenta ıntre valorile presiunii pe

extrados si intrados dispare si prin urmare

Γ(θ) = Γ(−θ) = 0. (3.9)

Daca C este o curba ınchisa, continuitatea circulatiei de-a lungul liniei portante impune ca

Γ(θ) = Γ(−θ). (3.10)

Admitem ca functia Γ este derivabila pe intervalul (-1,1), iar derivata sa este functie Holder

pe intervalul [-1,1].

Propozitia 6. In orice punct M0 ∈ C, componenta normala vitezei pe aripa este:

vn(η0) = limM 7→M0

v(M) · n0 =1

4π√y′2s (η0) + z′2s (η0)

∗θ∫−θ

Γ′ (η) k1(η0, η)dη, (3.11)

unde

k1(η0, η) =[ys(η0)− ys(η)]y′s(η0) + [zs(η0)− zs(η)]z′s(η0)

[ys(η0)− ys(η)]2 + [zs(η0)− zs(η)]2. (3.12)

Teorema 5. Ecuatia integrala a liniei portante nonplanare este

a(η0)Γ (η0)−c(η0)

2

∗θ∫−θ

Γ (η)K(η0, η)

R2s

dη = j(η0), −θ < η0 < θ, (3.13)

unde functia necunoscuta este Γ , iar functiile a, c si j sunt definite de formulele

c(η0) =x+(η0)− x−(η0)

2,

a(η0) =√y′2s (η0) + z′2s (η0),

j(η0) = 2

∫ x+(η0)

x−(η0)

h(ξ0, η0)

√ξ0 − x−(η0)

x+(η0)− ξ0dξ0,

h(ξ0, η0) =∂zε∂ξ

y′

s −∂yε∂ξ

z′

s,

(3.14)

25

Page 27: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

iar nucleul K este o functie regulata ın patratul [−θ, θ]2 si se exprima prin:

K(η0, η) = y′s(η)y′s(η0) + z′s(η)z′s(η0)−

− 2y′s(η)y′s(η0) [zs(η)− zs(η0)]2 + 2z′s(η)z′s(η0) [ys(η)− ys(η0)]2

R2s

+

+2 [ys (η)− ys (η0)] [zs (η)− zs (η0)] [y′s (η0) z

′s (η) + y′s (η) z′s (η0)]

R2s

, η 6= η0

K(η0, η0) = y′2s (η0) + z

′2s (η0) .

(3.15)

DeoareceK

R2s

=∂k1∂η

, urmeaza ca

Propozitia 7. O alta forma a ecuatiei liniei portante este

a(η0)Γ (η0) +c(η0)

2

∗θ∫−θ

Γ′ (η) k1(η0, η)dη = j(η0), −θ < η0 < θ. (3.16)

Solutia asimptotica a ecuatiei liniei portante nonplanare

Introducem noile variabile x1 si ξ1 prin ξ = εξ1, ξ0 = εξ01.

Facem conventia de notare ca relativ la noile variabile x1 si ξ1, functiile sa le notam cu

simbolul tilda deasupra.

Teorema 6. Daca circulatia Γ are un comportament asimptotic de forma Γ = Γ0 + εΓ1 +

O (ε2lnε) , atunci

Γ0 =1√

y′2s (η0) + z′2s (η0)j(η0),

Γ1 =c(η0)

2√y′2s (η0) + z′2s (η0)

∗θ∫−θ

Γ0 (η)K(η0, η)

R2s

dη.

(3.17)

Solutii exacte ale ecuatiei liniei portante nonplanare

i) Aripa pseudoeliptica

In cazul unei aripi pentru care imaginea curbei C este un arc de cerc de raza, dimensionala,a

26

Page 28: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

si deschidere 2θ, 0 < θ < π, daca functiile c si j care dau geometria aripii au forma

c (p)(1 + p2

)= c0

√1− p2

b20, j (p) = c0α0

√1− p2

b20, p = tg

η02, b0 = tg

2

), (3.18)

c0 si α0 fiind constante reale, ecuatia liniei portante nonplanare are solutia de forma

Γexact (p) = λ

√1− p2

b20, p = tg

η02, (3.19)

unde

λ =c0α0

1 +πc04b0

. (3.20)

Solutia asimptotica este:

Γasimptotic (p) = c0Γ0(η0) + c20Γ1(η0) = α0c0

√1− p2

b20

(1− πc0

4b0

), p = tg

η02. (3.21)

Coeficientul de eficienta aerodinamica este

Cef =1

cos4θ

4

> 1, daca 0 < θ ≤ π/2,

Cef =

(4 cos

θ

2tgθ

4

)2

, daca π/2 < θ < π.

(3.22)

Valoarea maxima coeficientului de eficienta este 1.44272 pentru o semideschidere de

θ = 1.80911 .

ii) Aripa inelara circulara

Pentru aripa inelara circulara de incidenta constanta α pentru care imaginea curbei C este

un cerc de raza (dimensionala) a iar lungimea (dimensionala) a semicorzii este constanta ac0.

Alegem a, raza cercului, ca lungime caracteristica. Atunci ecuatia liniei portante nonplanare

se scrie:

Γ(η0)−c08

∫ ∗π−π

Γ(η)

sin2 η − η02

dη = −2πc0α cos η0, −π ≤ η0 ≤ π. (3.23)

27

Page 29: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Solutia exacta ın cazul aripii inelare circulare este:

Γexact(η) = − 2πc0α

1 +πc02

cos η. (3.24)

Solutia asimptotica este:

Γasimptotic(η) = c0Γ(η) = −2πc0α cos η(

1− π

2c0

). (3.25)

Coeficientul de eficienta este Cef = 2.

Solutii numerice

i) Aripa redusa la un arc deschis de curba

Admitem ca C este un arc deschis de curba descris de ecuatiile parametrice:

C : y = y(s) z = z(s) − 1 ≤ s ≤ 1 (3.26)

Facem notatiile

a(s) =√y′2(s) + z′2(s), b(s) =

c(s)

2. (3.27)

Vom da cate o metoda numerica de rezolvare pentru fiecare din ecuatiile (3.13) si (3.16).

Metoda ıntai

Se aduce ecuatia (3.16) la forma

− 1

π

∗1∫−1

Γ′(σ)

σ − sdσ + A(s)Γ(s) +

1

π

1∫−1

H(s, σ)Γ(σ)dσ = J(s), −1 < s < 1, (3.28)

pentru care ın Prossdorf [14], pag.339-341, este data o metoda numerica de rezolvare, ınsotita

de o teorema de convergenta. Pentru discretizarea acestei ecuatii, urmand ideile din [14],

cautam solutii de forma

Γ(s) =√

1− s2C(s).

Aproximam functia C cu polinomul de interpolare Lagrange de gradul n cu nodurile sjnradacinile polinomului Cebasev de speta a doua Un. Astfel, se definesc

Γn(s) =√

1− s2n∑j=1

CjUn(s)

(s− sjn)U ′n(sjn), Cj = C(sj), (3.29)

28

Page 30: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Un(s) =sin(n+ 1)t

sin t, s = cos t, tjn =

n+ 1, sjn = cos tjn. (3.30)

Impunem ca ecuatia 3.16 sa fie satisfacuta ın nodurile Cebasev skn. Cu notatiile

bkj =(1− s2kn)H(skn, skn)

n+ 1, Ak = A(skn) Jk = J(skn)

reducem rezolvarea aproximativa ecuatiei (3.28) la rezolvarea unui sistem patratic de n

ecuatii cu n necunoscute.

AkCk +n∑j=1

(akj + bkj)Cj = Jk, k = 1, n (3.31)

Se introduce spatiul

X0 = {u ∈ C[−1, 1] ∩ C1(−1, 1)|u′ ∈ L2µ, u(−1) = u(1) = 0}

echipat cu norma ‖u‖0 = ‖u‖∞ + ‖u′‖µ si se noteaza cu W 21,µ completatul sau. Atunci are

loc urmatorul rezultat.

Teorema 7. [14] Daca ecuatia omogena asociata ecuatiei(3.28) are doar solutia nula, functiile

A,H si J sunt continue atunci sistemul liniar 3.31 are solutie unica pentru n suficient de

mare iar (Γn)n converge ın norma la solutia ecuatiei (3.28).

Metoda a doua

Pentru rezolvarea ecuatiei

a(s)Γ(s)− b(s)∫ ∗1−1

Γ(σ)N(s, σ)

(s− σ)2dσ = j(s),−1 < s < 1, (3.32)

cu conditiile Γ(−1) = Γ(1) = 0, cautam solutii de forma

Γ(s) =√

1− s2 Γ(s). (3.33)

Impunem ca ecuatia sa fie satisfacuta ın punctele sl = coslπ

n+ 1, l = 1, n, notam

Γl = Γ(sl), al = a(sl), bl = b(sl), jl = j(sl), Nlj = N(sl, sj),

29

Page 31: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

si obtinem sistemul liniar de n ecuatii cu n necunoscute {Γl}l=1,n :

al

√1− s2l Γl−

blπ

n+ 1

n∑j 6=l

[1− (−1)j+l

] 1− s2j(sj − sl)2

NljΓj−π(n+ 1)

2NllΓl = jl, l = 1, n. (3.34)

ii) Aripa inelara

Atunci cand imaginea curbei C este un cerc sau o elipsa, nucleul ecuatiei liniei portante

nonplanare se scrie sub forma

K

R2s

=

(K

R2s

)s

+

(K

R2s

)r

, (3.35)

unde (K

R2s

)s

=1

4 sin2 s− t2

(3.36)

reprezinta partea singulara a nucleului iar

(K

R2s

)r

reprezinta partea regulata a nucleului.

Ecuatia liniei portante nonplanare (3.13) se scrie atunci:

a(t)Γ(t)− c(t)

8

∫ ∗π−π

Γ(s)

sin2 s− t2

ds− c(t)

2

∫ π

−πΓ(s)

(K

R2s

)r

ds = j(t), −π ≤ t ≤ π, (3.37)

cu conditia Γ(−π) = Γ(π). Pentru a discretiza integrala singulara, folosim formula de

cuadratura:∫ ∗π−π

f(t)

sin2 t− tk2

dt ' 4π

2n+ 1

n∑i=−ni 6=k

f(ti)w′i(tk)−

4πn(n+ 1)

2n+ 1f(tk), k = −n, n (3.38)

Impunem ca ecuatia sa fie satisfacuta ın punctele tk =2k

2n+ 1π, k = −n, n, folosim formulele

de cuadratura (3.38) pentru integrala singulara, facem notatiile

Γk = Γ(tk), ak = a(tk), ck = c(tk), jk = j(tk), Rki =

(K

R2s

)r

(tk, ti)

30

Page 32: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

si obtinem sistemul liniar patratic de ordinul 2n+ 1 cu necunoscutele {Γk}k=−n,n :

akΓk−πck

2(2n+ 1)

n∑i=−ni 6=k

w′i(tk)Γi−πck

2n+ 1

n∑i=−n

RkiΓi+πn(n+ 1)ck2(2n+ 1)

Γk = jk, k = −n, n (3.39)

Aripa de rezistenta minima

Problema aripii de rezistenta minima poate fi formulata astfel:

Dintre toate aripile care au aceeasi portanta, sa se determine distributia optima Γopt a

circulatiei vitezei si geometria aripii pentru care rezistenta la ınaintare este minima.

Teorema 8. (Demasi et al. [5]) Atunci cand curba C este simetrica fata de axa OZ,

distributia optima circulatiei vitezei verifica ecuatia

1

∫ θ

−θΓ(η)k1(η, η0)dη = −λys(η0), η0 ∈ (−θ, θ). (3.40)

Ecuatia (3.40) trebuie integrata cu conditiile

Γ(θ) = Γ(−θ) = 0, daca C este un arc deschis de curba,

Γ(θ) = Γ(−θ), daca C este o curba ınchisa.(3.41)

In literatura (Demasi[4]) se cunosc repartitiile optime ale circulatiei viteze pentru aripa

inelara circulara si aripa inelara eliptica de semiaxe aw si bw (atat cazul aw > bw cat si aw < bw). In toate aceste trei situatii distributia optima circulatiei este de forma Γopt(η) = k cos η, k

fiind o constanta reala.

Teorema 9. i) Repartitia optima a circulatiei vitezei ın cazul aripii ın forma de arc de cerc

de deschidere 2θ este

Γopt(p) = κ

√b2 − p2

1 + p2, p = tg

η

2, b = tg

θ

2(3.42)

unde κ este o constanta reala.

31

Page 33: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

ii) In acest caz coeficientul de eficienta aerodinamica este

Cef (θ) = 0.5tg2(θ

2

)+ 1, daca 0 < θ < π/2,

Cef (θ) = 2

[1− cos4

2

)], daca π/2 < θ < π.

(3.43)

Conjectura (Demasi) Daca C este o curba inchisa, Cε este curba obtinuta din C prin

indepartarea unei mici portiuni, iar Γ, respectiv Γε sunt distributiile optime ale circulatiei

vitezei pentru care se obtine aripa de rezistenta minima, atunci limε→0 Γε = Γ + constanta.

Demasi[5] justifica numeric conjectura pentru mai multe tipuri de aripa. Vom valida analitic

conjectura pentru cazul aripii nonplanare ın forma de arc de cerc.

Teorema 10. Fie o aripa ın forma de arc de cerc, cu lungimea semicorzii o constanta c0,

descrisa ın planul Y OZ de un arc de cerc de deschidere 2θ. Fie Γopt(θ) distributia optima a

circulatiei la un coeficient de portanta CL dat. Atunci limθ→π

Γopt(θ) este o distributie optima

a circulatiei vitezei pentru aripa inelara circulara avand aceelasi coeficient de portanta.

32

Page 34: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Concluzii

Rezultatele din aceasta teza pot fi sintetizate astfel:

• Metoda solutiilor fundamentale s-a dovedit, atat ın cazul liniei portante curbe cat si

ın cazul liniei portante nonplanare, un instrument matematic de studiu util pentru deter-

minarea miscarii perturbate. In ambele situatii, metoda solutiilor fundamentale conduce la

aceleasi rezultate ca modelul lui Prandtl. Zona de curgere nepotentiala din aval de aripa,

nu mai este postulata, ca ın modelele clasice, ea este o consecinta a ecuatiilor de miscare.

• In cazul aripii portante curbe, relatia geometrica ıntre unghiul geometric de atac,

unghiul efectiv de atac si deflectia curentului de aer, care conduce, ın modelul lui Prandtl,

la ecuatia liniei portante, trebuie reinterpretata.

• Modelul clasic al lui Prandtl pentru linia portanta dreapta poate fi adaptat la cazul

aripii nonplanare de anvergura mare, cu doua ajustari:

- aripa se reduce la un segment curb, situat ıntr-un plan perpendicular pe directia curentului

uniform;

- zona din spatele aripii, ın care curgerea nu este potentiala, este o suprafata cilindrica,

avand generatoarele paralele cu directia curentului uniform.

• Ecuatia liniei portante nonplanare este o aproximare liniara a ecuatiei suprafetei por-

tante nonplanare, ın raport cu inversul anvergurii aripii. Este o ecuatie cu acelasi tip de

singularitate ca ecuatia lui Prandtl pentru linia portanta dreapta. S-au pus ın evidenta trei

metode numerice de rezolvare, cu o rata de convergenta foarte buna.

• La fel ca ın cazul liniei portante drepte, solutia asimptotica este divergenta la capetele

aripii, in cazul ın care bordurile laterale ale aripii sunt segmente de dreapta.

33

Page 35: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Ca directii de cercetare ulterioare propunem:

• Pentru suprafata planara si pentru linia portanta nonplanara am demonstrat riguros

trecerea la limita de la viteza ın fluid la viteza pe aripa. Pentru suprafata portanta nonpla-

nara si pentru linia portanta curba, problema a ramas nerezolvata.

•Metoda numerica, descrisa ın Capitolul 2 are o rata de convergenta ridicata. Incercarea

de extindere a metodei la aripi simetrice fata de axa OX, de forma arbitrara, conduce la o

schema numerica instabila. De asemenea, se impune elaborarea unei metode numerice efi-

ciente de rezolvare a ecuatiei suprafetei portante nonplanare.

• Studiul efectelor de sol si de tunel de vant ın teoria liniei portante nonplanare.

• Studiul miscarii nestationare ın jurul aripii portante nonplanare.

• O noua abordare a ecuatiei liniei portante nonplanare, prin reducerea acesteia la o

ecuatie Fredholm, folosind un procedeu de regularizare de tip Vekua-Carleman, a operatorilor

integrali singulari.

34

Page 36: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

Bibliografie selectiva

[1] H.Ashley, M.Landahl, Aerodynamics of wings and bodies, Dover Publications (1985)

[2] S.M.Belotserkovsky, I.K.Lifanov, Method of discrete vortices, CRC Press (1993)

[3] C.Cone, The theory of induced lift and minimum induced drag of nonplanar lifting sys-

tems, NASA TR R-139 (1962)

[4] L.Demasi, A theoretical investigation on the conditions of minimum induced drag of

closed wing systems and C-Wings , Journal of Aircraft, 44(1), 81-99 (2007)

[5] L.Demasi, A.Dipace, G.Monegato, R.Cavallaro, Invariant formulation for the minimum

induced drag conditions of nonplanar wing systems, AIAA Journal, 52, 2223-2240 (2014)

[6] L.Dragos, Method of fundamental solutions. A novel theory of lifting surface in subsonic

flow, Arch. of Mechanics, 35(5-6), 579-590 (1983)

[7] L.Dragos, Metode matematice ın aerodinamica, Editura Academiei Romane

(2000),Kluwer Academic Publishers(2003)

[8] J.L.Guermond, Generalized lifting-line theory for curved and swept wings, Journal of

Fluid Mechanics, 211, 497-513 (1990)

[9] D.Homentcovschi, Aerodinamique stationaire linearisee I.(Subsonique), Archives of Me-

chanics, 27, 325-341 (1975).

[10] G.Iosilevskii, Asymptotic theory of high-aspect-ratio arched wings in steady incompress-

ible flow, Journal of Fluid Dynamics, 303, 367-377 (1995)

[11] T.Kida, Y.Miyai, An alternative treatment of lifting-line theory as a perturbation prob-

lem, ZAMP, 29(4), 591-607 (1978)

[12] H.Multhopp, Methods for calculating the lift distribution of wings(subsonic lifting sur-

face theory), ARC R&M, No2884 (1955)

35

Page 37: METODE ASIMPTOTICE S˘I NUMERICE ^IN STUDIUL ARIPII …fmi.unibuc.ro/ro/pdf/2016/doctorat/rezumatSTOICA.pdf · 2016-08-30 · dec^at aripile planare, cu alte cuvinte la anverguri

[13] A.Plotkin, J.Katz, Low-speed aerodynamics, McGraw-Hill (1991)

[14] S.Prossdorf, B.Silbermann, Numerical analysis for integral and related operator equa-

tions, Akademie Verlag 1990

[15] A.Stoica, A Numerical solution for the lifting surface integral equation, INCAS Bul-

letin,6(1), 159-164 (2014)

[16] A.Stoica, A rigorous derivation of certain equations arising in the lifting wing theory,

Bulletin Mathematique de la Societe des Sciences Mathematiques de Roumanie (trimis

spre publicare)

[17] A.Stoica, An exact solution for the minimum induced drag of an circular arched wing

and the validation of a conjecture, Journal of Engineering Mathematics(trimis spre pub-

licare)

[18] A.Stoica, A singular integral equation for curved wings of high aspect ratio, Proc. Ro.

Acad. Series A (trimis spre publicare)

[19] E.O.Tuck, Some accurate solutions of the lifting surface integral equation,

J.Austral.Math.Soc.Ser.B, 35, 127-144 (1993)

[20] J.Weissinger, The lift distribution of swept-back wings, ARC R&M, No1120 (1947)

36