La Gubu et al//Paradigma, Vol. 15 No .2, Oktober 2011, hlm. 69-78
METODE BEDA HINGGA UNTUK ANALISIS PROSES TRANSFER MASSA La Gubu1, Edi Cahyono1, dan La Hamimu2 1
Program Studi Matematika, Jurusan Matematika, FMIPA, Universitas Haluoleo, Kendari 93231 2 Program Studi Fisika, Jurusan Matematika, FMIPA, Universitas Haluoleo, Kendari 93231
ABSTRAK Dalam tulisan akan dibicarakan proses transfer massa dengan menggunakan metode beda hingga yang diimplementasikan dengan software Matlab. Daerah observasi transfer massa dibagi menjadi sejumlah grid 6x6 dengan dimensi grid ∆x = ∆y = 0,1 . Selanjutnya dalam skema numerik metode beda hingga rapat massa pada setiap titik dalam daerah observasi dapat diketahui. Dengan menggunakan input parameter ∆t = 10 −3 dan kecepatan difusi massa c=1 diperoleh hasil bahwa seiring dengan berjalannya waktu t, massa merambat pada semua titik dalam daerah observasi. Ditinjau 2 buah sumber yang melepaskan massa secara kontinu. Hasil simulasi numerik menunjukkan bahwa semakin jauh titik pengamatan dari sumber massa maka waktu yang diperlukan juga semakin lama. Rapat massa di sekitar sumber memiliki nilai yang lebih besar bila dibandingkan dengan titik yang jauh dari sumber. Rapat massa sebagai fungsi waktu diperoleh hasil bahwa nilai rapat massa pada suatu titik tertentu dalam daerah observasi akan semakin bertambah seiring bertambahnya waktu t. Tetapi untuk selang waktu yang mendekati orde ribuan maka rapat massa u menjadi konstan (stabil). Secara alamiah hal ini bisa terjadi karena adanya syarat batas yang diberlakukan serta kriteria parameter input yang diberikan. Kata kunci: metode beda hingga, difusi massa, rapat massa ABSTRACT In this paper we discuss mass transfer process using finite difference method which is implimeted in Matlab. The observed domain is partioned into 6x6 grids with stepsize ∆x = ∆y = 0,1 . Moreover, in the numerical scheme of the finite difference method the mass density at each point in the domain is able to be obtained. With the parameter input ∆t =10 −3 and mass diffusion c=1, the results show that the mass propagates to whole observed area as time t increases. Considering 2 points as continuous mass sources, the numerical results show that as the distance of observed points and the sources is getting large, the time needed to propagate is getting larger. The mass density is greater around the sources compared to the density far from the sources. As a function of time, the mass density at a certain point in the observed area increases as time t increases. However, for the time interval of thousandth order the mass density u become stable or constant. Scientifically, it could happen as the boundary condition as well as certain criteria for the given parameter are provided. Key words: finite difference method, mass transver, mass density Diterima: 1 Juni 2011 Disetujui untuk dipublikasikan: 1 Agustus 2011
70
Metode Beda Hingga Untuk Analisis Proses Transfer Massa
1. PENDAHULUAN Proses transfer massa, disadari maupun tidak, sering dijumpai dalam kehidupan sehari-hari. Proses larutnya gula di dalam segelas air merupakan transfer massa gula ke seluruh air di dalam gelas tersebut. Mula-mula gula berada di dasar gelas, gula ini akan larut di dalam air, makin lama massa gula ini akan menyebar pada seluruh air di dalam gelas tersebut. Bila air dalam gelas tersebut diaduk, proses larutnya gula ini semakin cepat. Pada masalah yang lebih kompleks, transfer massa dapat terjadi pada proses menyebarnya tumpahan minyak pada peristiwa pecahnya kapal tangker. Titik pada saat mana kapal tangker tersebut pecah dianggap sebagai titik sumber tumpahan minyak. Pada peristiwa ini titik sumber dianggap melepaskan massa hanya sekali. Sumber yang hanya melepaskan massa sekali saja disebut sumber dengan denyut tunggal. Solusi yang dihasilkan dari denyut tunggal ini merupakan solusi fundamental persamaan transfer massa. Gambar 1 menunjukkan difusi massa dengan denyut tunggal berdasarkan solusi fundamental persamaan difusi.
a
b
c
d
Gambar 1. Difusi massa dengan denyut tunggal dengan urutan waktu dari a ke d. Pada saat baru terjadi pelepasan massa, kerapatan massa masih terkonsentrasi di sekitar sumber, massa belum merambat jauh. Hal ini terlihat dari daerah hitam yang terkonsentrasi di sekitar sumber dan terang di titik-titik yang jauh dari sumber (Gambar 1a). Berikutnya, pada Gambar 1b daerah hitam di sekitar sumber mulai membesar, menunjukkan massa telah mulai menyebar. Penyebaran massa ini terus terjadi dan semakin tampak jelas pada Gambar 1c yaitu dengan membesarnya daerah hitam dan juga memudarnya warna hitam di sekitar sumber. Pada Gambar 1d daerah di sekitar sumber
71
La Gubu et al//Paradigma, Vol. 15 No.2, Oktober 2011, hlm. 69-78
bahkan tidak hitam lagi, yang menunjukkan penurunan kerapatan massa yang cukup berarti di titik sumber. Penurunan kerapatan massa di sekitar sumber ini semata-mata karena massa telah menyebar dan secara total tidak terjadi pengurangan ataupun penambahan massa. Contoh yang lain proses difusi massa dapat berupa bocornya pipa pada instalasi pipa di dasar laut atau dalam tanah. Instalasi pipa demikian sering digunakan pada industri kimia. Dalam industri yang demikian, instalasi pipa-pipa sering dipasang di dalam tanah di bawah pabrik. Kebocoran pada instalasi pipa ini akan mengakibatkan transfer massa zat kimia ke daerah sekitar kompleks industri kimia tersebut. Hal ini akan menyebabkan polusi yang sangat berbahaya. Dalam contoh ini titik sumbernya adalah tempat bocornya instalasi pipa tersebut. Cahyono (2004) menggunakan pendekatan sumber berdenyut untuk kasus ini. Pendekatan ini dimaksudkan untuk tujuan ke depannya dalam mencari titik sumber. Pada peristiwa kebocoran seperti ini, titik sumber melepaskan massa secara kontinu, dan bukan berdenyut. 2. PERSAMAAN DIFUSI TRANSFER MASSA Transfer massa yang dimodelkan dalam persamaan difusi banyak diterapkan dalam industri seperti perambatan air dari dalam kayu menuju ke permukaan kayu pada proses pengeringan kayu, Cahyono etal (2003). Persamaan difusi ini juga merupakan model perambatan panas dan telah digunakan untuk mempelajari perambatan panas pada proses sterilisasi makanan dalam kaleng, seperti pada Banga etal (1991) dan Wibowo & Ponidi (2002). Persamaan difusi dalam bentuk yang paling sederhana adalah pada kasus dimensi satu dan diberikan oleh
u t = c 2 u xx
(1)
dengan x dan t adalah variabel jarak dan waktu, sedangkan u=u(x,t) menunjukkan densitas massa di titik x pada saat t. Konstanta c merupakan difusifitas yang bergantung pada sifat fisik material dimana massa merambat. Penurunan persamaan (1) dapat dilihat pada beberapa buku, sala satunya adalah Kreyszig (1994).
72
Metode Beda Hingga Untuk Analisis Proses Transfer Massa
Solusi persamaan difusi untuk proses transfer massa dengan sumber tunggal diberikan oleh
u ( x, t ) =
1 2 exp − x . 2π t 2t 1
(2)
Ungkapan ini dikenal sebagai suatu fungsi densitas (kerapatan massa) peluang normal dengan rata-rata nol dan varians t. Aplikasi persamaan difusi untuk gerak Brown pada pemodelan stokastik seperti yang dibahas pada Taylor & Karlin (1998), menggunakan solusi
u ( y, t | x) =
1 exp − ( y − x) 2 . 2π t 2t 1
(3)
Persamaan (3) menyatakan peluang posisi partikel berada di titik y setelah pengamatan berjalan selama t satuan waktu. Peluang ini yang terdistribusi normal dengan rata-rata posisinya berada di titik awal x, sedangkan variansnya berupa lamanya pengamatan t. Pada dasarnya pembahasan solusi berdenyut untuk persamaan difusi melibatkan singularitas di titik sumber (Cahyono dkk, 2004).
u ( x, y , t ) =
1
N
∑a π n =0
n
x2 + y2 1 exp − 4(t − t n ) 4(t − t n )
(4)
untuk t N < t ≤ t N +1 . Solusi pada persamaan (4) bukan merupakan pendekatan masalah bila sumber melepaskan massa secara kontinu. Untuk masalah dengan sumber melepaskan massa secara kontinu dapat diselesaikan dengan menggunakan metode numerik. Salah satu metode numerik untuk menyelesaikan masalah yang demikian adalah metode beda hingga. Metode beda hingga seperti disajikan pada Morton & Mayers (1996). Turunan pertama u terhadap t dan turunan kedua u terhadap x secara berturut-turut didekati sebagai berikut
u t ( x = x n , t = t p +1 ) =
U np +1 − U np ∆t
(5)
73
La Gubu et al//Paradigma, Vol. 15 No.2, Oktober 2011, hlm. 69-78
u xx ( x = x n , t = t p ) =
U np+1 − 2U np + U np−1
(∆x )2
.
(6)
Dengan demikian persamaan (1) dapat didekati dengan
U np +1 − U np U np+1 − 2U np + U np−1 = ∆t (∆x )2
(7)
yang secara eksplisit terhadap U np +1 dapat dituliskan dengan
[
U np +1 = U np + U np+1 − 2U np + U np−1
]
∆t . (∆x )2
(8)
3. SKEMA NUMERIK METODE BEDA HINGGA Skema numerik metode beda hingga untuk kasus 2-Dimensi (bidang x-y) dapat diilustrasikan dengan gambar berikut:
y ui,j+1 ui-1, j
ui,j
ui+1,
ui,j-1 x Gambar 2 Skema numerik metode beda ingga
Variabel u i , j pada skema di atas menyatakan densitas massa untuk titik x=i dan y=j pada saat tertentu t=p. Bentuk persamaan difusi untuk media dua dimensi yang homogen dapat dituliskan menjadi:
u t = u xx + u yy
(9)
dimana u menyatakan densitas massa yang merupakan fungsi dalam variabel ruang (x,y) dan waktu t, yaitu u=u(x,y,t). Misalkan u ( xi , y j , t p ) = u ip, j , maka untuk selang waktu ∆t akan diperoleh persamaan:
74
Metode Beda Hingga Untuk Analisis Proses Transfer Massa
u ip, +j 1 − u ip, j ∆t
= u t ( xi , y j , t p )
(10)
Turunan kedua u terhadap x dan y secara berturut-turut didekati sebagai berikut:
u xx ( xi , y j , t p ) = u yy ( xi , y j , t p ) =
u ip+1, j − 2u ip, j + u ip−1, j
(11)
(∆x )2
u ip, j +1 − 2u ip, j + u ip, j −1
(12)
(∆y )2
Jika media yang diamati dibuat dengan ukuran grid yang sama yakni ∆x = ∆y , maka pada selang waktu ∆t densitas massa u i , j dapat dituliskan menjadi:
u ip, j+1 = u iP, j +
∆t (∆x )2
{(u
p i +1, j
) (
− 2u ip, j + u ip−1, j + u ip, j +1 − 2u ip, j + u ip, j −1
)}
(13)
4. SIMULASI NUMERIK DENGAN SOFTWARE MATLAB Langkah awal dalam melakukan simulasi adalah membuat grid dari media atau luasan yang hendak diamati. Pada penelitian ini media yang diamati digrid seperti pada gambar berikut ini:
y
∆x
7
∆y
6 5 4 3 2 1
2
3
4
5
6
7
x Gambar 3. Grid media pengamatan
La Gubu et al//Paradigma, Vol. 15 No.2, Oktober 2011, hlm. 69-78
75
Media yang diobservasi digrid dengan ukuran grid 6x6 dan ∆x = ∆y . Parameter input dalam program Matlab meliputi: nx = 6 (jumlah grid dalam arah x); ny = 6 (jumlah grid dalam arah y); ∆x = ∆y =1/10; dt = 10-3; dan kecepatan difusi c=1. Koefisien laju difusi didefiniskan dengan alfa = c 2
dt . ∆x 2
Syarat batas pada x=1, x=7, y=1 dan y=7 ditentukan oleh densitas massa terdekatnya dan dikali dengan koefisien laju difusi (alfa). Selanjutnya diasumsikan bahwa pada titik (3,5) dan titik (5,3) terdapat sumber massa yang melepaskan massa secara kontinu dengan densitas u(3,5) = u(5,3) = 1 satuan. Berikut ini adalah beberapa contoh hasil simulasi numerik metode beda hingga dengan implementasi program Matlab.
Gambar 4. Transver massa untuk t=20, t=40, t=60, t=80 Gambar 4 di atas memperlihatkan proses terjadinya transfer massa yang dihasilkan oleh 2 buah sumber . Berdasarkan hasil simulasi numerik pada berbagai interval waktu pengamatan yakni untuk
t=20, t=40, t=60 sampai t=80 tampak terlihat jelas proses
penyebaran massa yang ditandai oleh meningkatnya nilai rapat massa u di sekitar titik sumber. Seiring dengan berjalannya waktu t, proses transfer massa juga menyebar ke seluruh daerah observasi. Gambar 5 berikut ini memperlihatkan proses perambatan massa
Metode Beda Hingga Untuk Analisis Proses Transfer Massa
76
untuk interval waktu dalam orde ratusan, berturut-turut untuk t=200, t=400, t=600, dan t=800. Pada Gambar 5 ini tampak terlihat jelas adanya proses transfer massa. Hal ini ditandai dengan bertambahnya nilai rapat massa u yang cukup signifikan dan massa telah menyebar dihampir seluruh wilayah observasi.
Gambar 5. Transver massa untuk t=200, t=400, t=600, t=800
Gambar 6. Transver massa untuk t=2000, t=4000, t=6000, t=8000 Hasil simulasi numerik untuk interval waktu dalam orde ribuan diperlihatkan oleh Gambar 6 yakni untuk t=2000, t=4000, t=6000, dan t=8000. Pada Gambar 6 ini tampak terlihat bahwa laju perubahan nilai rapat massa u lebih kecil bila dibandingkan dengan laju perubahan nilai rapat massa u untuk interval waktu dalam orde puluhan maupun ratusan.
La Gubu et al//Paradigma, Vol. 15 No.2, Oktober 2011, hlm. 69-78
77
Untuk interval waktu t=4000, t=6000, dan t=8000 telah memperlihatkan adanya nilai rapat massa u yang konstan terhadap waktu. Dalam hal ini nilai rapat massa u telah mencapai tingkat saturasi (kestabilan). Secara alamiah hal ini bisa terjadi karena adanya syarat batas yang diberlakukan serta jumlah grid yang terbatas. Jika jumlah gridnya diperbesar maka waktu yang diperlukan dalam proses transfer massa untuk semua daerah observasi juga akan semakin lama. Selain itu cepat lambatnya proses transfer massa juga dipengaruhi oleh besarnya besar kecilnya laju difusi, dimensi grid ( ∆x , ∆y ) serta ∆t . Pada Gambar 7 berikut ini diperlihatkan nilai rapat massa u sebagai fungsi waktu (t) di empat titik berbeda yaitu u(4,4), u(3,2), u(5,6) dan u(1,1). Pada waktu
t >0 rapat
massa di titik u(4,4) mempunyai nilai yang lebih besar dibandingkan dengan tiga titik lainnya. Hal ini disebabkan karena titik u(4,4) berdekatan dengan sumber massa yakni di titik (3,5) dan titik (5,3). Pada dasarnya grafik rapat massa sebagai fungsi waktu (t) untuk ke empat titik di atas mempunyai pola yang sama. Terlihat jelas bahwa untuk t. Mendekati seribu maka nilai rapat massa mulai stabil (tidak berubah terhadap waktu).
Gambar 7. Rapat massa u sebagai fungsi dari waktu t 5. SIMPULAN Berdasarkan uraian dan hasil simulasi numerik yang diperoleh, dapat ditarik beberapa kesimpulan sebagai berikut : 1. Metode beda hingga melalui implementasi skema numerik dapat diterapkan untuk
Metode Beda Hingga Untuk Analisis Proses Transfer Massa
78
menyelesaikan masalah transfer massa dengan sumber melepaskan massa secara kontinu. 2. Hasil simulasi numerik transfer massa memperlihatkan bahwa untuk daerah observasi −3 dengan ukuran grid 6x6, dimensi grid ∆x = ∆y = 0.1 dan parameter input ∆t = 10
membutuhkan waktu dalam orde ribuan untuk mencapai semua daerah observasi. 3. Semakin dekat titik pengamatan dengan sumber maka waktu yang dibutuhkan oleh transfer massa untuk mencapai titik tersebut akan semakin cepat dan begitu pula sebaliknya. Nilai rapat massa di sekitar sumber selalu lebih besar dari pada yang jauh dari sumber. 4. Dari grafik rapat massa sebagai fungsi waktu diperoleh bahwa rapat massa pada suatu titik akan bertambah seiring dengan bertambahnya waktu. Tetapi pada interval waktu dengan orde ribuan rapat massa u di titik tersebut mulai stabil (tersaturasi). Secara alamiah hal ini bisa terjadi karena adanya syarat batas yang diberlakukan serta kriteria parameter input yang diberikan.
DAFTAR PUSTAKA [1] Banga, J. R., Perez-Martin, R. I., Gallardo, J. M. & Casares, J. J. (1991) Optimization of the thermal processing of conduction-heated canned foods: study of several objective functions. Journal of Food Engineering 14, pp. 25-51. [2] Belolipetski, A.A & Krikorov,A.M.(1984) Fundamental Solutions of the Non Linear Equations of Heat Condution,USSR Computational Mathematics and Mathematical Physics. 22(3),141-149. [3] Cahyono, E & La Gubu,(2004).Wood drying process: modeling for linear diffusivity with respect to the water content”, Proceeding SEMINAR MIPA IV ITB, 112114. [4] Edi Cahyono & Yudi Soeharyadi, “On the interaction of Barenblatt solutions to the porous medium equation”, International Journal of Mathematical Sciences and Engineering Applications, Vol. 4, No. II, 2010. [5] Kreyszig, E. (1994) Advance Engineering Mathematics, John Wiley & Son, Singapore. [6] Morton, K. W. & Mayers, D. F. (1996) Numerical Solution of Partial Differential Equations, Cambridge University Press, Cambridge. [7] Wibowo, A. & Ponidi. (2002) Dinamika panas pada pusat kaleng dalam proses sterilisasi makan kaleng yang berbentuk tabung. MATEMATIKA: Jurnal Matematika dan Pembelajarannya, ISSN: 0852-7792, Th VIII, Edisi Khusus.
La Gubu et al//Paradigma, Vol. 15 No.2, Oktober 2011, hlm. 69-78
79