JMA, VOL. 2, NO.1, JULI, 2003, 37-44
37
SOLUSI NUMERIK PERSAMAAN BOLTZMANN ENDAR H. NUGRAHANI Departemen Matematika, Fakultas Matematika dan Imu Pengetahuan Alam, Institut Pertanian Bogor Jln. Meranti, Kampus IPB Dramaga, Bogor 16680, Indonesia
Abstract : Particle methods are well known tools to solve the kinetic Boltzmann equation numerically. The usual procedure of such method is the direct simulation Monte-Carlo, which directly simulate the collision between particles. The second method of interest is the stochastic weighted particle method, which is developed to improve the previous method. The main idea of the second method is to use random weight transfer between particles during collisions. In order to reduce the stochastic fluctuations, this method provides a way to increase the number of particles. But if the additional particles cannot be compensated in some natural way, then this number should be reduced. Several reduction procedures have been proposed. Some numerical results using both methods are presented. It is shown that the second method gives some better results in some ways. Katakunci: persamaan Boltzmann, metode Monte-Carlo, metode partikel terbobot stokastik.
1. PENDAHULUAN Persamaan Boltzmann adalah persamaan kinetik yang mendeskripsikan perubahan yang terjadi pada partikel gas dengan berlalunya waktu karena terjadinya fenomena tumbukan antar partikel. Misalkan f ( x, v, t ) : 3 3
adalah fungsi kepekatan peluang suatu partikel gas yang berada pada posisi x, bergerak dengan kecepatan v, serta diamati pada waktu t. Bentuk dari persamaan Boltzmann adalah: f v, x f Q( f , f ) t
(1.1)
f ( x, v, t 0) f 0 ( x, v) .
(1.2)
dengan nilai awal
38
ENDAR H. NUGRAHANI
Ruas kanan dari persamaan (1.1) adalah integral lipat berdimensi tinggi yang mendeskripsikan perubahan dari f(x,v,t) sebagai akibat dari tumbukan antar partikel, yaitu dengan bentuk berikut: Q( f , f )
B(v, w, e) f (v' ) f (w' ) f (v) f (w) dw de
3 S 2
dengan B(v,w,e) adalah fungsi kernel tumbukan partikel. Solusi analitik dari persamaan Boltzmann (1.1) tersebut tidaklah mudah untuk diperoleh. Beberapa hasil solusi analitik eksak yang telah didapatkan hanya bisa diperoleh untuk distribusi nilai awal tertentu, di antaranya adalah untuk distribusi awal Maxwell oleh Bobylev(2) dan beberapa solusi self-similar oleh Bobylev dan Cercignani (3). Solusi numerik dari persamaan Boltzmann harus memperhatikan diskretisasi dari ruang berdimensi tinggi serta melibatkan perhitungan integral lipat berdimensi tinggi pula pada setiap titik diskretisasinya. Dengan demikian penggu-naan metode yang umum digunakan pada solusi numerik persamaan diferensial parsial, misalnya metode elemen-hingga ataupun metode bedahingga, akan sulit diimplementasikan. Oleh karena itu penggunaan metode partikel akan menawarkan jalan keluar bagi penyelesaian numerik persamaan kinetik yang berupa persamaan diferensial–integral berdimensi tinggi tersebut. 2. METODE PARTIKEL Metode partikel untuk penyelesaian numerik persamaan kinetik berdimensi tinggi pertama kali dikembangkan adalah simulasi langsung dengan metode Monte-Carlo (Direct Simulation Monte-Carlo – DSMC) oleh G. A. Bird (1). Pada perkembangannya dirumuskan metode partikel terbobot stokastik (Stochastic Weighted Particle Method – SWPM) oleh S. Rjasanow dan W. Wagner (6). Ide dasar metode partikel adalah pendugaan barisan ukuran berikut f t k , x, v dxdv, t k kt , k 0,1, , t 0
dengan sistem N
t k , dx, dv g i t k xi tk ,vi tk dx, dv i 1
di mana
g i t k , xi t k , vi t k iN1 ,
k 0,1,
melambangkan sekumpulan partikel, dengan k adalah indeks tahap waktu. 2.1 Metode Simulasi Langsung Monte-Carlo (DSMC). Metode DSMC yang dikembangkan oleh G. A. Bird(1) memiliki algoritma berikut. Prosedur simulasi diawali dengan pendefinisian evolusi waktu dari suatu sistem sejumlah N partikel
g i (t ), xi (t ), vi (t ), dengan nilai awal
i 1,, N , t t 0
JMA, VOL. 2, NO.1, JULI, 2003, 37-44
g i (t 0 ), xi (t 0 ), vi (t 0 ),
39
i 1, , N
yang dibangkitkan menurut suatu sebaran peluang tertentu berdasarkan (1.2). Asumsi pada nilai awal adalah bobot partikel konstan g i (t 0 ) g , i 1,, N .
(2.1)
Diskretisasi waktu didefinisikan sebagai berikut: t k kt , k 0,1, t k 1 t k t , t 0.
Simulasi akan dilakukan pada setiap interval waktu (tk, tk+1] dengan membedakannya dalam dua fase: fase gerakan bebas dan fase tumbukan partikel. Pada fase pertama hanya akan disimulasikan gerakan partikel secara bebas tanpa memperhatikan terjadinya tumbukan antar partikel. Sebaliknya pada fase kedua hanya disimulasikan tumbukan antar partikel tanpa memperhatikan gerakan partikel itu sendiri. Fase gerakan bebas menghasilkan perubahan posisi partikel setiap waktu hanya dipengaruhi oleh kecepatannya, yaitu xi (t ) xi (t k ) (t t k )vi (t k ), t t k , t k 1 .
Bobot partikel diasumsikan tidak mengalami perubahan selama periode simulasi. Apabila daerah definisi adalah berbatas, maka kecepatan partikel akan berubah menurut syarat batasnya. Fase tumbukan partikel disimulasikan sebagai berikut. Pertama-tama posisi partikel dipartisi ke dalam c buah sel
c
. 1
Simulasi tumbukan dilakukan mendefinisikan proses stokastik
di
dalam
masing-masing
sel
dengan
Z (t ) g i (t ), xi (t ), vi (t ), i 1,, n
dengan t tk dan n = nl(tk). Dari proses pada keadaan z, misalkan partikel ke – i dan j akan bertumbukan dengan peluang tertentu. Pembangkit proses ditentukan berdasarkan fungsi loncatan proses Jk(z,i,j,e) berikut g , x k , v k , k i, j J k (.) g , xi , v'i , k i g,x , v' , k j j j
(2.2)
di mana v'i vi e e, vi v j v' j v j e e, vi v j
(2.3)
adalah kecepatan partikel setelah tumbukan. Dengan algoritma DSMC tersebut jumlah partikel di dalam sistem tidak akan mengalami perubahan. Dengan demikian juga akan berlaku konservasi berbagai besaran makroskopik dari sistem, yaitu diantaranya besaran massa, kerapatan, momentum, energi, dan sebagainya. Sedangkan kekonvergenen
40
ENDAR H. NUGRAHANI
pendugaan numerik dengan metode DSMC ini telah dibuktikan dalam Wagner (8) . 2.2 Metode Partikel Terbobot Stokastik (SWPM). Simulasi numerik dengan menggunakan metode partikel terbobot stokastik memiliki latar belakang ide yang serupa dengan DSMC, dengan letak perbedaan pada pendefinisian fungsi lompatan proses sebagai alternatif dari (2.2) berikut: g k , x k , v k , k i, j g G (.), x , v , k i i i i J k . g j G (.), x j , v j , k j G (.), x , v' , k n 1 i i G (.), x j , v' j , k n 2
(2.4)
di mana v’i dan v’j adalah kecepatan partikel ke – i dan j setelah tumbukan partikel. Besaran G(z,j,i,e) adalah fungsi perpindahan bobot, G z , i, j , e
min g i , g j 1 ( z , i , j , e)
(2.5)
yang akan memungkinkan terjadinya penambah-an jumlah partikel dalam sistem bila dipilih parameter 0. Prosedur yang memungkinkan penambahan jumlah partikel ini terutama akan sangat menguntungkan untuk simulasi sistem dengan kerapatan rendah untuk mengurangi keragaman stokastik(6). Akan tetapi penambahan jumlah partikel dalam simulasi numerik dapat meningkatkan beban komputasi secara signifikan, sehingga sampai tahap tertentu akan diperlukan prosedur tambahan untuk mengurangi kembali jumlah partikel tersebut. Prosedur pengurangan jumlah partikel yang telah dirumuskan antara lain diberikan dalam Rjasanow dan Wagner (6) sebagai prosedur reduksi deterministik, sedangkan dalam Nugrahani dan Rjasanow (5) dirumuskan prosedur reduksi stokastik. Prosedur Reduksi Deterministik Prosedur untuk mengurangi jumlah partikel pada simulasi dengan metode SWPM ini dilakukan sedemikian sehingga sistem partikel tersebut tetap dapat menduga penyelesaian numerik persamaan Boltzmann. Dengan demikian tetap harus diperhatikan terpenuhinya sifat konservasi masa, kerapatan, momentum dan energi, serta dapat dikontrolnya tambahan galat pendugaan yang ditimbulkan. Ada dua langkah yang dirumuskan dalam prosedur reduksi deterministik ini. Langkah pertama adalah pengelompokan partikel dalam sejumlah cluster berdasarkan kecepatannya:
g
ij
, xij , vij , i 1,, nˆ, j 1,, ni
di mana nˆ adalah jumlah cluster yang dibentuk, dengan i adalah indeks cluster dan j indeks partikel dalam cluster. Kemudian dalam masing-masing cluster dengan ni 3 akan dilakukan pengurangan jumlah partikel, yaitu dengan menggantikan sejumlah 2 partikel berikut:
JMA, VOL. 2, NO.1, JULI, 2003, 37-44
vi1 V (i )
41
g (i ) g i1 g i 2 2 (i ) e, vi 2 V (i ) (i ) e, e S 2
di mana g(i) adalah massa cluster, V(i) adalah rataan impuls cluster, serta (i) adalah suatu besaran yang berhubungan dengan aliran energi dalam cluster. Bukti kekonvergenan serta jaminan kontrol atas galat tambahan dari prosedur reduksi ini diberikan oleh Rjasanow dan Wagner (6). Prosedur Reduksi Stokastik Prosedur ini ditujukan secara khusus untuk mempertahankan partikel pada daerah ’luar’
v 3 ; || v U || R, R . Sedangkan partikel pada daerah 3 akan dipilih menurut kriteria stokastik tertentu. Dengan memperhatikan partikel ke – i:
g i , xi , vi ,
i 1,, n
prosedur reduksi dilakukan dengan cara berikut:
Dengan peluang p
gi , tetapkan g i g . Hal ini berarti partikel tersebut g
dipertahankan dan diberi bobot standar.
Dengan peluang 1 p
g gi , tetapkan g
gi = 0. Hal ini berarti partikel tersebut dikeluarkan dari sistem, sehingga jumlah partikel akan berkurang. Prosedur reduksi ini diusulkan dalam Nugrahani dan Rjasanow (5) sebagai bentuk khusus dari prosedur yang lebih umum yang dirumuskan dan dibuktikan kekonvergenannya oleh Matheis dan Wagner (4). 3. HASIL SIMULASI NUMERIK 3.1 Simulasi Extra Momen. Studi numerik yang dilakukan adalah penyelesaian persamaan Boltzmann (1.1) pada ruang posisi yang telah dipartisi sehingga fungsi kepekatan peluang f(v,t) adalah homogen dalam ruang. Kondisi awal mengikuti distribusi molekul pseudo-Maxwell campuran berikut f (v )
| v V1 | 2 | v V2 | 2 (1 ) exp exp 2T1 2T2 2T1 3 / 2 2T2 3 / 2
Simulasi numerik dilakukan terhadap besaran extra momen berikut (t )
f (v, t )dv
3
dengan mempergunakan metode DSMC dan SWPM pada prosedur reduksi stokastik. Hasil simulasi tersebut disajikan pada Gambar 1 untuk hasil simulasi metode DSMC, serta Gambar 2 untuk metode SWPM.
42
ENDAR H. NUGRAHANI
Gambar 1. Simulasi extra momen dengan DSMC
Gambar 2. Simulasi extra momen dengan SWPM
Dari kedua grafik tersebut dapat diamati bahwa metode DSMC memberikan hasil pendugaan yang tidak cukup memuaskan, terutama setelah suatu tahap waktu tertentu metode ini tidak dapat lagi memberikan interval kepercayaan pendugaan yang positif. Sedangkan metode SWPM memberikan hasil pendugaan yang lebih baik dalam hal masih dapat diberikannya interval kepercayaan yang positif (5). Lebih lanjut disajikan pula berbagai besaran hasil simulasi extra momen menggunakan metode SWPM dengan prosedur reduksi stokastik. Gambar 3 menyajikan jumlah partikel dalam sistem yang memperlihatkan bahwa prosedur reduksi stokastik memperkecil fluktuasi jumlah partikel. Gambar 4 memberikan bobot maximal partikel dalam keseluruhan proses simulasi, yang ternyata benilai relatif cukup kecil. Serta Gambar 5 memberikan dugaan nilai kepadatan partikel dalam sistem.
Gambar 3. Jumlah partikel pada simulasi extra momen dengan SWPM
JMA, VOL. 2, NO.1, JULI, 2003, 37-44
43
Gambar 4. Bobot maximal partikel pada simulasi extra momen dengan SWPM
Gambar 5. Kepadatan partikel pada simulasi extra momen
3.2 Simulasi Self-similar dengan SWPM. Bagian berikutnya dari studi simulasi numerik diberikan untuk memberikan penyelesaian bagi besaran α – momen berikut [|| v || ](t ) || v || f (v, t )dv 3
dengan f(v,t) adalah penyelesaian self-similar menurut Bobylev dan Cercignani(3) dengan distribusi awal berikut: f (v )
4
ds
1 s
2 2
0
1 M || v ||, 2 s
dengan M(.) adalah distribusi Maxwell. Keunikan dari model self-similar adalah bentuk distribusinya yang divergen, sehingga secara umum konservasi energi tidak dapat dipertahankan(9). Sehingga solusi numerik α –momen tersebut hanya dapat dilakukan pada nilai α tertentu. Pada contoh ini dipilih α = 0,75 dan metode yang dipergunakan adalah SWPM dengan prosedur reduksi deterministik. Hasil simulasi disajikan pada Gambar 6, yang memberi gambaran bahwa metode SWPM terbukti cukup mampu untuk memberikan hasil numerik bagi pendugaan fungsi distribusi yang bukan maxwellian. Akan tetapi dapat dilihat bahwa kedekatan hasil simulasi numerik ini akan sangat tergantung pada jumlah partikel awal yang dipergunakan dalam simulasi, semakin kecil jumlah partikel maka pendugaan akan semakin jauh dari solusi analitik eksak.
44
ENDAR H. NUGRAHANI
Gambar 6. Simulasi self-similar pada α – momen dengan metode SWPM.
4. KESIMPULAN Penyelesaian numerik persamaan Boltzmann yang merupakan persamaan kinetik berdimensi tinggi dapat diperoleh dengan menggunakan metode partikel. Metode partikel dalam bentuk simulasi langsung Monte-Carlo (DSMC) adalah metode yang umum diterapkan, akan tetapi ternyata tidak cukup mampu memberikan hasil yang memuaskan pada daerah dengan kepadatan rendah. Di lain pihak simulasi dengan metode partikel terbobot stokastik (SWPM) terbukti mampu memberikan hasil yang lebih memuaskan. Metode ini dapat memberikan pendugaan yang baik di daerah dengan kepadatan partikel rendah. Lebih lanjut metode ini juga terbukti mampu memberikan dugaan numerik bagi fungsi distribusi awal yang divergen, asalkan banyaknya partikel yang disimulasikan adalah cukup besar. DAFTAR PUSTAKA 1. 2. 3. 4. 5.
6. 7.
8. 9.
G. A. Bird. Molecular Gas Dynamics. Clarendon Press, Oxford, 1976. A. V. Bobylev. Exact solutions of the Boltzmann equation. Sov. Phys. Dokl., 20(12): 822 – 824, 1975. A. V. Bobylev and C. Cercignani. Self-similar solutions of the Boltzmann equation and their applications. J. Statist. Phys., 106(5-6): 1039 – 1071, 2002. Matheis and W. Wagner. Convergence of the stochastic weighted particle method for the Boltzmann equation. Preprint 739, WIAS, Berlin, 2002. E. H. Nugrahani and S. Rjasanow. On the stochastic weighted particle method. In M. Griebel and M. A. Schweitzer, editors, Meshfree methods for partial differential equations, Volume 26 of Lecture notes in computational science and engineering, pages 319 – 326. Springer, Berlin, 2003. S. Rjasanow and W. Wagner. A Stochastic Weighted Particle Method for the Boltzmann Equation. J. Comput. Phys., 124: 243 – 253, 1996. S. Rjasanow, T. Schreiber and W. Wagner. Reduction of the number of particles in the stochastic weighted particle method for the Boltzmann equation. J. Comput. Phys., 145(1): 382 – 405, 1998. Wagner, W. A convergence proof for Bird’s direct simulation Monte-Carlo method for the Boltzmann equation. J. Stat. Phys., 66(3-4): 1011-1044, 1992. E. H. Nugrahani. Beiträge zur Numerik der Boltzmann Gleichung. Dissertation. University of Saarland, Germany, 2003.
