Gli studi fluidodinamici di

laboratorio e la tecnica PIV

Massimiliano Manfrin - Turlab, Dip. Fisica, Università di Torino - manfrin@ph.unito.it

• Campi applicativi che spaziano dalla fisica dei fluidi

all’ingegneria (meccanica, aeronautica e aerospaziale,

civile, ...), dalla medicina alla meteorologia, dalla

biologia all’oceanografia……

Fluidodinamica: studio delle forze che agiscono sui fluidi

(gas e liquidi).

L’equazione «omnicomprensiva» è quella di Navier-Stokes.

2

Si ricordi che :

• La fisica studia le equazioni del moto.

• Se si determinano i campi di velocità, e lo si riesce a fare nel tempo, il problema dello studio del fenomeno è potenzialmente risolto.

• Posso ottenere le accelerazioni e quindi le forze; posso descrivere il moto.

3

Necessità per gli studi in laboratorio di avere anche

apparecchiature rotanti.

Perché una vasca rotante ?

Sulla Terra i processi atmosferici che

coinvolgono i fluidi (atmosfera e oceano)

si svolgono in un sistema di riferimento

rotante.

Presenza della forza di Coriolis: agisce

solo su una particella in moto perturban-

done la direzione del vettore velocità, ma

non il modulo.

Forza particolarmente importante su

grandi scale e per alte velocità.

Finestre per

permettere

di

fotografare il

flusso

Tubo per

l’aggiunta di

acqua

Collettore

strisciante

per la

trasmissione

dei dati

Aspetto unico: canali di vetro sul fondo che permettono la

creazione di lame verticali potenti e ben definite.

Laser verde a diodi allo stato solido: potenza massima 15 W

trasmessa direttamente (senza fibra ottica).

DiametroDiametroDiametroDiametro vascavascavascavasca 5 5 5 5 mmmm

AltezzaAltezzaAltezzaAltezza massimamassimamassimamassima acquaacquaacquaacqua 0.8 0.8 0.8 0.8 mmmm

VelocitVelocitVelocitVelocità rotazionerotazionerotazionerotazione 0000----20 20 20 20 rpmrpmrpmrpm

5

I flussi nella vasca possono essere generati con

pompe oppure, con, in generale migliori risultati,

variando bruscamente il periodo di rotazione

dello strumento dopo che si è ottenuta la

condizione di rotazione solida

Tini > Tfin spin-up Tini < Tfin spin-down

Il flusso risultante ha una velocità che decresce

esponenzialmente e il cui valore massimo è dato da:

Umax=2π(1/Tfin-1/Tini)r

(r distanza dal centro )

Il moto della vasca è generato

da un motore elettrico ed il

sistema di rotazione si basa

sulla tecnologia del velo

d’olio ad alta pressione che

riduce gli attriti e le

vibrazioni.

6

1 Mikrotron MC1310:

CMOS

1-500 Hz a piena

risoluzione 1280x1024;

501- 1000 Hz at 1/2 ris …

8 bit

2 Dalsa 1M60P: 1-60 Hz

CCD 1024x1024

12-10-8 bit

2 Dalsa 4M60P: 1-60 Hz

CMOS 2352x1728

10-8 bit

1 Flare2MC320:

CMOS,

1-280 Hz,

a colori,

2048x1088,

8 bit

A causa dell’attrito sul fondo il flusso nello

strato limite diviene

cross-isobarico

Dovendo conservarsi la massa si

sviluppa una circolazione

secondaria (Ekman pumping)

Spin-down Spin-up

Convergenza (L) Divergenza (H)

Tipi di Circolazione

8

E’ la nota Spirale di Ekman in

cui il vettore velocità è ruotato

verso sinistra avvicinandosi al

suolo e contemporaneamente il

modulo del vettore diminuisce.

Andamento della velocità nel fluido in moto in una vasca in rotazione.

9

Importanza di poter riprodurre i fenomeni

naturali in laboratorio:

- Maggiore possibilità di controllare le condizioni

sperimentali e al contorno;

- Possibilità di eliminare o ridurre l’importanza di alcune

delle forze (analisi di scala) in azione potendo così scindere il

problema in sotto-problemi di più semplice interpretazione;

- Possibilità di variare alcuni parametri per verificarne

l’influenza/importanza;

- Fornire misure sperimentali per validare i modelli numerici;

- Ripetibilità dell’esperimento;

• Esperimenti in laboratorio devono quindi riprodurre opportunamente le

condizioni naturali.

• Le tecniche di misura non dovrebbero essere intrusive né perturbative.

• La misura di per sé comporta una perturbazione che deve essere

minimizzata per «non» influenzare il fenomeno osservato.

10

Fenomeni naturali (prototipo) e riproduzioni in

laboratorio (modello): rappresentazione in

laboratorio di fenomeni naturali a scale molto

diverse: dalla microscala alla scala sinottica, quindi

dai cm alle migliaia di km.

Come far sì che lo studio di laboratorio sia

rappresentativo del prototipo naturale?

11

• Si dice che esiste similarità fra due quantità

quando è possibile metterle in relazione

attraverso una costante adimensionale

(indipendente dal sistema di coordinate)

Teoria della similitudine

La tecnica della Similitudine

Nella riproduzione dinamica su modello di unfenomeno naturale occorre che il pattern del flussonel prototipo sia una copia ingrandita e non distortadel “pattern” del flusso nel modello.

12

Si ha similitudine geometrica fra due sistemi quando è

possibile far coincidere uno di essi con l’altro mediante un

opportuno cambiamento di unità di misura delle lunghezze

(è la meno importante fra le similitudini possibili).

Si dice che due “oggetti” sono cinematicamente simili se è possibile trovare due SR, tali che gli “oggetti”siano geometricamente simili per

ogni t e che il rapporto di similitudine geometrica sia indipendente da t.

Se le rispettive masse

stanno fra loro in un

rapporto costante μ,

allora essi sono anche

materialmente simili

13

è necessario che, in corrispondenza di punti

geometricamente similari, le forze agenti sulle

particelle siano in rapporti fissi ad istanti omologhi

di tempo.

Similarità dinamica

Se esiste similarità dinamica fra dueprocessi fisici, è possibile trasferire irisultati ottenuti per il primo al secondo, eviceversa («più facile a dirsi che non afarsi»……)

14

Questo si traduce nella conservazione di alcune quantità

adimensionali che vengono individuate da opportuni rapporti di forze

(le variabili sono variabili di scala).

Numero di Rossby: importanza della rotazione [U/(f*L)];

Numero di Reynolds: importanza della viscosità [UL/υ]

Ad esempio, nel moto dei grandi sistemi atmosferici, anticicloni apportatori di

bel tempo e depressioni apportatrici di vento e pioggia, si ha, alle nostre

latitudini,

U= 10 m/s L= 1000 Km f=10-4

Ro=0.1

Se si vuole riprodurre su scale di laboratorio tale moto, riducendo la velocità a

0.1 m/s e la lunghezza a 10 m, f deve essere 0.1.

Come posso determinare i campi di velocità ?

Uso di traccianti passivi: determinazione del moto nel

tempo

15

Se voglio fare delle valutazioni quantitative ho bisogno di un

approccio che mi registri qualcosa permettendomi una

analisi successiva, in qualsiasi momento, indipendente

dall’acquisizione e dall’acquisitore.

Conviene riprendere con una

foto/telecamera e poi analizzare

quanto catturato.

Devo fare anche in modo che il

risultato sia scientificamente

indipendente da chi effettua le

misure o le analisi.

16

Due approcci:

Euleriano LagrangianoSeleziono un tracciante ,

«etichettandolo» e lo seguo nelle

sue evoluzioni.

Seleziono un’area di mio interesse

ed eseguo le misure in quell’area

valutando sempre nuovi traccianti

che entrano ed escono dal campo

di misura

Ho un’ottima visione d’insieme

nell’area di misura, che però è

limitata:

Ho un’ottima valutazione del moto

di quel particolare elemento nello

spazio e nel tempo ma perdo la

visione d’insieme:

PIV/LIF

Particle Image

Velocimetry / Laser

Induced Fluorescency

PTV

Particle Tracking

VelocimetryMeno usata perché richiede zone di

misura molto ampie e/o sistemi di

inseguimento complessi.

Fornisce informazioni solo in

pochissimi punti: è limitante!

17

Permettono misure quantitative

-Densità =1.03 g/cm3

-Diametro: 5-1000 micron

-Necessitano di un agente

bagnante per renderle

sospensibili

-Buona riflettività

-Forma simil-sferica

-Economiche

Approcci sperimentali

Le più piccole (5-10 micron) rimangono in

sospensione, in acqua dolce, per alcune

decine di ore; quelle di grandi dimensioni

(1000 micron) per qualche minuto.

ParticelleInchiostri

Permettono misure qualitative

-Gestione molto più

semplice e veloce

-Normali o fluorescenti

-Utili per studi di diffusione

PIV

(Particle Image Velocimetry)

Abilità di definire campi di flusso istantanei su aree e volumi che

possono spaziare di diversi ordini di grandezza.

Possibile identificare strutture spaziali in flussi non stazionari e studiarli

quantitativamente.

Tecnica recente, nata negli anni 70.

18

Le riprese possono essere bidimensionali (una telecamera; setup

molto più «semplici») o tridimensionali (due o più telecamere; setup

molto più complessi).

Come devono essere i traccianti ?

PICCOLI !!!

RIFLETTENTI (in realtà è uno scattering di Mie)

Possibilmente in equilibrio idrostatico con il fluido

Rappresentazione «fedele» di un elemento fluido

Quando vengono illuminati da una opportuna luce

diventano visibili (come la polvere in controluce)19

- è facile da maneggiare (comporta qualche rischio però);

- può avere potenze ben regolabili ed elevate;

- l’assorbimento da parte dell’acqua è minimizzato;

- è fra le lunghezze d’onda più facilmente creabili;

Nelle riprese, bisogna stare molto attenti alle aberrazioni:

- introdotte dagli obiettivi e/o dal passaggio di mezzi

diversi (rifrazione),

- dalle deformazioni prospettiche,

- dalla luminosità di background (le particelle devono

essere ben distinguibili le une dalle altre, risaltare rispetto

al bg, avere una distribuzione di luminosità simil-gaussiana

con bordi definiti ma sfumati (tecniche di blurring).

Utilizzata luce verde laser:

20

E’ poi molto importante definire bene a priori (deve essere

studiata in anticipo la fisica che si vuole riprodurre) i tempi di

acquisizione e quindi la differenza temporale tra due

immagini successive; e definire i tempi di esposizione per

evitare «l’effetto cometa» (sovraesposizione).

Ci deve essere un buon seeding: non poche (pochi punti, radi

e poco «legati»), non troppe (punti sufficienti, ma rischio di

eccessiva luminosità e sovrapposizione).

21

I tempi devono

essere grandi

abbastanza da

avere un moto

«sensibile», ma

piccoli a

sufficienza da non

perder di vista le

particelle

inquadrate nella

prima immagine.

Due possibili approcci di acquisizione:

- In continua (vantaggio di poter «scegliere»

l’accoppiamento dei campi e di non dover sapere con

troppa precisione la fisica in anticipo; svantaggio di avere

una grande mole di dati);

- A bursts: si scelgono 2-7 tempi di acquisizione

consecutivi con intervalli che consentano la massima

flessibilità temporale; poi pausa di consistente lunghezza.

Inizialmente la ripresa veniva effettuata con camere analogiche su

rullino (solo con il secondo approccio) .

Grandissima risoluzione ma problemi nell’analisi che doveva

essere fatta a mano su una o due coppie di immagini

considerando «particella a particella»; tempi di valutazione di

qualche giorno, fino anche a qualche settimana !!22

Con l’avvento delle telecamere digitali è aumentato

enormemente il numero di immagini acquisite (sebbene

almeno all’inizio a prezzo di risoluzioni e sensibilità molto

inferiori).

I tempi sarebbe esplosi se non si fosse provveduto ad

elaborare tecniche digitali di analisi.

In questo modo si possono decidere i parametri inizialmente

e poi i computer si occupano di iterare l’analisi lungo tutta la

serie (le condizioni fisiche non devono mutare troppo !!).

Come si può mappare lo spazio delle immagini (in pixel) nello

spazio fisico reale?

E come si può valutare lo spostamento dei traccianti?23

Problema: Le equazioni del moto sono

espresse in unità fisiche e non in pixel

24

Prima di passare all’implementazione via software della PIV è necessario mappare le posizioni delle particelle nel piano dell’immagine (pixel) al piano della lama laser (cm).

Il primo step per effettuare la

calibrazione è l’allineamento

della lama laser con un

bersaglio dalle dimensioni

note (ideale struttura

quadrettata o crocettata),

assicurandosi che sia ben

centrato.

La calibrazione

L’informazione data

dall’immagine è in pixel, ma il

numero di questi dipende

strettamente dal tipo di obiettivo.

Acquisita un’immagine è necessario localizzare il centro delle croci. I software in uso sfruttano un procedimento simile alla tecnica usata nella PIV.Sono però imprecisi, essendo automatici.

La localizzazione

manuale del centro di

tutte le croci reca

ottimi risultati in

termini di calibrazione,

ma richiede tempi

lunghissimi.

25

26

I centri celle croci opportunamente mappati nella realtà fisica vengono

poi fittati da una superficie curva che minimizza gli scarti e fornisce alcuni

parametri del tipo di ripresa (distanza focale, centro ottico,….)

Il risultato permette di

avere una funzione

definita in tutto lo

spazio inquadrato e di

stimare le possibili

deformazioni dovute a

difetti delle lenti, alla

loro curvatura, a

problemi di parallasse,

…..

27

28

La disposizione angolare ha un inconveniente: produce campi di immagini

stirati in modo opposto.

Essendo la distorsione delle

due telecamere diversa, è

impossibile sovrapporre

direttamente le due aree di

interrogazione senza aver

prima interpolato i dati su

una griglia comune!

Calibrazione Stereo

29

Bersaglio «raddrizzato» mediante la

calibrazione.

Bersaglio come osservato dalla

telecamera di destra in un approccio

stereo

PIV

Analisi di coppie di immagini legate temporalmente (se si è sbagliata la

fase sperimentale, a questo livello non resta che rifare tutto….).

Tra le due immagini è necessario che le particelle si spostino di alcuni pixels (best tra 5 e 15).

x x

z

Lo spostamento osservato dipende da frame rate di acquisizione, dalla velocità del fluido, dalla risoluzione e dalle dimensioni del campo.

30

Partendo dalle intensità luminose (a maggiore risoluzione in bit dell’immagine corrisponde maggior accuratezza del risultato) relative a ciascun pixel (in toni di grigio), viene calcolato il coefficiente di cross-correlazione tra le pattern box dell’immagine a e b.

Esistono diverse possibili implementazioni della PIV. Quella adottata è chiamata CIV (Correlation Image Velocimetry) e si basa sul calcolo diretto (più lento, ma più preciso del calcolo tramite FFT) della cross-correlazione (normalizzata):

31

32

Si valuta il valore di tale funzione per ogni pattern-box

possibile all’interno della search box della seconda immagine,

generando una superficie topografica di cross-correlazione, di

cui viene individuato il massimo.

33

X

Origine del

vettore

velocitàSearch box

Pattern box

Il vettore velocità è assegnato temporalmente all’immagine a della coppia, mentre il suo punto di applicazione è posto nel punto intermedio del segmento congiungente il punto corrispondente al

centro della pattern box nella prima immagine

ed il

centro della pattern box della seconda immagine

dtW

yy

dtU

xx

a

a

2

2

+=

+=

Dal massimo di correlazione si ricava lo spostamento più probabile e,

conseguentemente, il vettore della velocità (si conosce il tempo trascorso tra

le due immagini).

Sono poi utilizzati algoritmi gaussiani di interpolazione subpixel che migliorano

la risoluzione a 0.1 px. 34

Lo spostamento identificato in questo modo ha risoluzioni di 1 pixel (px).

t

t +∆t

L’analisi viene poi ripetuta per tutte le possibili pattern box

identificabili all’interno delle immagini.

Quante sono le possibili pattern box?

Dipende da chi fa le analisi.

Si può avere overlapping da 0 e fino al 50 %; oltre, i punti

sono così correlati da non presentare più informazione e

proporre invece problemi in fase di elaborazione.

35

Nel corso del calcolo PIV vengono computati anche alcuni

flag associati ad ogni vettore e che ne caratterizzano la

qualità: soprattutto valore del coefficiente di correlazione e

localizzazione del massimo di correlazione (un massimo sui

bordi può non essere un massimo).

E’ una analisi statistica: lo spostamento indicato non è

verità assoluta, ma il numero che ha la massima

probabilità di rappresentare lo spostamento vero.

L’errore finale è, nel migliore dei casi (immagini

«perfette»), dell’ordine di qualche decimo di px (in

termini «fisici» l’errore dipende dalla risoluzione).

Quali sono gli errori di questo tipo di analisi?

Non si può valutare analiticamente, ma solo stimare da

medie di ensemble o metodi montecarlo (creazione di

immagini sintetiche dalla fisica nota, processing,

valutazione delle differenze).

Si ricordi che:

36

37

- Definire l’area media (pattern box) contenente almeno cinque

particelle: ad ogni area corrisponderà un vettore del campo di

velocità; dimensioni possono essere diverse lungo le due direzioni; le

box possono essere rettangolari.

- Stabilire un grigliato con su cui centrare le finestre;

- Definire l’area (search box; indicativamente circa 30 px più grande

della pattern box) entro la quale si muoverà la finestra per calcolare il

coefficiente di correlazione (bidimensionale);

- Dimensioni delle box sempre dispari in modo da assegnare il centro ad

un px reale e non interpolato.

Analisi PIV

38

Tutto questo processo può (e «deve») venire iterato più

volte portando così ad un miglioramento del risultato:

Il passaggio precedente viene utilizzato per imporre uno

spostamento piccolo nel suo intorno e non muoversi più su

tutta l’area identificata dalla search box.

Al secondo passaggio vengono introdotti inoltre algoritmi di

stretching e di convect flow.

Non è però possibile utilizzare le analisi direttamente a

cascata a causa della presenza di vettori errati.

39

Interfaccia grafica matlab per analisi, filtraggio e patch, sviluppata presso Coriolis-Legi - Grenoble

40

1) Utilizzo di alcuni flag computati in fase di

analisi.

Necessità di eliminare tali vettori errati.

Il campo calcolato contiene vettori errati dovuti a:

possibili problemi di illuminazione, seeding non uniforme,

agglomerati di particelle, poche particelle, parallassi, non buona

messa a fuoco, ….

Scartare i vettori ad «occhio» è relativamente semplice e molto

accurato (soprattutto nel caso di flussi ben definiti), ma farlo per

5000 vettori per campo e in 1000, 10000 campi può essere una

faccenda lunghetta…….

Si usano quindi software da impostare e poi applicare

iterativamente: è un metodo molto meno accurato ma molto più

veloce (…. e meno alienante…).

Un approccio molto rapido è:

- valutazione dei valori di correlazione (valori bassi

caratterizzano vettori spesso errati) calcolati in fase di

analisi;

- localizzazione dei massimi (sui bordi della finestra i vettori

possono essere errati, ma solo qualche volta).

41

42

Questo tipo di filtraggio può essere svolto attraverso l’interfaccia

precedente o quella qui presente che prepara ad un approccio più

raffinato

43

Scatter plot (o istogramma) nel piano

(u,v): permette di identificare (alcune

volte) aree non fisiche e quindi di

eliminare i punti che vi cadono.

2) Creazione di uno scatter plot

44

45

peaklocking

I valori

corrispondenti a

pixel interi sono

più probabili.

46

3) Filtraggio su base locale, utilizzando intorni

di opportuna grandezza di cui vengono valutati

media/mediana, std dev, abs dev (possibili

anche filtraggi tramite FFT)

Confronto della consistenza di ogni vettore con un suo

intorno (magari definito di dimensioni diverse, sia nello

spazio che nel tempo).

Approccio più raffinato ma decisamente più lungo e

soggetto alla scelta, soggettiva, di numerosi parametri.

Calcolo di valori medi, mediani e std,

massimi assoluti.

47

Interfaccia per la scelta dei parametri nel

filtraggio con l’approccio locale

48

L’applicazione del filtraggio solo tramite flag, nel caso sotto

riportato, pone dei problemi.

49

Invece l’applicazione dei filtri più sofisticati è molto più sensibile

I vettori filtrati vengono scartati (non rimossi dai

campi). Il campo risultante non è più omogeneo e non

è neanche su un grigliato regolare (questo già

dall’inizio per l’utilizzo del punto di applicazione).

Per evitare problemi in fase dell’applicazione delle

equazioni, si procede alla sostituzione dei vettori

scartati e alla regolarizzazione della griglia con

processi di interpolazione eseguiti su sottodomini (in

modo da preservare al massimo eventuali

disomogeneità del flusso) effettuati mediante spline

(che garantiscono anche la continuità delle derivate).

50