URČENÍ PARAMETRŮ FLUKTUAČNÍCH RYCHLOSTÍ TURBULENCE V MÍCHANÉ NÁDOBĚ Ing. Vít Pešava Školitel: Prof. Ing. Pavel Ditl, DrSc. České vysoké učení technické v Praze, Fakulta strojní, Ústav procesní a zpracovatelské techniky
Abstrakt Při návrhu mnoha procesních zařízení je třeba znát sedimentační rychlost částic, která je ovlivněna turbulencí vsádky. Okrajovou podmínkou při výpočtu této sedimentační rychlosti jsou časové závislosti fluktuačních rychlostí turbulence. Cílem této práce je stanovit parametry potřebné k popisu těchto závislostí. Záznam fluktuačních rychlostí obsahuje velké množství dat. Vysoký počet koeficientů, kterými lze plně aproximovat naměřený časový průběh je při dalších výpočtech nevhodný. Modelování sedimentace částic je pak časově extrémně náročné. Proto byl navrhnut a ověřen postup redukce počtu koeficientů aproximujících fluktuace axiálních složek rychlosti kapaliny. Nejjednodušším, ale nejméně přesným případem je aproximace fluktuačních rychlostí turbulence sinusovým průběhem. Cílem bylo vyvinout metodiku, pomocí které lze z grafů bezrozměrných rychlostí určit vstupní hodnoty modelu sedimentace pro různé geometrie míchaného systému a různé otáčky. Klíčová slova šesti lopatkové míchadlo, LDA, FFT, turbulence, sedimentační rychlost 1. Úvod Metodou Laser Doppler Anemometry (LDA) byly získány časové záznamy axiální rychlosti uz v určitých místech turbulentního proudového pole v míchané nádobě. Tato místa jsou reprezentována objemem měřící sondy tvořené křížícími se paprsky laseru. Sonda má tvar rotačního elipsoidu o rozměrech přibližně 2x1 mm. Ukázka záznamu je na obrázku 1. __
Podle [1] lze rychlost uz v určitém časovém okamžiku rozložit na střední hodnotu u z a její fluktuační složku u z' . Okamžitá rychlost bude potom rovna součtu těchto složek: __
u z = u z + u z'
(1)
Střední hodnota rychlosti je v míchané nádobě při nezměněných otáčkách konstantní. Můžeme ji tedy určit dle vztahu: _
uz =
t +T
ò t
u z dt =
1 N
N
åu 1
zi
(2)
Časovou závislost fluktuačních rychlostí pak určíme ze vzorce (1), přičemž střední hodnota fluktuačních rychlostí je nulová. Střední hodnota čtverce fluktuační rychlosti je charakterizována tzv. urms (root mean square velocity), kterou vypočítáme dle vztahu: ____
u rms = u z' 2 =
1 N
N
å (u i =1
' 2 zi
)
(3)
Vztahy (2) a (3) platí pro diskrétní hodnoty za předpokladu konstantního časového kroku __
mezi jednotlivými hodnotami. Výše zmíněné charakteristiky u z a u z' je výhodné vyhodnotit v bezrozměrném tvaru, kdy jednotlivé hodnoty podělíme obvodovou rychlostí míchadla v = p × n × D . Výsledkem jsou univerzální charakteristiky společné pro různé otáčky a geometricky podobné nádoby různé velikosti. Další charakteristikou turbulence je Intenzita turbulence. V případě izotropní turbulence můžeme její intenzitu vypočítat z následujícího vztahu: I=
u rms _
(4)
u
Výše zmíněné závislosti budou využity pro určování parametrů sinusové vlny – viz dále. Například podle [2] je rozdíl mezi sedimentační rychlostí částice v klidné kapalině a v turbulentním poli zapříčiněn změnou součinitele odporu prostředí CD v důsledku interakce mezi fluktuacemi kapaliny a částicí. Proto byl vyšetřován časový průběh fluktuačních rychlostí u‘ = u‘(t), který je okrajovou podmínkou do matematického modelu sedimentace částice v kapalině. Cílem této práce je: a) Aproximovat časový záznam naměřených fluktuací rychlosti pomocí vhodné matematické funkce tak, aby koeficienty této funkce mohli být použity v modelu sedimentace. b) Vhodným způsobem snížit počet koeficientů potřebných pro aproximaci původního záznamu. c) Určit parametry sinusové vlny tak, aby tato vlna dokázala aproximovat s co nejnižší chybou původní záznam. Důvodem pro úkoly uvedené v bodech b, a c, je vysoká výpočetní náročnost při použití plné aproximace původního záznamu. Vhodnou matematickou metodou pro aproximaci záznamu turbulentních fluktuací je Fourierova transformace. Tato metoda převádí časovou závislost rychlostí na závislost frekvenční. Tedy v případě diskrétních hodnot (tzn. použití diskrétní Fourierovy transformace) na součet sinusovek příslušejících jednotlivým frekvencím. Tímto způsobem je možné s vysokou přesností sestrojit s pomocí parametrů jednotlivých sinusovek původní záznam. Další výhodou vyjádření záznamu ve frekvenční oblasti je možnost jeho analýzy s využitím Kolmogorovovy teorie turbulence. Tato vlastnost je výhodná pro snížení počtu koeficientů pro aproximaci původního záznamu i jeho nahrazení sinusovým průběhem. Její využití bude popsáno později. Podle [3] lze převést funkci času na funkci frekvence prostřednictvím dopředné Fourierovy transformace dle následujícího předpisu:
~
+¥
u ( f ) = ò u (t )e 2pift dt
(5)
-¥
V této práci byla použita rychlá Fourierova transformace, která je součástí Matlabu (blíže viz [4]). Nutnou podmínkou správné funkce Fourierovy transformace jsou konstantně vzorkovaná vstupní data (stejný časový krok mezi jednotlivými hodnotami). Této podmínce však záznam pořízený metodou LDA nevyhovuje (obr. 1). Proto bylo třeba upravit vzorkování dat. V tomto případě byla použita metoda „Sample and Hold“. Při použití této metody nejprve definujeme vlastní časový vektor s konstantním krokem. Na pozice příslušných rychlostí dosadíme hodnoty rychlostí odpovídající nejbližšímu nižšímu času původního časového vektoru [5]. __
Dále byla vypočtena střední rychlost záznamu u z . Odečtením střední rychlosti od upraveného vektoru rychlostí byl vypočten vektor fluktuačních rychlostí u´, který byl dále zpracováván.
Obr. 1. Záznam měření okamžitých rychlostí turbulence pomocí LDA. Takto upravená data již můžeme pomocí FFT převést do frekvenční oblasti. Výsledkem je soubor n sinusovek. Jejich počet odpovídá Nyquistově frekvenci (tzn. n = polovina počtu členů jednotlivých vektorů upravených dat). Vyšetřovaný záznam lze velice přesně aproximovat jako součet těchto sinusovek. Jednotlivé sinusovky jsou zde jednoznačně určeny pomocí tří parametrů: amplitudy C, frekvence f a počáteční fáze φo. Z výše uvedených dat byla dále odfiltrována část záznamu odpovídající makronestabilitám, které neovlivňují sedimentaci částice. Odfiltrovány byly s použitím výsledků v [6] frekvence v rozsahu 0 – 2 Hz. Všechny výše zmíněné úpravy byly prováděny pomocí Matlabu. 2. Experimenty Záznamy rychlostí proudění v nádobě byly pro šesti lopatkové míchadlo se skloněnými lopatkami naměřeny pomocí metody LDA. Většina předkládaných výsledků pochází z těchto měření. Dále byly provedeny měření pomocí metody Particle Image Velocimetry (PIV) pro míchanou nádobu se šesti lopatkovým turbinovým míchadlem s dělícím kotoučem (Rushtonova turbína). Tato měření proběhla v letošním roce a stále probíhá vyhodnocování naměřených hodnot. Výsledky měření budou použity ke stejnému účelu jako výsledky LDA.
2.1 LDA měření Experimenty byly prováděny ve skleněné válcové nádobě s průměrem D = 300 mm s plochým dnem. Nádoba byla vybavena čtyřmi radiálními narážkami šířky b/D = 0,1. Míchadlo bylo šesti lopatkové se šikmými lopatkami skloněnými pod úhlem 45° dle normy CVS 691020 o průměru d/D = 1/3. Nádoba byla naplněna destilovanou vodou do výšky H/T = 1. Míchaná nádoba byla uložena ve vnější skleněné nádobě čtvercového půdorysu. Prostor mezi nádobami byl vyplněn destilovanou vodou z důvodu minimalizace vlivu zakřivení válcové nádoby na paprsky laseru.
b
H d
h2
Oblast měření
D
Obr. 2. Schéma měřené nádoby Rovina, ve které byly měřeny rychlosti turbulentního proudění je znázorněna na obrázku 2. Všechny záznamy byly naměřeny v turbulentní oblasti. Hodnoty Reynoldsových čísel, při kterých byla prováděna jednotlivá měření a příslušné otáčky míchadla jsou uvedeny v tabulce 1. Tabulka 1. – Otáčky a Reynoldsova čísla během měření. Otáčky míchadla [1/s] 300 450 600
Re [-] 49 900 74 850 99 800
Měření byla prováděna argonovým laserem CASIX LDC 1500 s maximálním výkonem 100 mW a vlnovou délkou 532 nm. Přenosová optika DANTEC 57x, přijímací optika s fotonásobičem DANTEC 57x08 a procesor DANTEC BSA. Jako značkovací částice byly použity S-HGS (Silver coated – Hollow Glass Spheres) se středním průměrem 10 µm a hustotou 1,1 . 2.2 PIV měření Experimenty byly prováděny ve skleněné válcové nádobě s průměrem D = 300 mm s plochým dnem (obr. 3). Nádoba byla vybavena čtyřmi radiálními narážkami šířky b/T = 0,1. Míchadlo bylo šesti lopatkové turbinové s dělícím kotoučem CVS 691021 o průměru d/D = 1/3. Výška spodní hrany míchadla nade dnem byla h2/D = 1/4. Nádoba byla naplněna v první části měření destilovanou vodou a v druhé části měření roztokem etylenglykolu s vodou do výšky H/D = 1. Míchaná nádoba byla uložena ve vnější skleněné nádobě čtvercového půdorysu. Prostor mezi nádobami byl vyplněn destilovanou vodou z důvodu minimalizace vlivu zakřivení válcové nádoby na paprsky laseru. Jednotlivé varianty měření jsou uvedeny v tabulkách 2 a 3. b
oblast1/1a
H
r1 d
s1
v1 h1
v1‘ h1‘ h2
r1‘ oblast2
s1‘
D
Obr. 3. Parametry míchané nádoby a měřené oblasti.
Tabulka 2. – Seznam měření Otáčky Oblast Víko [1/min] 300, 450, 1 ano 600
Kapalina*/
Záznamy (pro každé otáčky)
Teplota kapaliny v nádobě [°C]
voda
2x600ms, 1x5000ms
25
1
ne
300
voda
2x600ms, 1x5000ms
25
2
ne
300
voda
1x5000ms
23
2
ano
300, 450, 600
voda
3x5000ms
23
2
ano
600
etylenglykol
3x5000ms
20,3; 21,6; 22,6
1a
ano
600
etylenglykol
3x5000ms
23,0; 23,3; 23,5
*/ Voda = z kohoutku upravená ultrazvukem (zbavena bublin). Etylenglykol = roztok etylenglykolu s vodou. Kinematická viskozita n = 2 × 10 -6 m 2 / s při 20 °C. Měřící aparatura se skládala z laseru Litron Class 4 Laser Product, čočky Dantec Light Sheet, kamery Dantec SpeedSense 611 a objektivu Sigma MacroDg. Pro zlepšení rozlišení kamery byly použity fluorescenční značkovací částice o velikosti 11,95 µm (střední odchylka 0,25 µm) a optický filtr propouštějící světlo o vlnové délce 570 nm. Fluorescenční částice osvětlené laserem při vlnové délce 532 nm emitují světlo o vlnové délce 570 nm. Poloha těchto částic je během měření zaznamenávána kamerou vybavenou optickým filtrem propouštějícím světlo o vlnové délce 570 nm. Tím dojde k odfiltrování nežádoucích částic ze záznamu (nečistoty, bublinky). Tabulka 3. – Umístění a velikosti oblastí měření. Oblast
r1 (r1‘) [mm]
h1 (h1‘)[mm]
s1(s1‘) [mm]
v1(v1‘) [mm]
1
55
88
45
30
1a
55
88
46
30
2
10
55
43
27
3. Výsledky V míchané nádobě byly naměřeny lokální rychlosti kapaliny na celém poloměru nádoby ve výšce 15 mm pod spodní hranou míchadla CVS 691020. Vyhodnocení průběhů střední rychlosti a střední hodnoty čtverce odchylek fluktuační rychlosti s použitím těchto dat již dříve provedli a publikovali autoři Fořt, Kysela a Jirout v [7]. Průběhy výše uvedených rychlostí v této práci jsou rekonstrukcí výsledků v [7]. Rychlosti byly vyjádřeny v závislosti na radiální poloze v bezrozměrném tvaru. Obrázek 4 znázorňuje bezrozměrný profil střední axiální rychlosti dle vzorce (6) v závislosti na radiální poloze. Pro představu je v levém horním rohu obrázku znázorněna lopatka míchadla.
__ __
u*=
u p ×n×D
(6)
Obr. 4. Střední rychlost v poloze z = -15mm. Z obrázku 4 je zřejmé, že průběh bezrozměrných středních axiálních rychlostí dle vztahu (6) je pro všechny otáčky téměř identický. Profil středních axiálních rychlostí lze tedy popsat bezrozměrnou závislostí dle vztahu (6), což souhlasí s dříve publikovanými výsledky v [7]. Na obrázku 4 je také patrná změna směru proudění na bezrozměrném poloměru R = 0,75. Nalevo od tohoto poloměru proudí kapalina od míchadla směrem ke dnu. Vpravo je pak vzestupný proud od dna k hladině. V levé části obrázku můžeme vidět kuželovou oblast reverzního vzestupného proudění pod míchadlem. Tato oblast se vyznačuje nízkou rychlostí proudění opačného směru než proud v oblasti výtoku z míchadla. Rozhraní mezi oběma oblastmi je přibližně na bezrozměrném poloměru R = 0,08.
Obr. 5. Střední hodnota čtverce fluktuační rychlosti v poloze z = -15mm.
Obrázek 5 znázorňuje bezrozměrný profil střední hodnoty čtverce fluktuační axiální rychlosti dle vzorce (7) v závislosti na radiální poloze. u rms * =
u rms p ×n× D
(7)
Z obrázku 5 je zřejmé, že rovněž průběh bezrozměrných středních hodnot čtverce axiálních fluktuačních rychlostí dle vztahu (7) je pro všechny otáčky téměř identický. Profil středních hodnot čtverce axiálních fluktuačních rychlostí lze popsat bezrozměrnou závislostí dle vztahu (7). Také tento výsledek souhlasí se [7]. Oblast přibližně konstantní velikosti urms na obrázku 5 leží v rozsahu bezrozměrných poloměrů R = 0,45 – 0,95. Tato oblast navazuje na zónu vysokých fluktuací pod míchadlem a končí zónou, která je ovlivněna přítomností stěny nádoby. Velikost urms* se v této oblasti pohybuje v blízkosti hodnoty 0,08.
Obrázek 6 znázorňuje profil intenzity turbulence dle vzorce (4) v závislosti na radiální poloze. Tento profil je výrazně ovlivněn změnami směru proudění. Vhodným bodem pro určení parametrů turbulence ve vzestupném proudu je bezrozměrný poloměr R = 0,94. Rychlosti __
u z a urms a také intenzita turbulence jsou zde konstantní. Tabulka 4 ukazuje tyto parametry v daném bodě pro různé otáčky míchadla.
Obr. 6. Intenzita turbulence v poloze z = -15mm. Tabulka 4. – Parametry turbulence v bodě R = 0,94. _
n [1/min] ɛ [W/m3] u ax [m/s] 300 450 600
100 338 800
0,489 0,859 1,141
_
urms [m/s] 0,1122 0,1776 0,2602
urms /
u ax
0,229 0,207 0,228
Výkonová spektrální hustota (PSD) dat naměřených metodou LDA umožnuje využití Kolmogorovovy teorie pro určení vlastností turbulence. Pro výpočet PSD byla spočtena Fourierova transformace signálu a její konvoluce byla zobrazena v závislosti na vlnovém čísle (8) v logaritmických souřadnicích. 2 ×p × f k = __ (8) u ax
Pro vyhlazení průběhu byl signál rozdělen na části a spektrum získané na částech signálu bylo zprůměrováno. Obrázek 7 ukazuje spektrum získané ze 128 vzorků po 2048 bodech. Zobrazený průběh je průměrnou hodnotou z těchto 128 spekter.
Obr. 7. Spektrum fluktuací lokálních axiálních rychlostí kapaliny v míchané nádobě.
Z obrázku 7 je patrné že část závislosti se chová dle vztahu (9). Tuto část proto můžeme považovat za inerciální oblast. PSD = E (k i ) » k i-5 / 3
(9)
Výsledky Fourierovy transformace časové závislosti fluktuací axiálních rychlostí kapaliny budou dále použity pro aproximaci této závislosti.
4. Závěr V míchané nádobě byly naměřeny lokální axiální rychlosti kapaliny na celém poloměru nádoby ve výšce 15 mm pod spodní hranou míchadla CVS 691020. Z naměřených hodnot byla ověřena bezrozměrná závislost střední axiální rychlosti na radiální poloze dle vzorce (6). Z obrázku 4 je zřejmé, že tato závislost je pro všechny otáčky téměř identická. Profil středních axiálních rychlostí tedy lze popsat bezrozměrnou závislostí dle vztahu (6). Profil středních rychlostí rozděluje poloměr nádoby na tři části. Dvě oblasti vzestupného proudu po okrajích poloměru a oblast proudění od míchadla ke dnu uprostřed poloměru. Současně bylo zjištěno, že rovněž profil středních hodnot čtverce axiálních fluktuačních rychlostí urms lze popsat bezrozměrnou závislostí dle vztahu (7). Přičemž je vidět, že tento profil nejdříve v zóně míchadla roste s poloměrem. Na konci této zóny profil prudce klesá na hodnotu 0,08 a dále je přibližně konstantní až do oblasti ovlivněné stěnou nádoby. Profily rychlostí jsou v souladu s výsledky publikovanými dříve jinými autory v [7]. Spektra fluktuací lokálních axiálních rychlostí kapaliny byla vyhodnocena pomocí Fourierovy transformace (FT). Z grafu jejich závislosti na vlnovém čísle byla identifikována inerciální oblast, ve které energie fluktuací závisí pouze na vlnovém čísle dle vztahu (9). Rozklad časové závislosti fluktuačních rychlostí kapaliny na sinusové vlny pomocí FT dále umožní aproximaci původní časové závislosti. Tato aproximace bude součtem jednotlivých sinusovek. Dominantní frekvence, určené z energií jednotlivých vln budou využity k vytvoření aproximace se sníženým počtem koeficientů.
Poděkování Tato práce byla podpořena Grantovou agenturou České republiky prostřednictvím grantu č. P101-12-2274. Tato práce byla podpořena grantem Studentské grantové soutěže ČVUT číslo SGS13/066/OHK2/1T/12. Autor děkuje Ing. Bohuši Kyselovi, Ph.D. z Ústavu pro hydrodynamiku AV ČR za poskytnutí dat z LDA měření a Ing. Michalu Kotkovi, Ph.D. a Ing. Darině Jašíkové, Ph.D. z Fakulty mechatroniky, informatiky a mezioborových studií TU Liberec za poskytnutí dat z PIV měření. Seznam symbolů b šířka narážek C amplituda fluktuací rychlosti D průměr nádoby d průměr míchadla E výkonová spektrální hustota signálu f frekvence H výška hladiny h2 výška míchadla nad dnem nádoby I intenzita turbulence N počet vzorků záznamu n otáčky míchadla t čas urms root mean square velocity uz axiální rychlost
(m) (m/s) (m) (m) (m2/s2) (Hz) (m) (m) (1) (1) (s-1) (s) (m/s) (m/s)
_
uz
střední hodnota axiální rychlosti
(m/s)
uz ‘
fluktuace axiální rychlosti
(m/s)
Řecká písmena φo
k r ν
počáteční fáze vlnové číslo hustota kinematická viskozita
(rad) (1/m) (kg·m-3) (m2/s)
Seznam použité literatury [1] [2] [3] [4] [5]
[6]
[7]
ŠESTÁK, Jiří a František RIEGER. Přenos hybnosti, tepla a hmoty. Vyd. 3. Praha: ČVUT, 2004. 299 s. ISBN 80-01-02933-6. Brucato, A., Grisafi, F., Montante, G.: Particle drag coefficients in turbulent fluids, Chemical Engineering Science, 1998, vol. 53, Issue: 18, p. 3295-3314. Žitný, Rudolf. NAZ. Přednáška 3. Praha: ČVUT, 2012 ZAPLATÍLEK, Karel a Bohuslav DOŇAR. MATLAB: začínáme se signály. 1. vyd. Praha: BEN - technická literatura, 2006. 271 s. ISBN 80-7300-200-0. Benedict, L.H., Nobach, H., Tropea, C.: Estimation of turbulent velocity spectra from laser Doppler data, MEASUREMENT SCIENCE & TECHNOLOGY, 2000, vol. 11, Issue: 8, p. 1089-1104. Roussinova,V., Kresta, S., Weetman, R.. Low frequency macroinstabilities in a stirred tank: scale-up and prediction based on large eddy simulations, Chemical Engineering Science, Volume 58, Issue 11, June 2003, Pages 2297-2311, ISSN 00092509, 10.1016/S0009-2509(03)00097-6. Fort, I., Kysela, B., Jirout, T.: Flow characteristics of axial high speed impellers, Chemical and Process Engineering, 2010, vol. 31, pp 661-679