Moti di fluido ideale irrotazionale:
metodologie di soluzione
Cesare Corrado
Marzo 2009
In copertina: linee di corrente e campo di velocità per un profilo
NACA0012 risolte con la teoria del potenziale
Indice
Introduzione
1
1 Alcuni richiami sulle equazioni dei fluidi
3
1.1
Le equazioni dei fluidi: ipotesi semplificative ed equazioni del
potenziale cinetico . . . . . . . . . . . . . . . . . . . . . . . .
3
1.2
La funzione di Corrente ψ . . . . . . . . . . . . . . . . . . . . 11
1.3
Connessione del dominio . . . . . . . . . . . . . . . . . . . . . 12
1.4
Condizioni di unicità per il potenziale . . . . . . . . . . . . . . 13
1.5
Condizioni al contorno per la scia . . . . . . . . . . . . . . . . 15
2 I metodi a singolarità virtuali
2.1
17
Il metodo della funzione di Green . . . . . . . . . . . . . . . . 20
2.1.1
Applicazione della formula di Green al problema aerodinamico
. . . . . . . . . . . . . . . . . . . . . . . . . 22
3 I metodi a Pannelli
25
3.1
Il metodo di Hess e Smith . . . . . . . . . . . . . . . . . . . . 27
3.2
Il metodo di Morino nel caso stazionario . . . . . . . . . . . . 30
3.2.1
Alcune osservazioni sulle matrici D e Z . . . . . . . . . 33
4 Altri metodi numerici per l’aerodinamica
35
4.1
La linea portante-schema di Weissinger . . . . . . . . . . . . . 36
4.2
La superficie portante . . . . . . . . . . . . . . . . . . . . . . . 37
4.2.1
Il Vortex-Lattice Method . . . . . . . . . . . . . . . . . 37
iii
INDICE
4.2.2
iv
La Piston Theory . . . . . . . . . . . . . . . . . . . . . 38
Introduzione
La risoluzione delle equazioni che regolano il campo di moto di un fluido è
un problema non risolvibile per via analitica e praticamente impossibile da
risolvere per via numerica a causa dell’elevata separazione delle scale spaziotemporali che regolano il fenomeno; inoltre, per le applicazioni pratiche la
soluzione conterrebbe un maggior numero di informazioni di quelle necessarie. Di qui la necessità di creare modelli semplificati in grado di fornire
un’approssimazione della soluzione sufficientemente accurata per i fini pratici e non troppo onerosa dal punto di vista computazionale.
Per quanto concerne le applicazioni aeronautiche, l’informazione necessaria
per un’analisi preliminare riguarda la distribuzione del carico aerodinamico
sulle superfici bagnate dal fluido, carico fornito in gran parte dalla portanza,
essendo le forze resistive inferiori di diversi ordini di grandezza: di qui la
possibilità di assumere modelli a fluido ideale. Se il campo di moto del fluido
ideale è anche di tipo isoentropico ne conseguirà che sarà anche irrotazionale,1
permettendo di introdurre una grandezza detta potenziale cinetico che ridurrà il problema ad un’equazione scalare. Un’ulteriore semplificazione potrà
essere ottenuta se si riuscirà ad esprimere le grandezze aerodinamiche significative senza dover ricorrere alla soluzione di tutto il campo di moto, ma
risolvendo un problema sulle sole superfici aerodinamiche.
Per la risoluzione numerica esistono vari metodi, tipicamente:
1
In realtà per un moto in cui siano presenti onde d’urto deboli, ovvero onde d’urto
approssimabili con un raggio di curvatura infinito, è possibile ancora applicare la teoria
irrotazionale
1
Introduzione
• Elementi finiti
• Volumi finiti
• Differenze finite
• Elementi di contorno
Nella trattazione seguente verrà esposta la sola categoria dei metodi ad
elementi di contorno, poiché questi trovano ampio impiego nelle applicazioni
aeronautiche; di questi inoltre verrà trattato il solo caso stazionario subsonico, in quanto è l’unico che non richiede una risoluzione nel dominio della
trasformata di Fourier; verrà tuttavia introdotto anche un metodo per regimi supersonici, sotto le cui condizioni di applicabilità è possibile trovare una
soluzione senza dover ricorrere al dominio delle frequenze.
Queste note sono organizzate nel modo seguente:
• Nel primo capitolo vengono fatti alcuni richiami sulla teoria dei fluidi,
in particolar modo per quanto riguarda l’ottenimento dell’equazione del
potenziale e le ipotesi di validità;
• Nel secondo capitolo vengono introdotti i metodi a singolarità virtuali
e la funzione di Green,in particlare per il caso Aerodinamico;
• Nel terzo capitolo vengono illustrati i metodi a pannelli tipicamente più
utilizzati: il metodo di Hess e Smith ed il metodo di Morino;
• Nel quarto capitolo vengono illustrate tecniche di soluzione per profili
di spessore trascurabile: in particolar modo metodi detti “a superficie
portante”. Viene inoltre illustrato un metodo applicabile al caso supersonico, sotto l’ipotesi che la corrente asintotica abbia un numero di
Mach superiore a 2.
2
Capitolo 1
Alcuni richiami sulle equazioni
dei fluidi
1.1
Le equazioni dei fluidi: ipotesi semplificative ed equazioni del potenziale cinetico
Per descrivere il moto di un fluido comprimibile occorre risolvere un sistema
di equazioni algebriche e differenziali alle derivate parziali che ne traduce
le leggi di conservazione; nel caso più generale di un fluido omogeneo, le
equazioni in gioco tradurranno le seguenti leggi:
• Conservazione della massa
• Conservazione della quantità di moto
• Conservazione dell’energia
• Equazione termica di stato (p = p(ρ, T ))
• Equazione calorica di stato (e = e(ρ, T ) )
• Relazione costitutiva (σ = σ(p, ∇v))
3
1. Alcuni richiami sulle equazioni dei fluidi
Considerando anche l’equazione del momento della quantità di moto, si
ottiene che il tensore degli sforzi fluidi è simmetrico e quindi, nel caso tridimensionale, sarà costituito da 6 incognite, legate al campo di moto ed alla
pressione attraverso le relazioni costitutive. Si tratta dunque di risolvere un
sistema complesso nelle incognite p, v, T, ρ, e, in tutte le sue scale spaziali e
temporali.
Per gli scopi aerodinamici, le informazioni fornite dalla risoluzione di un
tale sistema sono sovrabbondanti: è allora lecito domandarsi se un modello
semplificato fornisca le informazioni necessarie con un costo computazionale
nettamente inferiore. Si introducono le seguenti ipotesi semplificative:
1. flusso inviscido (µ = 0)
Il campo di moto è descritto dalle equazioni di Eulero:1
∂ρ
+ ∇ · (ρv) = 0
∂t
∂ρv
+ ∇ · (ρv ⊗ v) + ∇p = 0
∂t
p
= cost
ργ
(1.1)
(1.2)
(1.3)
avendo indicato con ⊗ il prodotto tensoriale fra due vettori; l’ultima
relazione deriva dal fatto che il flusso è isoentropico (s = cost).
2. Flusso irrotazionale (∇ × v = 0)
Il campo di moto è descritto dalle equazioni di Eulero semplificate, che
assumono la seguente forma:
∂ρ
+ ∇ · (ρv) = 0
∂t
2
|v|
∂v
ρ
+ ρ∇
+ ∇p = 0
∂t
2
p
= cost
ργ
1
(1.4)
(1.5)
(1.6)
Queste vengono scritte in forma conservativa; esplicitando le derivate ed utilizzando
l’equazione di continuità è possibile trovare anche la forma non conservativa, tuttavia molti
metodi numerici si basano sulla prima
4
1. Alcuni richiami sulle equazioni dei fluidi
Occorre ora fornire le condizioni al contorno; considerando un problema
esterno e descrivendo la superficie del corpo immerso nel fluido con una
superficie S, la condizione al contorno di non penetrazione si traduce
nella:
v·n=0
su S
(1.7)
Ricordando che ogni campo vettoriale irrotazionale in un dominio semplicemente connesso2 può essere espresso come il gradiente di una funzione potenziale, la (1.5) può essere riscritta nella forma:
∂φ |v|2
1
∇
+
+ ∇p = 0
∂t
2
ρ
(1.8)
Si consideri ora una curva proveniente dall’infinito e terminante in un
punto P del campo di moto e si vada ad integrare la (1.8) lungo tale
cammino, ciò porta al teorema di Bernoulli per flussi irrotazionali:
Z P
Z P ∂φ |v|2
1
d
+
dp = 0
+
∂t
2
−∞
−∞ ρ
Z P
1
∂φ |v|2 − |v∞ |2
+
+
dp = 0
(1.9)
=⇒
∂t
2
−∞ ρ
avendo assunto un flusso stazionario a monte. Utilizzando ora l’equazione dell’energia:
p
p∞
p γ
= γ ⇒ ργ =
ρ
γ
ρ
ρ∞
p∞ ∞
(1.10)
e ricordando la definizione della velocità del suono:
r
p
c= γ
ρ
2
(1.11)
Si fa osservare che per il problema in questione il dominio non è in genere semplice-
mente connesso; come si vedrà in seguito, questa ipotesi sarà soddisfatta per un corpo
aerodinamico introducendo la scia.
5
1. Alcuni richiami sulle equazioni dei fluidi
Si ottiene:
p∞
p∞
p = ργ γ → dp = γ γ ργ−1 dρ
ρ∞
ρ∞
P
Z P
Z P
p∞ γ−2
1
γ p∞ γ−1
γ γ ρ dρ =
dp =
ρ
γ − 1 ργ∞
−∞ ρ∞
−∞ ρ
−∞
γ−1
γ1
γ
p
p
⇒ ργ−1 = ργ−1
ρ = ρ∞
∞
p∞
p∞

P
1
Z P
γ
γ−1 p∞
1
γ

dp = 
p γ
ρ
γ
−
1
ρ
∞
−∞
−∞
1
γ
1
γ
p∞
p
=
ρ∞
ρ
Z P
1
c2 − c2∞
dp =
=⇒
γ−1
−∞ ρ
da cui:
c2
|v∞ |2
c2
∂φ |v|2
+
+
=
+ ∞
∂t
2
γ−1
2
γ−1
(1.12)
che rappresenta il teorema di Bernoulli per un flusso:
• comprimibile
• isoentropico
• irrotazionale
• non stazionario
Considerando il coefficiente di pressione :
cp =
p − p∞
1
ρ |v |2
2 ∞ ∞
(1.13)
Dalle (1.10), (1.11) e (1.12) si può scrivere:
p
γρc2
=
=
p∞
γρ∞ c2∞
6
p
p∞
γ1
c2
→
c2∞
p
p∞
γ
=
p
p∞
c2
c2∞
γ
1. Alcuni richiami sulle equazioni dei fluidi
p
⇒
p∞
c2
c2∞
c2
c2∞
γ
γ−1
γ−1
=1− 2
c∞
ρ
=
ρ∞
c2
c2∞
1
γ−1
(1.14)
∂φ 1
+ (|v|2 − |v∞ |2 )
∂t
2
(1.15)
Si ottiene:
cp =
p
−1
p∞
1
2
γM∞
2
(1.16)
Si perviene così al problema del potenziale completo, più conveniente
rispetto al problema di Eulero in quanto si sono ridotte le incognite del
campo di veloictà da tre ad una (scalare):
∂ρ
+ ∇ · (ρ∇φ) = 0
∂t
1
2 γ−1
c
ρ(φ) = ρ∞ 2
c∞
2
c
γ − 1 ∂φ 1
2
2
+ (|v| − |v∞ | )
=1− 2
c2∞
c∞
∂t
2
∂φ
=0
∂n
(1.17)
(1.18)
(1.19)
su S
(1.20)
Osservazione 1 Il modello del potenziale completo non può essere applicato in presenza di onde d’urto, in quanto il flusso che si viene a creare non è
più isoentropico e di conseguenza cade l’ipotesi di irrotazionalità;3 tuttavia,
nel caso in cui si abbiano onde d’urto deboli (e di conseguenza rettilinee) il
modello può ancora essere applicato, poiché il salto entropico da luogo ad un
rotore della velocità praticamente nullo.
3. Approssimazione delle piccole perturbazioni
In questo caso il potenziale è costituito dalla sovrapposizione lineare del
3
Come si vede dall’equazione di Crocco
7
1. Alcuni richiami sulle equazioni dei fluidi
potenziale corrispondente al flusso indisturbato e di una perturbazione
(teoria linearizzata):
φ = V∞ (x + ϕ)
(1.21)
Si passa dunque a linearizzare le grandezze in gioco mediante lo sviluppo:
c = c∞ + ∆c
p = p∞ + ∆p;
Per la velocità, poichè i termini quadratici sono trascurabili:
∂ϕ
∂ϕ
∂ϕ
v = ∇φ =
1+
i+
j+
k V∞
∂x
∂y
∂z
∂ϕ
2
v · v ' V∞ 1 + 2
∂x
(1.22)
(1.23)
(1.24)
Trascurando i termini quadratici, dalla (1.19) si ottiene:
c2
c2∞ + 2c∞ ∆c
'
c2∞
c2∞
γ−1
∂ϕ 1
∂ϕ
2
2
+
− V∞
=1− 2
V∞
V∞ + 2V∞
c∞
∂t
2
∂x
γ−1
Dϕ
= 1 − 2 V∞
c∞
Dt
(1.25)
Per la (1.18), sostituendo e sviluppando in serie di Taylor nell’intorno
di c2∞ , arrestandosi al primo ordine:
1
2
ρ∞ + ∆ρ
c∞ + 2c∞ ∆c γ−1
'
ρ∞
c2∞
2
1
2
2−γ
c∞ + 2c∞ ∆c γ−1
1
c∞ + 2c∞ ∆c γ−1 2∆c
'
+
c2∞
γ−1
c2∞
∆c=0
∆c=0 c∞
2∆c
'1+
(1.26)
c∞ (γ − 1)
∆ρ
1
Dϕ
⇒
= − 2 V∞
(1.27)
ρ∞
c∞
Dt
8
1. Alcuni richiami sulle equazioni dei fluidi
dove:
Dϕ
Dt
(1.28)
è la derivata materiale (Lagrangiana) del potenziale linearizzata.
Dall’equazione di continuità si ottiene la linearizzazione:
∂∆ρ
∂ϕ
∂ϕ
∂ϕ
+ ∇ · (ρ∞ + ∆ρ) V∞
1+
i +
j+
k
=0
∂t
∂x
∂y
∂z
∂∆ρ
∂∆ρ
⇒
+ ρ∞ V∞ ∇ · ∇ϕ + v∞
'0
∂t
∂x
(1.29)
Andando a sostituire la (1.27), ricordando la definizione della (1.28) e
dopo alcuni passaggi:
∇ · ∇ϕ =
1 D2 ϕ
c2∞ Dt2
(1.30)
Per ottenere il coefficiente di pressione dalla (1.13), per prima cosa si
linearizza la pressione con uno sviluppo di Taylor arrestato al primo
ordine nell’intorno di c2∞ :
2
γ
p
p∞ + ∆p
c∞ + 2c∞ ∆C γ−1
=
=
p∞
p∞
c2∞
#
"
2
1
c∞ + 2c∞ ∆c γ−1 2c∞
γ
'1+
γ−1
c2∞
c2∞
∆p
2∆c γ
V∞ Dϕ
⇒
=
= −γ 2
p∞
c∞ γ − 1
c∞ Dt
∆c
(1.31)
∆c=0
(1.32)
andando ora a sostituire l’espressione appena trovata nella (1.16) si
ottiene:
cp = −
2 Dϕ
V∞ Dt
(1.33)
9
1. Alcuni richiami sulle equazioni dei fluidi
particolare importanza merita il caso stazionario, per cui la (1.30)
diventa:
2
(1 − M∞
)
∂ 2ϕ ∂ 2ϕ ∂ 2ϕ
+ 2 + 2 =0
∂x2
∂y
∂z
(1.34)
dove si vede che per M∞ < 1 (caso subsonico) l’equazione è di tipo
ellittico, mentre per M∞ > 1 l’equazione diventa di tipo iperbolico.
Per quanto riguarda le condizioni al contorno, ci si metta in un sistema
di riferimento solidale col corpo aerodinamico, di velocità VB nota; si
ottiene:
− VB · n = ∇φ · n
(1.35)
Applicando ora la teoria lineare per il potenziale:
∂ϕ
∂ϕ
∂ϕ
j+
k · n = V∞ i · n + V∞ ∇ϕ · n
i+
∇φ · n = V∞ 1 +
∂x
∂y
∂z
⇒
VtotBN
∂ϕ
(V∞ i + VB ) · n
=−
=−
∂n
V∞
V∞
(1.36)
Avendo indicato con VtotBN la velocità assoluta della corrente in direzione normale al corpo che, per le ipotesi di piccole perturbazioni, indicherà la pendenza locale del corpo aerodinamico rispetto la corrente
assoluta; infatti, osservando la figura 1.1 si vede che si può scrivere:
tan α ' sin α ' α
(1.37)
|VtotB | = |V∞ i + VB | ' V∞
(1.38)
α ' sin α ' −
⇒ α=
10
∂ϕ
∂n
|VtotBN |
V∞
(1.39)
(1.40)
1. Alcuni richiami sulle equazioni dei fluidi
Figura 1.1: Velocità ed incidenza
Per la chiusura del problema si imporrà infine una condizione al contorno di Dirichlet per il potenziale all’infinito, data dalla seguente:
φ∞ = V∞ x ⇒ ϕ∞ = 0
1.2
(1.41)
La funzione di Corrente ψ
Si consideri il caso stazionario ed incomprimibile: l’equazione del potenziale
cinetico si riduce alla risoluzione dell’equazione di Laplace:
∇2 φ = 0
(1.42)
munita delle opportune condizioni al contorno. Sotto queste ipotesi, è possibile riformulare il problema di Laplace per una nuova grandezza detta funzione di corrente ψ, che fisicamente rappresenta una portata attraverso una
superficie, così definita nel caso bidimensionale:
Z
P2
v · ndc
ψ(P2 ) = ψ(P1 ) +
(1.43)
P1
dove P1 e P2 sono due punti che giacciono in una regione connessa del
fluido ed l è il generico percorso aperto che li unisce. Nel caso di un fluido
11
1. Alcuni richiami sulle equazioni dei fluidi
incomprimibile, è possibile dimostrare attraverso il teorema di Gauss che
il valore di ψ dipende solamente dalla posizione dei punti rispetto cui si
calcola, ma non dal percorso che li unisce. Per quanto riguarda la funzione
di corrente, in maniera duale a quanto detto per il potenziale, si avrà che
essa esisterà solo sotto la condizione di moto incomprimibile4 e nel caso di
moto irrotazionale soddisferà l’equazione di Laplace:
∇2 ψ = 0
(1.44)
Si può inoltre dimostrare che la funzione di corrente, il vettore velocità
ed il potenziale sono legati attraverso le seguenti relazioni:
∂φ
∂ψ
=u=
∂x
∂y
∂φ
∂ψ
=v=−
∂y
∂x
(1.45)
Se ne deduce che le linee equipotenziali e le linee in cui ψ = cost costituiscono un sistema di curve ortogonali.
1.3
Connessione del dominio
Nella sezione 1.1 si era accennato al fatto che affinchè il potenziale sia univocamente determinato, oltre all’ipotesi di irrotazionalità, occorre anche ipotizzare che il dominio sia semplicemente connesso. Si introduce allora la
seguente:
Definizione 1 Una regione dello spazio si dice semplicemente connessa (o
connessa con molteplicità uno) quando ogni coppia di cammini che uniscono
due punti del dominio e che rimangono in esso compresi, costituisce una curva
chiusa riducibile
Dalla definizione introdotta ne risulta che, ad esempio nel caso bidimensionale, se la curva in questione racchiude il copro aerodinamico, questa non
4
In realtà è possibile definirne una “modificata ” che può essere utilizzata anche nel
caso comprimibile
12
1. Alcuni richiami sulle equazioni dei fluidi
è riducibile:5 esiste quindi una regione dello spazio non semplicemente connessa ed in essa il potenziale sarà una funzione polidroma (che assume cioè
più di un valore).
Si consideri allora una generica curva chiusa C che racchiuda il copro inmmerso nella corrente fluida e si calcoli la circolazione di v lungo detto cammino;
ciò che si ottiene è la seguente:
I
v · tdc = Γ
(1.46)
C
Questo risultato mostra una prima interessante proprietà: per tutti i
cammini chiusi non riducibili e che circondano una volta la lacuna, la circolazione della velocità è indipendente dal percorso e pari ad un valore costante
Γ, detta costante ciclica: se infatti il cammino non riducibile si avvolgesse
per n, n ∈ Z volte attorno alla lacuna, si avrebbe che:
I
v · tdc = nΓ
(1.47)
C
Il potenziale allora viene ad essere definito come:
φ(P ) = φ∗ (P ) + nΓ
,n ∈ Z
(1.48)
avendo indicato con φ∗ (P ) il valore principale del potenziale. Si è ottenuto un
risultato di notevole interesse: essendo interessati al gradiente del potenziale,
la parte relativa alla polidromicità avrà derivata nulla, essendo Γ costante.
1.4
Condizioni di unicità per il potenziale
Si è visto che affinchè il sia una funzione monodroma (e quindi determinabile
in maiera univoca) è necessario o essere in un dominio semplicemente connesso o avere un valore di circolazione Γ nullo; si è visto inoltre che per poter
5
Nel caso tridimensionale, per un’ala di apertuta finita, il dominio è invece semplice-
mente connesso; un dominio tridimensionale non semplicemente connesso è ad esempio
un dominio attraversato da una retta, ovvero un dominio in cui i “buchi ” non sono di
dimensione finita.
13
1. Alcuni richiami sulle equazioni dei fluidi
determinare il potenziale in maniera univoca nel caso in cui la costante ciclica
non sia nulla occorre conoscere il valore di detta costante.
Si vuole ora applicare il teorema di Gauss nel caso di un dominio non semplicemente connesso e ricavare il potenziale in maniera univoca: il trucco
consiste nel rendere semplicemente connsesso il dominio mediante l’introduzione di un talgio, il cui compito è quello di non consentire percorsi chiusi
che circondino la lacuna (si veda la figura 1.2). Tuttavia, con questa operazione il potenziale presenterà una discontinuità attraverso detta linea: due
punti a cavallo del taglio infatti dal punto di vista topologico stanno agli estremi di una regione semplicemente connessa. Inoltre è possibile dimostrare
Figura 1.2: Taglio del dominio
attraverso il teorema di Gauss che il salto del potenziale a cavallo della scia
è pari proprio al valore della costante ciclica:
φ+ − φ− = ∆φ = Γ
(1.49)
In altre parole: la soluzione è determinata univocamente quando si conosce
il valore della derivata normale del potenziale sul contorno ed il valore della
costante ciclica.
Nelle applicazioni aerodinamiche relative agli aeromobili si è in presenza di
corpi aerodinamici:se l’angolo con cui la corrente investe il corpo non è troppo
14
1. Alcuni richiami sulle equazioni dei fluidi
elevato, la scia che si forma dietro di essi è caratterizzata da una dimensione
longitudinale prevalente rispetto allo spessore, inoltre si forma a partire da
un bordo di uscita aguzzo (il bordo di uscita del profilo). Per i problemi
di interesse allora i fenomeni legati al rotore della velocità saranno confinati
all’interno della scia, che costituirà il taglio del dominio precedentemente
introdotto: taglio di cui si conosce il punto di distacco, e che servirà quindi per descrivere una corrente che contiene anche della vorticità attraverso
l’equazione di Laplace.
1.5
Condizioni al contorno per la scia
Perchè l’equazione di Laplace possa essere risolta occorre fornire anche delle
condizoni al contorno sulla scia, la cui posizione inoltre è incognitae deve essere ricavata essa stessa come parte della soluzione. Per risolvere il problema
si ricorre alle leggi di conservazione della massa e della quantità di moto in
forma integrale.6
Prendendo due punti P1 e P2 a cavallo della scia e facendoli tendere a questa,
dall’equazione di conservazione della massa si ottiene:
∂φ ∂φ =
∂n 1
∂n 2
(1.50)
ovvero: la derivata del potenziale in direzione ortogonale alla scia deve essere continua; nel caso in cui il flusso sia anche stazionario, la discontinuità
costituisce anche una linea di corrente: ne segue allora che
∂φ =0
∂n (1.51)
non potendoci essere flusso di massa in direzione ortogonale ad una linea di
corrente.
Dalla continuità della derivata normale del potenziale e dall’equazione di
6
Il potenziale attraverso la scia è discontinuo e quindi non può essere utilizzata la forma
differenziale.
15
1. Alcuni richiami sulle equazioni dei fluidi
conservazione della quantità di moto inoltre si ritrova che anche la pressione
deve essere continua a cavallo della scia:
p1 = p 2
(1.52)
Applicando il teorema di Bernoulli dalla (1.52) è possibile ottenere una
condizione al contorno sul potenziale; indicando con ∆φ = φ1 − φ2 il salto di
potenziale a cavallo della scia si ottiene:
∂∆φ
+ Vs · ∇∆φ = 0
∂t
dove:
∇φ1 + ∇φ2
Vs =
2
(1.53)
è detta velocità sulla scia
(1.54)
In altre parole, si ha una condizione di trasporto per la discontinuità. Nel
caso particolare di moto stazionario questa condizione mostra che il salto di
potenziale ∆φ si mantiene costante sulla scia
16
Capitolo 2
I metodi a singolarità virtuali
Sono metodi utilizzati tipicamente per il caso stazionario incomprimibile e
consistono nel rappresentare il corpo tramite una distribuzione di singolarità, o una combinazione lineare delle stesse1 (sorgenti, pozzi dipoli, etc)
di intensità incognita, funzione però delle condizioni al contorno relative
al corpo aerodinamico. Le soluzioni semplici dell’equazione di Laplace che
determinano dette singolarità sono tipicamente le seguenti:
• La corrente uniforme dove, detta v∞ la velocità asintotica, si ottiene
(si veda rigura 2.1):
φ(x, y) = u∞ x + v∞ y
(2.1)
ψ(x, y) = u∞ y − v∞ x
(2.2)
• La sorgente dove, detta q la portata di massa emessa da una sorgente
posta nell’origine degli assi, r la distanza da essa e θ la relativa posizione
angolare, si ha (si veda figura 2.2):
q
log r
2π
q
ψ(r, θ) =
θ
2π
φ(r, θ) =
1
(2.3)
(2.4)
Avendo un’equazione alle derivate parziali lineare, si può applicare il principio di
sovrapposizione degli effetti.
17
2. I metodi a singolarità virtuali
• Il vortice che costituisce il campo di moto duale a quello della sorgente;
detta Γ l’intensità (circolazione) di un vortice posto nell’origine degli
assi, si ha:
Γ
θ
2π
Γ
ψ(r, θ) = − log r
2π
φ(r, θ) =
(2.5)
(2.6)
• La dipolo o dipolo, consistente nella sovrapposizione di una sorgente
e di un pozzo di eguale intensità q, posti ad una distanza d che viene
fatta tendere a zero in modo tale che il prodotto µ = qd, che costituisce
l’intensità della dipolo, si mantenga finito. Si ottiene (si veda figura
2.3):
µ cos θ
2π r
µ sin θ
ψ(r, θ) =
2π r
φ(r, θ) = −
(2.7)
(2.8)
Figura 2.1: Potenziale e linee di corrente per la corrente uniforme
Considerando un corpo aerodinamico, per quanto visto nella sezione 1.5,
affinchè la velocità del dorso e quella del ventre siano eguali al bordo di uscita,
è necessario introdurre una vorticità sul profilo, attraverso una singolarità di
tipo vortice: questo farà si che il corpo sia di tipo prtante; dovendosi inoltre
conservare la vorticità nel campo di moto e considerando di essere partiti
18
2. I metodi a singolarità virtuali
Figura 2.2: Potenziale e linee di corrente per la sorgente
Figura 2.3: Potenziale e linee di corrente per il dipolo
19
2. I metodi a singolarità virtuali
da uno stato di quiete in cui sicuramente la vorticità è nulla, dovrà esistere
fra il corpo e l’infinito a valle dello stesso una vorticità di intensità uguale
e contraria a quella presente sul corpo: questa vorticità è costituita dalla
scia. Detta ξ la distribuzione del vettore vorticità e Q il generico punto dove
è generata la vorticità, la velocità indotta da tale distribuzione può essere
calcolata mediante la legge di Biot-Savart:
Z Z Z
1
ξ × rP Q
V (P ) =
dV
3
4π
V krP Q k
(2.9)
che costituisce un integrale di convoluzione. Grazie alla linearità dell’operatore di Laplace è poi possibile costruire la soluzione desiderata semplicemente combinando le soluzioni fondamentali ed imponendo il soddisfacimento
delle condizioni al contorno.
Esempio 1 Si voglia ricavare l’andamento della corrente attorno ad un profilo alare portante: si distribuiranno sulla stessa delle sorgenti che forniranno l’andamento dello spessore e delle dipoli per far si che
il corpo risulti di forma chiusa. Le uniche incognite del problema sono costituite dalle intensità delle singolarità virtuali, incognite che saranno determinate dall’imposizione della condizione di non penetrazione
sul contorno del profilo.
2.1
Il metodo della funzione di Green
Si consideri l’equazione di Poisson, generalizzazione dell’equazione di Laplace:
∇2 φ = g(r)
(2.10)
Essendo un’equazione lineare, sarà possibile applicare il principio di sovrapposizione degli effetti ed esprimere la soluzione nel generico punto r0 come
somma di due contributi, uno relativo al termine noto g con condizioni al
contorno omogenee ed uno relativo alle condizioni al contorno con termine
noto nullo:
Z
G(r0 , r)g(r)dV + φ(r0 ) |BC
φ(r0 ) =
(2.11)
V
dove la funzione G, che compare nell’integrale di convoluzione è detta
funzione di Green e costituisce la funzione di influenza che esprime il valore
20
2. I metodi a singolarità virtuali
del potenziale in un punto in funzione della causa che lo ha generato in un
altro punto; tale funzione è per ora indeterminata.
Si considerino ora i seguenti problemi:
∇2 G = δ(r − r0 )
(2.12)
∇2 φ = g
(2.13)
Si moltiplichi la (2.12) per φ, e la si sottragga alla (2.13) moltiplicata per
G; integrando:
Z
V
G(r0 , r)∇2 φ(r) − φ(r)∇2 G(r0 , r) dV
Z
Z
=
G(r0 , r)g(r)dV −
φrδ(r − r0 )dV
V
(2.14)
V
poichè:
φ∇2 G = ∇ · (φ∇G) − ∇φ · ∇G
G∇2 φ = ∇ · (G∇φ) − ∇φ · ∇G
utilizzando il teorema della divergenza si ottiene infine:
Z
Z
G(r0 , r)g(r)dV +
Eφ(r0 ) =
∂V
V
∂G(r0 , r)
∂φ(r)
φ(r)
− G(r0 , r)
∂n
∂n
dS
(2.15)
Dove il primo termine a secondo membro rappresenta l’effetto dovuto al
termine g, mentre il secondo termine a secondo membro rappresenta il contributo legato alle condizioni al contorno del problema. Il termine E che
moltiplica il membro di destra è tipico dei problemi aerodinamici ed assume
0 all’interno del corpo
i seguenti valori: 1 nella zona in cui è presente il fluido
1
2
sul contorno del corpo
Osservazione 2 Nel termine di bordo:
Z ∂G(r0 , r)
∂φ(r)
φ(r)
− G(r0 , r)
dS
∂n
∂n
∂V
(2.16)
21
2. I metodi a singolarità virtuali
si conoscerà rispettivamente o il potenziale o la sua derivata dalle condizioni al
contorno: ne segue che uno dei due termini rimane indefinito. Se il contorno è di
forma semplice e descrivibile analiticamente è possibile selezionare una particolare
funzione di Green G che annulli il contributo incognito.
Considerando, per semlicità, il caso r0 = 0, si ottiene per la funzione di
Green:
ln(r)
2π
1
G=−
4πr
G=
2.1.1
per il caso 2D
(2.17)
per il caso 3D
(2.18)
Applicazione della formula di Green al problema
aerodinamico
La formula di green appena introdotta vale per un punto r0 interno al dominio di integrazione: sembrerebbe allora impossibiler applicare questa tecnica risolutiva al caso aerodinamico, dove il dominio è di dimensione infinita e caratterizzato dal campo di moto esterno al profilo aerodinamico con
condizioni al contorno sia sul corpo che all’infinito.
Figura 2.4: Dominio di integrazione per la funzione di Green (Da ....)
Si consideri il contorno rappresentato in figura 2.4 e vi si applichi la
formula di Green: si può osservare che il contorno l viene percorso due volte in
senso opposto, per cui non fornirà alcun contributo. L’integrale sul contorno
22
2. I metodi a singolarità virtuali
S∞ è non nullo, tuttavia può tendere ad un valore costante. Grazie alla
linearità dell’equazione del potenziale, si può definire una nuova grandezza,
detta potenziale di perturbazione:
ϕ(r0 ) = φ(r0 ) − V∞ · r0
(2.19)
che si annulla all’infinito (e quindi non dà contributo su S∞ ) e che è del tutto
equivalente a calcolare φ Si può così utilizzare la formula di Green, essendosi
ricondotti a soli integrali sul corpo aerodinamico:
I ∂G
∂ϕ
ϕ(r0 ) =
ϕ
−G
dS
∂n
∂n
S1
∂ϕ
= −V∞ · n
∂n
(2.20)
su S1
(2.21)
Alternativamente, considerando che all’infinito l’effetto del corpo è trascurabile per cui:
∂G
∂φ
φ
φ(r0 ) =
−G
dS = V∞ · r0
∂n
∂n
S∞
I
(2.22)
che rappresenta la soluzione esatta della corrente indisturbata, è possibile
esprimere la soluzione sul potenziale come:
I ∂φ
∂G
φ(r0 ) =
−G
dS + V∞ · r0
φ
∂n
∂n
S1
(2.23)
Si è interessati al caso di profili portanti, profili cioè che sono in grado di
generare una forza aerodinamica. Per la definizione di corpo aerodinamico,
la scia sarà di spessore infinitesimo, si staccherà dal bordo di uscita ed in
essa sarà concentrata tutta la vorticità attraverso la quale il corpo genera
portanza. La regione costituita dalla scia si trova dunque all’esterno del
dominio di itegrazione e la formula di Green, valutata lungo il contorno della
scia, darà contributo non nullo; tuttavia, essendo questo contorno percorso
due volte in senso opposto come accadeva prima per il contorno l, ed essendo
il potenziale discontinuo a cavallo, ne seguirà che:
I Z
∂G
∂ϕ
∂G
ϕ(r0 ) =
ϕ
−G
dS +
∆ϕ dS
∂n
∂n
∂n
S1
Sw
(2.24)
23
2. I metodi a singolarità virtuali
Avendo indicato con Sw il contorno costituito dalla scia. Nel caso 2D
stazionario l’espressione del potenziale assume un’espressione più semplice,
essendo ∆φ costante lungo la linea di corrente costituita dalla scia; introducendo la funzione, continua lungo la scia:
ϕ2 (r0 ) = φ(r0 ) − V∞ · r0 −
∆φ
θ
2π
(2.25)
è possibile ancora utilizzare la formula di Green estesa al solo contorno
S1 .
24
Capitolo 3
I metodi a Pannelli
La formula di Green risulta particolarmente utile nelle applicazioni numeriche
quando il contorno del corpo non è descrivibile analiticamente, come spesso
accade per i profili aerodinamici: è infatti possibile esprimere la soluzione
in funzione di integrali estesi alla sola frontiera del dominio di integrazione.
Sulla funzione di Green si basa una recente categoria di metodi numerici,
detti metodi a pannelli. Si cominci con l’osservare che la formula di Green
costituisce una generalizzazione dei metodi a singolarità virtuali; considerando infatti un corpo non portante, ad esempio nel caso tridimensionale, e la
funzione di Green definita nella (2.18) per un punto che non appartiene al
contorno del corpo (⇒ E = 1) si ha:
1
φ(r0 ) =
4π
Z
∂V
∂1
φ r dS −
∂n
Z
∂V
1 ∂φ
dS
r ∂n
(3.1)
Il primo termine rappresenta il potenziale generato da una dipolo di intensità φ, mentre il secondo termine rappresenta il potenziale generato da una
∂φ
sorgente di intensità ∂n
. Analizzando la (2.24), si può notare che nel caso
di un corpo portante è presente un termine aggiuntivo costituito da una distribuzione di dipoli lungo la superficie della scia di intensità pari a ∆φ; con
la formula di Green è dunque possibile rappresentare anche il campo di moto
attorno ad un corpo portante, a patto di distribuire lungo la scia una distribuzione di dipoli, equivalente ad uno strato vorticoso: si può così riuscire
25
3. I metodi a Pannelli
a costruire una formulazione di potenziale che contempla la presenza di uno
strato vorticoso di spessore infinitesimo all’interno del campo di moto.
Osservazione 3 Una rappresentazione alternativa dello strato vorticoso può essere ottenuta sostituendo la distribuzione di dipoli lungo la scia con una distribuzione di vortici; questa seconda tecnica sarà alla base del metodo di Morino,
che verrà introdotto nella sezione 3.2.
Si può pensare di discretizzare il corpo aerodinamico in N pannelli rettilinei e porre su questi delle sorgenti, dando localmente una distribuzione nota
(tipicamente costante per il caso subsonico), ma dipendende da un certo numero di parametri incogniti. Per ottenere un corpo portante, occorrerà infine
aggiungere delle incognite ulteriori costituite da dipoli o vortici, sempre di
intensità incognita. Ne segue che gli integrali che compaiono nella (2.24) potranno esere calcolati analiticamente e costituiranno la matrice di un sistema
lineare. Per chiudere il problema occorre determinare il vettore dei termini
noti per le N + 1 equazioni nelle N + 1 incognite:1 a questo fine si utilizzeranno le N equazioni che impongono la condizione di non penetrazione più la
condizione di regolarità al bordo di uscita (condizione di Kutta); supponendo
che la numerazione dei pannelli sia tale per cui il primo e l’ultimo pannello
sono quelli del bordo di uscita, il tutto si traduce nelle seguenti condizioni:
Vi · ni = 0
i = 1, .., N
V1 · t1 = VN · tN
(3.2)
(3.3)
Si introducono adesso i due metodi a pannelli più utilizzati nel campo
aerodinamico: il metodo di Hess e Smith, che costituisce un metodo inverso,
in quanto risolve un sistema le cui incognite sono le intensità delle singolarità
virtuali sui pannelli; il metodo di Morino che costituisce un metodo diretto,
ovvero un metodo la cui incognita primaria è proprio il potenziale.
1
Si sta considerando per semplicità il caso di una distribuzione di potenziale costante
sul pannello; nel caso di una distribuzione, ad esempio, lineare tale per cui le incognite
sono 2 per pannello, si potrà imporre ad esempio la continuità del potenziale.
26
3. I metodi a Pannelli
3.1
Il metodo di Hess e Smith
Il metodo di Hess e Smith è un metodo del tipo a pannelli basato sull’ipotesi
di assumere nullo il potenziale relativo ai termini di dipolo e di descrivere
l’effetto portante attraverso una distribuzione di vortici sul contorno del corpo, essendo chiaramente descritto l’effetto dello spessore mediante termini di
sorgente. Il potenziale in un generico punto del campo di moto potrà dunque
essere scomposto in tre contributi:
φ(r0 ) = φ∞ (r0 ) + φs (r0 ) + φv (r0 )
(3.4)
relativi rispettivamente alla corrente asintotica, alla distribuzione di sorgenti ed alla distribuzione di vortici. Detto α l’angolo di incidenza fra la corrente indisturbata ed il profilo e considerando per semplicità il caso bidimensionale, le espressioni per ciascun contributo del potenziale sono le seguenti
(cfr. cap.2):
φ∞ (r0 ) = V∞ (x0 cos α + y0 sin α)
I
1
φs (r0 ) =
q(s) log |r − r0 |ds
2π
I
1
γ(s)θds
φv (r0 ) = −
2π
(3.5)
(3.6)
(3.7)
dove le distribuzioni di intensità γ(s) e q(s) sono incognite e dipenderanno
dalle condizioni al contorno. Si consideri ora la suddivisione del profilo in
pannelli rettilinei come figura 3.1; si supponga inoltre per semplicità che su
ogni pannello l’intensità q(s) delle sorgenti sia costante e che su tutto il corpo
si abbia un’intensità di vortice costante e pari a γ; il potenziale nel generico
punto r0 del campo di moto sarà pari a:
N Z
1 X
φ(r0 ) = V∞ (x0 cos α + y0 sin α) +
(qj log |r − r0 | − γθ) ds
2π j=1 pann.j
(3.8)
Espressione valutabile analiticamente e funzione degli N + 1 valori qj e γ.
Per determinare le intensità incognite ci si serve delle condizioni al contorno;
27
3. I metodi a Pannelli
Figura 3.1: Pannellizzazione del profilo
scegliendo come punti di controllo i punti medi di ciascun pannello e detta ni
la relativa normale, si ottengono le N relazioni : Vi · ni = 0; per l’incognita
N + 1 si utilizza invece la condizione di Kutta che, nel caso in cui il primo
e l’ultimo pannello vadano a cadere sul bordo di uscita, porta alla seguente:
V1 · n1 = VN · nN dove ti è la tangente all’i−esimo pannello.
Si scomponga ora la velocità nel punto di controllo del generico pannello i
nelle sue componenti orizzontale e vericale:
ui = V∞ cos α +
vi = V∞ sin α +
N
X
qj usij + γ
N
X
j=1
j=1
N
X
N
X
qj vsij + γ
j=1
uvij
vvij
(3.9)
(3.10)
j=1
dove con le quntità fgij rappresentano i coefficienti di influenza, ovvero la
componente di velocità orizzontale\vericale prodotta nel punto di controllo
del pannello i da una distribuzione unitaria di singolarità posta sul pannello
j e dipendono esclusivamente dalla geometria del problema e possono essere
calcolati una volta per tutte. Per quanto riguarda il termine sorgente, si
ottiene per le componenti di velocità in coordinate cartesiane:
1
x
2
2π x + y 2
1
y
vs =
2π x2 + y 2
us =
(3.11)
(3.12)
Passando, per comodità, ad un sistema di coordinate locali per il pannello
28
3. I metodi a Pannelli
Figura 3.2: Pannello j-esimo visto dal punto di controllo del pannello i-esimo
j − esimo:
uj = u∗j cos θj − vj∗ sin θj
(3.13)
vj = u∗j sin θj + vj∗ cos θj
(3.14)
dove per il significato dei simboli si faccia riferimento alla figura 3.2, ed
indicando nel nuovo sistema di riferimento con (x∗ , y ∗ ) le coordinate del punto
di controllo del pannello i−esimo, si ottiene:
t=lj
Z lj
q
x∗ − t
1
1
2
∗
2
+
log
(x
−
t)
+
y
=
dt
=
−
2π 0 (x∗ − t)2 + y ∗2
2π
t=0
1
rij+1
= − log
(3.15)
2π
rij
t=l
Z lj
y∗
1
x∗ − t j
1
∗
dt = − arctan
vsij =
2π 0 (x∗ − t)2 + y ∗2
2π
y ∗ t=0
βij
=
(3.16)
2π
u∗sij
In maniera analoga per il vortice:
βij
2π
1
rij+1
= − log
2π
rij
u∗vij = −
(3.17)
∗
vvij
(3.18)
29
3. I metodi a Pannelli
Assemblando la matrice dei coefficienti ed il termine noto si ottiene infine:
"
#( ) ( )
Aqq Aqγ
q
bq
=
(3.19)
Aγq Aγγ
γ
bγ
Si possono ora determinare le prestazioni del profilo a patto di determinare
prima per ogni pannello la componente tangenziale di velocità:
Vi · ti = V∞ cos(θi − α)
N
X
qj
rij+1
+
sin(θj − θi )βij − cos(θi − θj ) log
2π
rij
j=1
N γ X
rij+1
+
sin(θj − θi ) log
+ cos(θi − θj )βij
2π j=1
rij
(3.20)
Il coefficiente di pressione, valutato sul singolo pannello nel punto di controllo
vale:
cp (x̄i , ȳi ) = 1 −
Vi · ti
V∞
2
(3.21)
da cui è possibile ricavare la forza aerodinamica, andando ad integrare
lungo il contorno.
Si può ripetere il calcolo per un’incidenza diversa semplicemente ricostruendo
il termine noto, rimanendo invariata la matrice dei coefficienti.
3.2
Il metodo di Morino nel caso stazionario
Come il metodo di Hess e Smith anche il metodo di Morino è un metodo del
tipo a pannelli in cui però si determina il potenziale φ sul corpo aerodinamico
risolvendo un sistema di equazioni in cui l’icognita è il potenziale stesso.
Anche il metodo di Morino è basato sulla formula di Green e consiste in un
passaggio al limite nel far tendere il generico punto r0 del campo di moto al
generico punto rc del contorno; diversamente dal metodo di Hess e Smith,
questo metodo è basato su singolarità del tipo sorgenti e dipoli. Si consideri
la (??) con il termine g = 0 e si scompongano gli integrali in una parte Sw
30
3. I metodi a Pannelli
relativa alla scia ed in una parte Sb relativa al corpo. Si ricordi inoltre che
la scia costituisce una linea di flusso, ovvero: V · n| = 0 ⇒ ∂φ = 0,
Sw
∂n Sw
avendo allineato la scia con la corrente asintotica V∞ , parallela all’asse x,
per semplificare le espressioni che seguiranno. Per quanto riguarda l’integrale
lungo la scia si ottiene:
∂G(r0 , r)
∂φ(r)
φ(r)
− G(r0 , r)
dS =
∂n
∂n
Sw
Z
∂G(r0 , r)
=
∆φ(r)
∂n
Sw
Z
(3.22)
Inoltre la scia non deve produrre portanza, il che si traduce nella seguente
relazione:
∆cp |Sw = 0
∂∆φ ⇒ −2
= 0 ⇒ ∆φ|Sw = f (y)
∂x Sw
∂φ
cp = −2
∂x
(3.23)
avendo indicato con y la coordinata in apertura. Come si vede, il salto di
potenziale sulla scia varia solo in apertura, essendo costante nella direzione
della corrente asinotica: ne segue che il salto sarà eguale a quello che si ha
sul bordo d’uscita: ∆φ|Sw = ∆φBU ; inoltre, per la circolazione:
I
Γ=
S
∂φ
ds =
∂s
I
dφ = ∆φBU
(3.24)
Osservazione 4 Nel caso non stazionario il salto di potenziale è retto da un’equazione evolutiva, facendo dipendere così la soluzione dalla storia passata; in particolare non varrà più la relazione ∆φ|Sw = ∆φBU . Tipicamente questo è affrontato
mediante la trasformata di Fourier per la variabile temporale; questo argomernto
però esula dalle presenti note.
Osservazione 5 Nel caso comprimibile stazionario l’equazione per il potenziale
assumerebbe la forma:
2
(1 − M∞
)
∂2φ ∂2φ ∂2φ
+ 2 + 2 =0
∂x2
∂y
∂z
(3.25)
31
3. I metodi a Pannelli
Nel caso subsonico è possibile ricondursi nuovamente all’equazione di Laplace attraverso la trasformazione di Prandtl-Glauert:
x
X=p
2
1 − M∞
Y =y
Z=z
(3.26)
ed utilizzare nuovamente la teoria sviluppata finora per il caso incomprimibile,
avendo sempre a che fare con un’equazione differenziale di tipo ellittico.
Si ottiene infine:
Z ∂G(r0 , r)
∂φ(r)
1
φ(r0 ) =
φ(r)
− G(r0 , r)
dS
2
∂n
∂n
SB
Z
∂G(r0 , r)
∆φBU (r)
+
dS
∂n
Sw
(3.27)
Applicando ora la stessa pannellizzazione di figura 3.1 la (3.27) diventa:
N Z
X
∂G(ri , r)
1
∂φk
φk
φi =
− G(ri , r)
dS
2
∂n
∂n
S
B
k
Z k=1
∂G(ri , r)
dS
(3.28)
+
∆φBU k
∂n
Swk
Considerando, per semplicità, il caso in cui il potenziale e la sua derivata
normale siano funzioni costanti su ogni pannello e pari al valore al centro del
pannello stesso, si possono portare fuori dal segno di integrale, ottenendo:
N
X
1
φk
φi =
2
k=1
Z
SBk
Z
Z
N
N
X
X
∂G
∂φk
∂G
GdS +
∆φBU k
dS −
dS
∂n
∂n SB
Swk ∂n
k=1
k=1
k
(3.29)
che porta al seguente sistema lineare:
1
{φ} = [D] {φ} + [Z]
2
∂φ
∂n
(3.30)
Si può notare, andando a sostituire la funzione di Green, che la matrice
[D] rappresenta dei termini di dipolo, mentre la matrice [Z] rappresenta dei
termini di sorgente. Raccogliendo i termini comuni si ottiene infine:
∂φ
(3.31)
[Y ] {φ} = [Z]
∂n
32
3. I metodi a Pannelli
dove il termine a secondo membro è determinato mediante le condizioni al
contorno sul potenziale relative alla non compenetrazione.
Per calcolare ora le azioni aerodinamiche, noto il potenziale, si introduce una matrice di differenziazione [d] lungo la direzione x dopodichè, dalla
formula per il coefficiente di portanza, si ottiene:
∂φ
−1
{cpi } = [d] [Y ] [Z]
∂n
(3.32)
ottenedo direttamente l’espressione per i carichi aerodinamici.
3.2.1
Alcune osservazioni sulle matrici D e Z
Si è detto che la matrice D rappresenta dei termini di dipolo, parte del corpo
e parte della scia. Si scomponga allora:
[D] = [DB ] + [DW ]
(3.33)
Per quanto riguarda il termine [DW ], questo dipenderà in generale da come
è pannellata la scia; poichè ∆φBU varia solo in apertura, si può pensare di
pannellizare la scia con pannelli molto allungati nella direzione del moto;
inoltre, nella matrice [D] il contributo relativo alla scia si fa sentire solo nei
pannelli corrispondenti al bordo di uscita.
Per quanto riguarda la matrice [Z], andando a sostituire ad esempio la funzione di Green per il caso tridimensionale, è possibile notare che facendo
tendere a zero lo spessore l’integrale che contribuisce alla generazione di detti termini tende ad annullarsi, portando ad un’equazione (e quindi ad un
sistema) singolare: ne segue che il metodo di Morino può essere applicato
solo per corpi di spessore finito, mentre nel presenta un comportamento singolare per corpi di spessore nullo come le superfici portanti. Si è tuttavia
notato che può comunque essere utilizzato per la simulazione di superfici
portanti, a patto di utilizzare uno spessore fittizio abbastanza sottile.
33
Capitolo 4
Altri metodi numerici per
l’aerodinamica
Per studiare l’aerodinamica dell’ala oltre ai metodi appena introdotti ne esistono di semplificati, che forniscono solamente la distribuzione della differenza del coefficiente di portanza ∆cp fra dorso e ventre, ma non la distribuzione
che lo genera. Qualora si consideri l’ala come una superficie bidimensionale si utilizzeranno metodi detti a superficie portante, mentre se questa è
schematizzata con una linea, il metodo sarà detto a linea portante.
Si consideri una generica alla e si faccia tendere a zero lo spessore: si
ottiene una superficie in cui le condizioni al contorno sopra e sotto detta
superficie differiranno solo per il segno. inoltre, considerando la formulazione del problema con la funzione di Green si otterrà per il caso stazionario
un’equazione integrale della forma:
Z
α(x, y) =
K(x, y, ξ, η, M∞ )∆cp (ξ, η)ds
(4.1)
S
dove la funzione K(x, y, ξ, η, M∞ ) costituisce il nucleo dell’equazione integrale (funzione di influenza) che lega l’incidenza alla distribuzione del carico
aerodinamico. Si osservi che in K è compreso anche il contributo dovuto alla
scia.
Per risolvere numericamente il problema è possibile utilizzare un metodo
a singolarità virtuali.
35
4. Altri metodi numerici per l’aerodinamica
4.1
La linea portante-schema di Weissinger
Si discretizzi l’ala attraverso una distribuzione discreta di vortici a staffa
come indicato in figura 4.1
Figura 4.1: Pannellizzazione per la linea portante
I vortici vanno posti al 25% della corda (centro aerodinamico) e si estenderanno all’infinito. Per ottenere la velocità indotta nel generico punto P
ad opera questi vortici, si utilizza la legge di Biot-Savart, che per un filetto
vorticoso di lunghezza l e intensità Γ assume la seguente forma:
dl × rP O
3
l krP O k
−−−−−→
= (P − O),
Z
v(P ) =
rP O
(4.2)
O∈l
(4.3)
Osservazione 6 Passando alle coordinate polari nella (4.2) si trova l’espressione
per la parte di vortice che si estende all’infinito.
Per chiudere il problema ci si serve delle condizioni al contorno, scegliendo
come punto di collocazione il punto posto al 75% della corda alare (punto di
incidenza equivalente):
v · n = (v + V∞ ) · n = 0
(4.4)
Osservazione 7 Nel caso in cui si abbia un’ala rettangolare si può utilizzare un’unico vortice trovando lo stesso risultato della teoria della linea portante di Prandtl,
mentre nel caso in cui si abbia un’ala a freccia questa teoria prende il nome di Teoria di Weissinger. In generale, pur avendo un’ala a pianta rettangolare, è possibile
utilizzare più di un vortice, ottenendo una distribuzione discreta di portanza.
36
4. Altri metodi numerici per l’aerodinamica
4.2
La superficie portante
Nel seguito verranno analizzati due schemi di superficie portante: il VortexLattice method, che costituisce un’estensione dello schema di Weissinger e la
Piston Theory, che trova applicazione in campo supersonico.
4.2.1
Il Vortex-Lattice Method
Si può pensare di estendere la teoria della sezione 4.1, considerando l’ala
cone costituita da pannelli sia in apertura sia in corda, su cui sono collocati
dei vortici a staffa, come rappresentato in figura 4.2, solo che questa volta
i vortici ed i punti di controllo saranno rispettivamente al 25% ed al 75%
del pannello stesso. Anche in questo caso il sistema risolvente si ottiene
imponendo v · n = 0 negli N punti di collocazione.
Figura 4.2: Pannellizzazione per la superficie portante
Si ottiene il seguente sistema lineare:
1
1
{vn } =
[GΓ · n] {Γ}
V∞
VZ
∞
1
dl × rP O
GΓ · n =
4π l krP O k3
(4.5)
(4.6)
Ricordando che su ogni pannello la vorticità Γi è costante si ottiene il
seguente legame fra la circolazione e la portanza prodotta:
1
L = ρV∞ Γb = ρV∞2 Si ∆cp
2 1
Si
⇒ Γ = V∞
∆cp
2
bi
(4.7)
(4.8)
37
4. Altri metodi numerici per l’aerodinamica
dove Si è la superficie del singolo pannello e b la larghezza del vortice a staffa.
Andando a sostituire nella (4.5) si ottiene infine il legame cercato:
{α} = [A] {∆cp }
1
S
[A] = [GΓ · n]
2
b ii
(4.9)
(4.10)
Considerando ora un’ala a bassa incidenza ed a bassa curvatura di linea
media, la condizione al contorno per il pannello i − esimo diventa:
!
∂z
=0
(4.11)
vi · ni = wi − vi tan Λ + V∞ α −
∂x lmi
∂z
rappresenta la pendenza
dove α rappresenta l’angolo di incidenza, ∂x
lmi
locale della linea media, Λ rappresenta l’angolo di diedro, vi la componente
di velocità in direzione dell’apertura e wi la velocità in direzione ortogonale
alla superficie. Dividendo la (4.11) per V∞ si ottiene il termine noto del
sistema:
vi
αi =
tan Λ +
V∞
4.2.2
α−
∂z
∂x
!
=0
(4.12)
lmi
La Piston Theory
La piston Theory è un metodo a superficie portante utilizzato per flussi
supersonici quando M∞ > 2 e consiste nel considerare una perturbazione di
pressione rispetto al valore asintotico del tipo:
p = ρ∞ c∞ vn
(4.13)
dove c∞ è il valore asintotico della velocità del suono. Questa teoria si
basa sul considerare l’ala come costituita da una serie di pistoni affiancati
posti in direzione normale alla superficie alare che non si influenzano fra loro;1
1
Quest’ultima ipotesi giustifica la condizione di applicabilità per M∞ > 2: perché i
pistoncini siano indipendenti si deve avere un cono di Mach di ampiezza nulla e per il
valore suddetto è abbastanza stretto da rendere accettabili i risultati su cui si basa questa
teoria.
38
4. Altri metodi numerici per l’aerodinamica
il coefficiente di portanza e l’incidenza avrnno la seguente forma:
cp =
α=
vn
V∞
p
1
ρ V2
2 ∞ ∞
=
2 vn
M∞ V∞
(4.14)
(4.15)
da cui:
∆cp = 2cp =
4
α
M∞
(4.16)
39