Software
2D Lid-Driven Cavity v.1.0.0
Strujanje u šupljini sa pokretnim poklopcem
Namjena
2D Lid-Driven Cavity je edukativni program za simulaciju strujanja u kvadratnoj ili pravougaonoj šupljini čiji gornji zid (poklopac) klizi konstantnom brzinom i pokreće fluid unutra. To je klasičan test-problem mehanike fluida: geometrija je jednostavna, a u strujanju se ipak javljaju glavni vrtlog, sekundarni vrtlozi u uglovima i tanki granični slojevi uz zidove.
Program rješava Navier-Stokesove jednačine metodom konačnih razlika i prikazuje strujnice, vrtložnost, brzinu i pritisak, kao i profile brzine kroz sredinu šupljine. Reynoldsov broj se zadaje direktno ili se računa iz izabranog fluida i dimenzija.
Rezultati se poklapaju sa referentnim rješenjem iz literature (Ghia, Ghia i Shin, 1982) za Re = 100, 400 i 1000. Program nije zamjena za komercijalne CFD pakete — namijenjen je učenju i demonstraciji numeričkih metoda.
Fizički problem
Šupljina je ispunjena nestišljivim fluidom. Poklopac klizi brzinom U, a fluid uz njega se zbog viskoznosti povlači za njim. Uz desni zid fluid skreće naniže, vraća se po dnu i penje uz lijevi zid — tako nastaje glavni vrtlog. U donjim uglovima formiraju se manji vrtlozi koji se okreću u suprotnom smjeru.
Cijelo strujanje određuje jedan broj, Reynoldsov:
Re = U · L / ν
gdje je L širina šupljine, a ν kinematička viskoznost fluida. Pri malom Re dominira viskoznost i vrtlog je pomjeren prema poklopcu; kako Re raste, vrtlog se spušta ka sredini, a vrtlozi u uglovima rastu. Program računa u bezdimenzionalnom obliku, a rezultate prikazuje i izvozi u SI jedinicama.
Kako program računa
Umjesto brzine i pritiska, program računa dvije pomoćne veličine: strujnu funkciju ψ, čije su izolinije strujnice, i vrtložnost ω, koja opisuje lokalno obrtanje fluida:
∇2ψ = −ω
∂ω/∂t + u ∂ω/∂x + v ∂ω/∂y = (1/Re) ∇2ω
Prva jednačina iz vrtložnosti daje strujnu funkciju, a iz nje brzine. Druga opisuje kako vrtložnost putuje sa fluidom i širi se zbog viskoznosti. Prednost ovog pristupa je što se pritisak u toku proračuna uopšte ne pojavljuje; računa se tek na kraju, iz gotovog polja brzine.
Vrtložnost na zidovima nije poznata unaprijed, već se računa iz strujne funkcije u susjednim čvorovima (Thomova formula). Upravo tu poklopac “ubacuje” vrtložnost u fluid.
Numeričke metode
Mreža je ravnomjerna, a izvodi se zamjenjuju razlikama između susjednih čvorova. Za član koji opisuje prenos vrtložnosti strujom bira se jedna od tri šeme:
| Šema | Tačnost | Kada je koristiti |
|---|---|---|
| Centralna | drugi red | Za tačne rezultate; traži manji vremenski korak pri većem Re |
| Upwind | prvi red | Uvijek stabilna, ali “razmazuje” rješenje — kao da je Re manji nego što jeste |
| Hybrid | između | Kombinacija: centralna gdje je mreža dovoljno gusta, upwind gdje nije |
Jednačina za strujnu funkciju rješava se u svakom vremenskom koraku. Podrazumijevano se koristi direktni rješavač: matrica se rastavi (LU) jednom na početku, pa je svaki naredni korak brz. Može se izabrati i iterativni SOR, na kojem se vidi kako relaksacioni faktor utiče na broj iteracija.
Vrijeme se integriše eksplicitno, pa vremenski korak ne smije biti prevelik. Program prije pokretanja računa najveći dozvoljeni korak iz tri uslova (difuzija, CFL i, za centralnu šemu, uslov ćelijskog Pecletovog broja) i prikazuje ih na ekranu. Ako je zadati korak veći, program ga sam smanji i o tome obavijesti. Simulacija se završava kada se strujanje ustali ili kada istekne zadato vrijeme.
Mogućnosti programa
- Reynoldsov broj — zadaje se direktno, ili se računa iz izabranog fluida (vazduh, voda, ulje SAE 30, glicerin), brzine poklopca i širine šupljine
- Pokretni zidovi — može se pokrenuti bilo koji zid, ne samo poklopac, u bilo kom smjeru
- Geometrija — kvadratna ili pravougaona šupljina, sa gustinom mreže po izboru
- Prikaz — strujna funkcija, vrtložnost, brzina ili pritisak, sa strujnicama, vektorima brzine i mrežom po izboru
- Profili brzine — kroz vertikalnu i horizontalnu sredinu šupljine, uz referentne tačke iz literature za Re = 100, 400 i 1000
- Očitavanje mišem — koordinate i brzina u tački ispod kursora
- Rad u pozadini — prozor ostaje upotrebljiv dok traje proračun, a simulacija se može prekinuti u svakom trenutku
- Izvoz — brzine, strujna funkcija, vrtložnost i pritisak u SI jedinicama, kao CSV ili NumPy (.npz)
Provjera rezultata
Kao referenca korišten je rad U. Ghia, K. N. Ghia i C. T. Shin, High-Re solutions for incompressible flow using the Navier-Stokes equations and a multigrid method, Journal of Computational Physics 48 (1982) 387–411. To je vjerovatno najčešće korišteno poređenje za ovaj problem: autori su ga riješili na mreži 129 × 129 i objavili brzine u tačkama kroz sredinu šupljine.
Na istoj mreži razlika između programa i Ghia rezultata, izražena u odnosu na brzinu poklopca, iznosi:
| Re | najveća razlika u brzini u | najveća razlika u brzini v |
|---|---|---|
| 100 | 0,4 % | 0,8 % |
| 400 | 0,5 % | 0,4 %* |
| 1000 | 1,4 % | 0,8 % |
* Bez jedne tačke iz originalne tabele (Re = 400, x = 0,9063), koja se vidi i na slici iznad. Za nju je objavljeno v = −0,238, a program daje −0,382. U susjednim tačkama slaganje je na trećoj decimali, a objavljena vrijednost ne leži na glatkoj krivoj između njih. Na finijoj mreži rezultat programa se od nje udaljava, a ne približava. Najvjerovatnije je riječ o štamparskoj grešci u originalnom radu; u programu je vrijednost ostavljena onakva kakva je objavljena.
Da li je rješenje tačno ili se samo poklapa?
Ni rješenje iz 1982. godine nije savršeno — i ono ima grešku koja zavisi od gustine mreže. Zato je za Re = 100 program pokrenut na četiri sve gušće mreže, a rezultati su ekstrapolirani na beskonačno gustu mrežu:
| Mreža | najmanja brzina u | najveća brzina v |
|---|---|---|
| 33 × 33 | −0,2042 | 0,1683 |
| 65 × 65 | −0,2116 | 0,1768 |
| 129 × 129 | −0,2135 | 0,1789 |
| 257 × 257 | −0,2139 | 0,1794 |
| ekstrapolirano | −0,2140 | 0,1796 |
Sa svakim prepolovljenjem koraka mreže greška pada četiri puta, kako i treba za metodu drugog reda tačnosti. Ghia vrijednosti se od rješenja na najgušćoj mreži razlikuju za 0,5 do 0,9 % brzine poklopca — isto koliko i rezultati programa na mreži 129 × 129. Drugim riječima, odstupanja u tabeli iznad nisu greška programa, već granica do koje ovakvo poređenje ima smisla.
Uticaj šeme
Pri Re = 1000 razlike između šema postaju očigledne. Centralna šema odstupa od reference oko 1,4 %, a upwind i hybrid oko 16 do 17 %: profil brzine je zaravnjen, kao da je viskoznost veća nego što jeste. Za kvantitativne rezultate treba koristiti centralnu šemu, a ostale dvije dobro služe da se pokaže šta je numerička difuzija.
Ograničenja
- Strujanje je dvodimenzionalno i laminarno, bez turbulencije
- Fluid je nestišljiv, sa konstantnim svojstvima; nema prenosa toplote
- Šupljina je zatvorena — nema ulaza ni izlaza fluida, a zidovi se kreću samo duž sebe
- Mreža je ravnomjerna, bez zgušnjavanja uz zidove
- Zbog eksplicitne vremenske šeme proračun za veći Re na gustoj mreži traje duže (Re = 1000 na mreži 129 × 129: nekoliko minuta)
- Pritisak se računa samo za ustaljeno strujanje
Tehnologije
Python 3.11 · PySide6 · NumPy · SciPy · Matplotlib
Grafički interfejs koristi Qt preko biblioteke PySide6, pod licencom GNU LGPL v3. Qt biblioteke se isporučuju kao zasebni fajlovi, pa ih korisnik može zamijeniti vlastitom verzijom, kako LGPL i zahtijeva. NumPy, SciPy i Matplotlib su pod BSD licencom.
Cijena i sugestije
2D Lid-Driven Cavity je potpuno besplatan, bez ograničenja u funkcionalnosti i bez naknadnih verzija koje se plaćaju. Sve sugestije su dobrodošle — prijedlozi, primjedbe i prijave grešaka.
Kontakt
Za pitanja, prijedloge ili prijave grešaka — direktno na mkozica@outlook.com.
Preuzimanje
Instalacioni program za Windows (64-bitni), veličine 67 MB. Sve potrebne biblioteke su uključene, pa se Python ne mora posebno instalirati. Program se može instalirati i bez administratorskih prava, samo za tekućeg korisnika.