DIFUSI PERMANEN SATU DIMENSI
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
Diskritisasi Persamaan Difusi Permanen Satu Dimensi dengan Metode Volume Hingga Istiarto – JTSL FT UGM
Persamaan Transpor Persamaan aliran turbulen (persamaan Reynolds), yang terdiri dari persamaan konservasi momentum dan konservasi massa (kekontinuan massa), adalah:
∂u ∂uu ∂uv ∂uw ∂ ⎛ ∂u ⎞ ∂ ⎛ ∂u ⎞ ∂ ⎛ ∂u ⎞ ⎟ − ⎜ ν t + + + − ⎜ ν t ⎟ − ⎜ ν t ⎟ = ∂t ∂x ∂y ∂z ∂x ⎝ ∂x ⎠ ∂y ⎜⎝ ∂y ⎟⎠ ∂z ⎝ ∂z ⎠ 1 ∂p 2 − − k + gx ρ ∂x 3 ∂v ∂vu ∂vv ∂vw ∂ ⎛ ∂v ⎞ ∂ ⎛ ∂v ⎞ ∂ ⎛ ∂v ⎞ + + + − ⎜ ν t ⎟ − ⎜ ν t ⎟ − ⎜ ν t ⎟ = ∂t ∂x ∂y ∂z ∂x ⎝ ∂x ⎠ ∂y ⎜⎝ ∂y ⎟⎠ ∂z ⎝ ∂z ⎠ 1 ∂p 2 − − k + gy ρ ∂y 3 ∂w ∂wu ∂wv ∂ww ∂ ⎛ ∂w ⎞ ∂ ⎛ ∂w ⎞ ∂ ⎛ ∂w ⎞ ⎟ − ⎜ ν t + + + − ⎜ ν t ⎟ − ⎜ ν t ⎟ = ∂t ∂x ∂y ∂z ∂x ⎝ ∂x ⎠ ∂y ⎜⎝ ∂y ⎟⎠ ∂z ⎝ ∂z ⎠ 1 ∂p 2 − − k + gz ρ ∂z 3 ∂u ∂v ∂w + + =0 ∂x ∂y ∂z
(1)
(2)
(3) (4)
Persamaan-persamaan di atas dapat dibaca sebagai persamaan transpor konvektifdifusif. Pada Pers. 1 s.d. 3, empat suku pertama di sisi kiri adalah suku-suku konveksi dan tiga suku berikutnya adalah suku-suku difusi. Suku-suku di sisi kanan dapat dibaca sebagai source. Persamaan-persamaan tersebut dapat dituliskan dalam bentuk persamaan transpor besaran skalar φ sebagai berikut: "# "# ∂φ ℓ + ∇ φ ℓ V − ∇ Γ ℓ ∇φ ℓ = R ℓ (5) ∂t !" Dalam persamaan di atas φ ℓ adalah suatu besaran skalar, V adalah vektor kecepatan
( ) (
)
aliran, Γ ℓ adalah koefisien difusi, dan Rℓ adalah source berupa besaran skalar. Definisi setiap suku pada persamaan transpor tersebut apabila dikaitkan dengan persamaan Reynolds disajikan pada Tabel 1.
1
Tabel 1. Definisi suku-suku pada persamaan transpor yang dikaitkan dengan persamaan aliran turbulen (Persamaan Reynolds).
"# "# ∂φ ℓ + ∇ φℓ V − ∇ Γℓ ∇ φℓ = R ℓ ∂t Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
( ) (
ℓ
φℓ
Γℓ
1
u
νt
2
v
νt
3
w
νt
4
1
0
)
(5)
Rℓ −
1 ∂p 2 − k + gx ρ ∂x 3
−
1 ∂p 2 − k + gy ρ ∂y 3
−
1 ∂p 2 − k + gz ρ ∂z 3 0
Difusi Permanen Satu Dimensi Proses transpor besaran skalar φ yang paling sederhana adalah difusi permanen (steadystate pure diffusion). Persamaan transpor difusi permanen dapat dengan mudah diperoleh dari persamaan umum transpor, Pers. 5, dengan menghilangkan suku derivatif terhadap waktu dan suku konvektif. Dengan cara ini, diperoleh: "# ∇ Γ ℓ ∇φ ℓ = R ℓ (6)
(
)
Pada persamaan di atas, tanda negatif di depan suku difusif telah dihilangkan untuk menyederhanakan penulisan. Dalam hal ini, suku source pada persamaan di atas berkebalikan tanda dengan suku source pada Pers. 5. Untuk difusi di bidang satu dimensi, persamaan di atas dapat dituliskan menjadi: d # dφ& %Γ (= R dx$ dx'
(7)
Dalam persamaan tersebut, subskrip ℓ tidak diperlukan sehingga tidak dituliskan. Aplikasi persamaan difusi permanen satu dimensi dipaparkan melalui contoh kasus di bawah ini. Suatu saluran (flume) tampang segi empat diisi air dan kedua ujung saluran ditutup sehingga tidak ada aliran dalam saluran tersebut. Panjang saluran 10 m, lebar saluran 0.40 m, dan kedalaman air 0.20 m. Air di ujung kiri saluran diberi temperatur 30°C dan di ujung kanan saluran diberi temperatur 25°C. Keduanya dipertahankan konstan. Temperatur di sepanjang saluran bervariasi antara kedua nilai tersebut, mengikuti gerak difusi temperatur. Koefisien difusi 5 m2/s. Temperatur di sepanjang saluran dapat dihitung dengan memakai Pers. 7; dalam hal ini, dipakai metode volume hingga (finite volume).
2
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
Coding dan hitungan dilakukan dengan bantuan spreadsheet. Para mahasiswa peserta kuliah diminta mengerjakan sendiri dengan mengacu pada langkah-langkah yang dijabarkan pada papan tulis dan tayangan yang ditunjukkan di ruang kuliah. Mahasiswa tidak disarankan untuk hanya membaca penjabaran langkah dan meng-copy file spreadsheet. Langkah #1: Pembuatan grid (diskritisasi domain model) Domain model, dalam hal ini saluran, dibagi menjadi sejumlah volume kontrol (control volume) diskrit, atau sering pula disebut dengan cell. Volume kontrol di dekat batas domain model diletakkan sedemikian hingga sisi volume kontrol berimpit dengan batas domain. Titik hitung (node) adalah titik tengah volume kontrol. Dalam kasus difusi di dalam saluran ini, domain model dibagi menjadi 10 volume kontrol berukuran seragam, panjang 1 m, agar hitungan menjadi sederhana.
Gambar 1. Diskritisasi domain model menjadi sejumlah volume kontrol. Langkah #2: Diskritisasi persamaan difusi Integrasi persamaan diferensial difusi, Pers. 7, untuk sebuah volume kontrol adalah: d # dφ&
(8)
∫∫∫ cv d x %$Γ d x (' d ϑ = ∫∫∫ cv Rd ϑ
Dengan memakai teorema divergensi Gauss, integral volume suku difusi di sisi kiri pada persamaan di atas dapat diubah menjadi integral luasan tertutup yang menyelimuti volume kontrol: dφ ! ! Γ (9) ∫∫ S d x n ⋅ d S = ∫∫∫ cv Rd ϑ " Dalam persamaan tersebut, S adalah vektor luas tegak lurus ke arah luar dari selimut volume kontrol. Untuk volume kontrol satu dimensi (Gambar 2), persamaan di atas dapat diubah menjadi: # dφ& # dφ& %Γ S ( − %Γ S ( = R Δϑ d x 'e $ d x 'w $
(10)
Dalam persamaan di atas, S adalah luas tampang sisi volume kontrol, Δϑ adalah volume volume-kontrol, dan R adalah source rata-rata di dalam volume kontrol. Titik hitung (node) suatu volume kontrol diidentifikasikan dengan simbol P yang memiliki volume kontrol tetangga di sisi kiri W (West) dan sisi kanan E (East). Sisi volume kontrol di kiri adalah w dan di kanan adalah e.
3
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
Gambar 2. Volume kontrol satu dimensi. Nilai koefisien difusi Γ di sisi-sisi volume kontrol, w dan e, dihitung dengan interpolasi linear dari nilai-nilai yang ada di dua node volume kontrol di kiri dan kanannya:
Γw =
ΔxwP ΓW + ΔxWw ΓP ΔxWw + ΔxwP
dan Γe =
ΔxeE ΓP + ΔxPe ΓE ΔxPe + ΔxeE
(11)
Cara interpolasi linear untuk menghitung nilai-nilai di sisi volume kontrol di atas dikenal pula dengan istilah beda tengah (central difference). Apabila ukuran volume kontrol seragam, maka:
Γw = 12 (ΓW + ΓP ) dan Γe = 12 (ΓP + ΓE )
(12)
Suku yang menunjukkan fluks difusi (suku-suku di sebelah kiri pada Pers. 10) dihitung dengan cara sebagai berikut:
#φ −φ & #φ −φ # dφ& # dφ& E P ( dan P W % Γ S = Γ S Γ S = Γ S % ( % ( ( e e% w w %% d x 'e d x 'w $ $ $ Δx PE ' $ ΔxWP
& (( '
(13)
Suku source di sebelah kanan pada Pers. 10, pada beberapa kasus, sering merupakan fungsi variabel φ itu sendiri. Dalam hal ini, suku source dihitung dengan cara linear sebagai berikut:
R Δϑ = Ru + R p φ P
(14)
Substitusi Pers. 13 dan 14 ke Pers. 10 menghasilkan persamaan berikut ini:
%φ −φ ( % P * − Γ S ' φ P − φW Γ e Se '' E * w w' & Δx PE ) & ΔxWP
( ** = Ru + R p φ P )
yang dengan pengelompokan koefisien menghasilkan persamaan berikut ini:
#Γ S & # Γ S Γ S & #Γ S & %% e e (( φ E + %%− e e − w w − R p (( ϕφ + %% w w (( φW = Ru $ Δx PE ' $ Δx PE ΔxWP ' $ ΔxWP '
(15)
Apabila koefisien-koefisien suku φE, φW, dan φP dinamai aE, aW, dan aP, maka persamaan di atas dapat dituliskan sebagai berikut: (16)
aW φW + a P φ P + a E φ E = Ru
Koefisien-koefisien dalam persamaan tersebut adalah:
4
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
aW =
Γw S w Γ S , a E = e e , dan a P = −aW − a P − R p ΔxWP ΔxPE
Pers. 16 merupakan persamaan diskrit di setiap volume kontrol. Untuk 10 volume kontrol (lihat Gambar 1), akan diperoleh 10 persamaan. Nilai variabel φ di setiap volume kontrol (φP) merupakan fungsi nilai-nilai φ di 2 volume kontrol tetangga (φE dan φW). Volume kontrol yang berada di kedua ujung saluran, yaitu volume kontrol 1 dan 10, hanya memiliki satu tetangga, sedangkan salah satu sisi berimpit dengan batas domain model. Di batas domain, nilai φ diketahui. Nilai ini dikenal sebagai syarat batas (boundary condition). Di bawah ini dipaparkan penyusunan kesepuluh persamaan tersebut. Volume kontrol 2 s.d. 9 Dengan 10 volume kontrol seragam, Δx = 1 m, koefisien difusi seragam, Γ = 5 m2/s, dan luas sisi volume kontrol seragam, S = 0.40×0.20 = 0.08 m2, maka koefisien-koefisien pada persamaan diskrit difusi adalah sebagai berikut:
aE =
Γe Se 5 × 0.08 = = 0.40 ΔxPE 1
aW =
Γw S w 5× 0.08 = = 0.40 ΔxWP 1
a P = −aW − a E − R p = −0.40 − 0.40 − 0 = −0.80 Perlu dicatat bahwa dalam kasus ini tidak ada source, sehingga Ru = 0 dan Rp φP = 0. Volume kontrol 1 Sisi kiri volume kontrol ini, w, berimpit dengan batas domain. Di sebelah kanan (timur), volume kontrol ini bertetangga dengan volume kontrol 2. Koefisien aE dihitung seperti hitungan pada volume kontrol 2 s.d. 9. Di sebelah kiri (barat), tidak ada volume kontrol tetangga, sehingga koefisien aW tidak ada. Koefisien aP dan suku source dihitung dengan cara berbeda dari cara sebelumnya untuk volume kontrol 2 s.d. 9. Persamaan diskrit transpor difusif di volume kontrol 1 adalah sebagai berikut (lihat Gambar 3):
Gambar 3. Volume kontrol di batas kiri domain model.
%φ −φ ( % ( P * − Γ S ' φP − φw * = 0 Γ e Se '' E * w w ' Δx * & Δx PE ) & wP ) Dalam hal ini, φw = φA = 30°C.
5
#Γ S & # Γ S Γ S & Γ S %% e e (( φ E + %%− e e − w w (( φ P = − w w φ w ΔxwP $ Δx PE ' $ Δx PE ΔxwP '
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
aW = 0 aE =
Γe Se 5 × 0.08 = = 0.40 ΔxPE 1
Rp =
Γ w S w 5× 0.08 = = 0.80 Δx wP 12
a P = −aW − a E − R p = −0 − 0.40 − 0.80 = −1.20 Ru = −
Γ w Sw 5× 0.08 φw = − 30 = −24 Δx wP 12
Volume kontrol 10 Perlakuan untuk volume kontrol 10 ini mirip dengan volume kontrol 1, hanya saja sisi yang berimpit dengan batas domain model adalah sisi kanan, e.
Gambar 4. Volume kontrol di batas kanan domain model. Persamaan diskrit transpor difusif di volume kontrol 10 adalah sebagai berikut (lihat Gambar 4):
%φ −φ ( % P * − Γ S ' φ P − φW Γ e Se '' e * w w' & Δx Pe ) & ΔxWP
( ** = 0 )
Dalam hal ini, φe = φB = 25°C.
#Γ S & %% w w (( φW $ ΔxWP ' aW =
# Γ S Γ S & Γ S + %%− w w − e e (( φ P = − e e φ e Δx Pe $ ΔxWP Δx Pe '
Γw S w 5 × 0.08 = = 0.40 ΔxWP 1
aE = 0 Rp =
Γ e Se 5× 0.08 = = 0.80 Δx Pe 12
6
a P = −aW − a E − R p = −0.40 − 0 − 0.80 = −1.20
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
Ru = −
Γe Se 5 × 0.08 φe = − 25 = −20 ΔxPe 12
Langkah #3: Penyelesaian persamaan (sistem persamaan linear) Kesepuluh persamaan yang diperoleh dari penjabaran persamaan diskrit transpor difusif di setiap volume kontrol pada langkah ke-2 merupakan satu sistem persamaan linear. Sistem persamaan ini dapat dituliskan dalam bentuk perkalian matriks sebagai berikut:
AΦ=R
(17)
A adalah matriks bujur sangkar berdimensi 10×10 yang elemennya adalah koefisien aE, aW, dan aP, Φ adalah vektor kolom yang elemennya adalah variabel φ, dan R adalah vektor kolom yang elemennya adalah konstanta/nilai source.
! # # # # # # # # # # # # # # # #"
aP
aE
0
0
0
0
0
0
0
aW
aP
aE
0
0
0
0
0
0
0
aW
aP
aE
0
0
0
0
0
0
0
aW
aP
aE
0
0
0
0
0
0
0
aW
aP
aE
0
0
0
0
0
0
0
aW
aP
aE
0
0
0
0
0
0
0
aW
aP
aE
0
0
0
0
0
0
0
aW
aP
aE
0
0
0
0
0
0
0
aW
aP
0
0
0
0
0
0
0
0
aW
0 $( &* 0 &* &* 0 &* &* 0 &* 0 & ** &) 0 &* &* 0 &* & 0 &* * aE & * &* a P &% *+
φ1 , ( * * φ2 * * * * φ3 * * * * φ4 * * φ5 ** ** - ) φ6 * = * * * φ7 * * φ8 * * * * φ9 * * * * φ10 *. *+
Ru1 , * Ru2 * * Ru3 * * Ru4 * Ru5 ** Ru6 * * Ru7 * Ru8 * * Ru9 * * Ru10 *.
!########## #"########### $! #"# $ !#"# $ A
Φ
R
Substitusi nilai-nilai koefisien yang telah diperoleh pada langkah ketiga ke dalam matriks di atas menghasilkan matriks sebagai berikut:
7
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
) + " −1.2 0.4 % + 0 0 0 0 0 0 0 0 $ '+ 0 0 0 0 0 0 0 '+ $ 0.4 −0.8 0.4 $ 0 0.4 −0.8 0.4 0 0 0 0 0 0 '+ $ ' 0 0.4 −0.8 0.4 0 0 0 0 0 '+ $ 0 + $ 0 0 0 0.4 −0.8 0.4 0 0 0 0 '+ * $ 0 0 0 0 0.4 −0.8 0.4 0 0 0 '+ $ ' 0 0 0 0 0.4 −0.8 0.4 0 0 '+ $ 0 + $ 0 0 0 0 0 0 0.4 −0.8 0.4 0 '+ $ ' 0 0 0 0 0 0 0.4 −0.8 0.4 ' + $ 0 $# 0 0 0 0 0 0 0 0 0.4 −1.2 '& + + +,
φ1 + φ 2 + ) −24 + + + φ3 + + 0 + + + 0 + φ4 + + + + 0 + + φ5 + + 0 + .=* . φ6 + + 0 + + + + φ7 + + 0 + 0 + φ8 + + + + 0 + φ9 + +, −20 +/ + φ10 +/
Matriks A adalah matriks diagonal dengan lebar (bandwidth) 3. Jika dimensi matriks besar, penyelesaian persamaan di atas umumnya dilakukan dengan tri-diagonal matriks algorithm (TDMA). Jika dimensi matriks kecil, penyelesaian persamaan tersebut dapat dilakukan dengan cara langsung, misal dengan matriks inversi atau cara iterasi. Didalam MSExcel, ada sejumlah fungsi yang telah disediakan untuk melakukan operasi matriks, misal untuk menghitung determinan, mencari inversi matriks, atau melakukan perkalian matriks. Perkalian inversi matriks A (A 1) dengan suku di sisi kiri dan kanan Pers. 17 menghasilkan persamaan berikut: −
Φ = A −1 R
(18)
Langkah penyelesaian persamaan di atas dengan bantuan MSExcel dipaparkan pada paragraf-paragraf di bawah ini. 1. Misalkan elemen-elemen matriks A dituliskan dalam sel B2:K11 dan elemen matriks R dalam sel N2:N11. 2. Matriks A 1 dicari dan disimpan menjadi sel B14:K23 dengan langkah sebagai berikut: −
− pilih sel B14:K23 − tuliskan fungsi untuk menghitung invers matriks =MINVERSE(B2:K11) − tekan tombol CNTRL+SHIFT+ENTER bersama-sama − sel B14:K23 berisi A 1 seperti di bawah ini −
1.1875 1.0625 0.9375 0.8125 0.6875 0.5625 0.4375 0.3125 0.1875 0.0625
1.0625 3.1875 2.8125 2.4375 2.0625 1.6875 1.3125 0.9375 0.5625 0.1875
0.9375 2.8125 4.6875 4.0625 3.4375 2.8125 2.1875 1.5625 0.9375 0.3125
0.8125 2.4375 4.0625 5.6875 4.8125 3.9375 3.0625 2.1875 1.3125 0.4375
0.6875 2.0625 3.4375 4.8125 6.1875 5.0625 3.9375 2.8125 1.6875 0.5625
8
0.5625 1.6875 2.8125 3.9375 5.0625 6.1875 4.8125 3.4375 2.0625 0.6875
0.4375 1.3125 2.1875 3.0625 3.9375 4.8125 5.6875 4.0625 2.4375 0.8125
0.3125 0.9375 1.5625 2.1875 2.8125 3.4375 4.0625 4.6875 2.8125 0.9375
0.1875 0.5625 0.9375 1.3125 1.6875 2.0625 2.4375 2.8125 3.1875 1.0625
0.0625 0.1875 0.3125 0.4375 0.5625 0.6875 0.8125 0.9375 1.0625 1.1875
3. Vektor kolom Φ dihitung dengan fungsi yang telah disediakan dalam MSExcel dan disimpan menjadi sel N26:N35 dengan langkah sebagai berikut: − pilih sel N26:N35
− tekan tombol CNTRL+SHIFT+ENTER bersama-sama − sel N26:N35 berisi Φ seperti di bawah ini 29.75 29.25 28.75 28.25 27.75 27.25 26.75 26.25 25.75 25.25 31 30 Temperatur [° C]
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
− tuliskan fungsi untuk menghitung perkalian dua matriks =MMULT(B14:K23,N26:N35)
29 28 27 26 25 24 0
1
2
3
4
5
6
7
8
9
10
Jarak [m]
Gambar 5. Profil temperatur air di flume yang diperoleh dari penyelesaian numeris persamaan difusi.
Penutup Sampai di sini, para pembaca atau mahasiswa diharapkan telah memahami proses dan mekanisme transpor difusif suatu besaran skalar serta mampu melakukan hitungan untuk kasus sederhana, yaitu transpor difusif permanen satu dimensi. Satu hal yang perlu dicamkan oleh para mahasiswa adalah bahwa sangat penting untuk melakukan sendiri penjabaran dan diskritisasi persamaan transpor difusif dan selanjutnya
9
melakukan sendiri coding langkah-langkah penyelesaian aplikasi transpor difusif pada contoh kasus yang telah diberikan.
Istiarto Jurusan Teknik Sipil dan Lingkungan FT UGM http://istiarto.staff.ugm.ac.id email:
[email protected]
Topik berikutnya adalah transpor konvektif-difusif untuk kasus yang juga masih sederhana, yaitu transpor konvektif-difusif permanen satu dimensi.
Referensi Ferziger, J. H., and Peric, M., 1997, Computational Methods for Fluid Dynamics, Springer-Verlag, Berlin, Germany. Istiarto, I., 2001, Flow Around A Cylinder In A Scoured channel Bed, Doctoral Dissertation, EPFL, Switzerland. Versteeg, H. K., and Malalasekera, W., 1995, An Introduction to Computational Fluid Dynamics, The Finite Volume Method, Longman Group, Essex, England.
-o0o-
10