Research and Development on Nanotechnology in Indonesia, Vol.2, No.2, 2015, pp. 68-79 ISSN : 2356-3303
SIMULASI DEPOSISI BERDASARKAN PENUMBUHAN BUTIRAN MENGGUNAKAN DINAMIKA MOLEKULAR Sparisoma Viridi1,a, Suprijadi2,b, Veinardi Suendo3,c 1
Fisika Nuklir dan Biofisika, FMIPA, Institut Teknologi Bandung, Bandung 40132, Indonesia 2 Fisika Teoretik Energi Tinggi dan Instrumentasi, FMIPA, Institut Teknologi Bandung, Bandung 40132, Indonesia 3 Kimia Fisik, FMIPA, Institut Teknologi Bandung, Bandung 40132, Indonesia * e-mail:
[email protected],
[email protected], 3
[email protected] Received : 16 January 2015 Accepted : 18 February 2015
ABSTRAK Penumbuhan butiran berbentuk cakram dua-dimensi dapat digunakan untuk melakukan simulasi deposisi bahan. Metode dinamika molekuler digunakan dengan menerapkan gaya tolak berupa gaya normal berdasarkan model pegas linier-dashpot, gaya gravitasi, dan gaya ikat buatan yang bergantung dari besar overlap sebelum ikatan terjadi. Laju jatuhnya butiran menentukan bagaimana struktur yang terbentuk pada lapisan tempat deposisi terjadi. Variasi ukuran butiran dilakukan untuk menunjukkan efek dari konsentrasi. Kata kunci: deposisi, dinamika molekuler, simulasi.
PENDAHULUAN Peristiwa deposisi merupakan suatu hal menarik untuk dikaji, yang dalam simulasi telah dilakukan menggunakan grid pada kasus etching dan deposisi [Hsiau et al., 1997], mengunakan random walk pada kasus deposisi partikel akibat aliran turbulen [Nagy et al., 1996], dan menggunakan dinamika molekuler pada kasus deposisi partikel butiran pada butiran lebih besar [Sustini et al., 2011]. Dalam tulisan ini akan dibahas deposisi pada suatu lapisan 2-d yang juga dibentuk dari butiran (butiran alas). Peran dari gaya luar selain gaya gravitasi bumi akan dibahas untuk melihat pengaruhnya
68| CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
terhadap situs-situs deposisi yang terbentuk, yang diharapkan dapat memperoleh bentuk-bentuk seperti dilaporkan dalam beberapa eksperimen [Wu et al., 2002; Yu et al., 2013; Sagu et al. 2014], setidaknya untuk versi potongan dalam 2-d. SIMULASI Deposisi partikel dalam bentuk 2-d dimodelkan dengan menggunakan partikel butiran i bermassa m dan berdiameter D yang dijatuhkan dari posisi tertentu (1) ri (0) = H + LU [0,1] di mana L adalah lebar sistem 2-d, H adalah ketinggian tetap, dan U[0,1] adalah distribusi uniform yang membangkitkan bilangan acak antara 0 dan 1 dengan fungsi kerapatan probabilititasnya adalah
0 ≤ x ≤ 1, 1, f (x ) = 0, x < 0 ∨ x > 1.
(2)
Kecepatan awal partikel i selalu dipilih sama dengan nol saat mulai dijatuhkan dari ri (0 ) . Ilustrasi sistem yang dimaksud diberikan dalam Gambar 1. Selain percepatan gravitasi g , jatuhnya partikel juga dipengaruhi oleh gaya luar F (t ) yang dapat berbentuk sembarang dan juga dapat merupakan fungsi dari waktu t. Dinamika partikel diperoleh dengan menggunakan hukum kedua Newton tentang gerak yang memberikan 1 ai (t ) = F (t ) + g . (3) m Selanjutnya, dengan menggunakan metode improved Euler yang sudah sering digunakan dalam kasus terkait interaksi butiran [Haris et al., 2014], dapat diperoleh kecepatan dan posisi partikel butiran saat t + Δt melalui (4) vi (t + ∆t ) = vi (t ) + a (t )∆t , (5) ri (t + ∆t ) = ri (t ) + v (t )∆t . Persamaan (3)-(5), yang umumnya dilakukan dengan metoda integrasi dengan ketelitian lebih tinggi, dikenal dengan metode dinamika molekuler. Dalam simulasi ini interaksi antar butiran yang umumnya menggunakan model pegas linier-dashpot tidak digunakan, akan tetapi definisi dari overlapnya tetap digunakan untuk syarat terminasi gerak butiran yang dijatuhkan bebas dalam pengaruh g dan F (t ) , sebagaimana telah digunakan pada kasus deposisi beberapa butiran kecil pada sebuah butiran besar [Sustini et al., 2011].
69 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
ri (0) vi (0) = 0 m, D g
F(t) l h B L
Gambar 1. Model deposisi 2-d berbasiskan granular dengan setiap butiran i bermassa m dan berdiameter D yang akan tumbuh pada butiran alas B. Overlap antara antara butir i dan j yang masing-masing terletak pada posisi ri dan r j dan berdiameter D i dan D j adalah [Schäfer et al., 1996]
ξ ij = max(0, Ri + R j − rij ) ,
(6)
di mana
1 Di , 2 rij = rij = ri − r j . Ri =
(7) (8)
Fungsi max(x, y) dalam Persamaan (6) didefinisikan sebagai x, x > y , max(x, y ) = . (9) y, x < y. Saat suatu butiran pertama kali dijatuhkan nilai overlapnya dengan butiran alas B adalah nol. Kemudian, saat butiran tersebut mencapai butiran alas, nilai overlap yang diberikan oleh Persamaan (6) menjadi bernilai, yang mengindikasikan bahwa kecepatan butiran tersebut secara mendadak dibuat nol dan posisinya ditetapkan di sana. Dengan demikian, bila suatu saat suatu butiran baru dijatuhkan telah terdapat N B butiran alas dan telah terdapat N D butiran, maka interaksi butiran yang baru ini perlu diperiksa terhadap N B + N D butiran sebelumnya yang telah ada. Sebagai konsekuensi, waktu
70 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
perhitungan akan semakin lama dengan bertambahnya butiran yang telah terdeposisi pada butiran basis. Syarat batas periodik juga diterapkan pada butiran yang jatuh bila gaya luar F (t ) yang diterapkan membuat butiran melewati rentang posisi x > L atau x < 0 dengan cara x (10) x = L, L di mana notasi menyatakan fungsi floor(x) atau
(11) x = floor(x ) = max{m ∈ Ζ | m ≤ x}, di mana x adalah bilangan riil, m bilangan bulat, dan Ζ adalah himpunan bilangan bulat (positif, negatif, dan nol).
KARAKTERISASI SITUS DEPOSISI Suatu kotak berukuran l × h digunakan untuk menghitung kerapatan situs deposisi yang tumbuh yang dinyatakan dengan ϕ N πD 2 φ = lh , (12) 4lh dimana N lh adalah jumlah partikel butiran dalam luas l × h tersebut. Selain itu jarak antar masing-masing situs juga akan dibahas terkait dengan gaya luar yang diberikan. Bila terdapat beberapa nilai maka sebaran nilai atau distribusi dari λ akan didisusikan. Ilustrasi mengenai jarak antar situs deposisi diperlihatkan dalam Gambar 2.
λ1
λ2
B Gambar 2. Tiga situs deposisi yang terbentuk dengan jarak pisahnya masing-masing adalah λ 1 dan λ 2 . HASIL DAN DISKUSI Parameter simulasi yang digunakan adalah L = 1, H = 1, g = 1, Δt = 10-3, D = 10-2, m = 10-3, dan N B = 100, sedangkan parameter lain diubah dengan nilai dijelaskan pada gambar atau teks terkait.
71 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
(a)
(b)
(c)
(d)
Gambar 3. Hasil deposisi dengan N D = 100 dan jumlah situs N S : (a) 10, (b) 20, (c) 50, dan (d) 100. Hasil pertama yang dibahas adalah situs tempat terjadinya deposisi ditentukan jumlahnya N S dan jarak antara situs dibuat tetap. Apabila jumlah partikel yang terdeposisi dibuat tetap N D = 100 maka dengan semakin besar nilai N S akan semakin rendah situs deposisi yang tumbuh seperti diberikan dalam Gambar 3. Tinggi situs diberikan oleh paremeter l. Dengan menggunakan N S = 20 dan N D = 500 ditambah dengan adanya gaya luar F ext yang berarah ke kiri akan mempengaruhi arah pertumbuhan butiran pada situs-situs deposisi sebagaimana digambarkan dalam Gambar 4 berikut. Jumlah situs deposisi yang mungkin apabila dikaitkan dengan ukuran partikel dan lebar sistem dapat menentukan apakah pertumbuhan situs akan menghasilkan rongga-rongga ataupun tidak. Gambar 5 sekilas memperlihatkan bahwa rongga akan muncul saat N S ≈ N B atau lebih besar. Pertumbuhan seperti bunga karang mulai muncul saat N S > N B .
(a)
72 | CAS – Center for Advanced Sciences
(b)
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
(b) (d) Gambar 4. Hasil deposisi dengan nilai parameter N D = 500, jumlah situs N S = 20, dan gaya luar F ext : (a) -2×10-4, (b) -1×10-4, (c) -5×10-5, dan (d) -1×105. Kemiringan situs deposisi yang tumbuh bergantung dari nilai gaya luar yang diberikan. Semakin besar gaya luar, semakin miring situs butiran yang tumbuh.
(a)
(b)
(c)
(d)
Gambar 5. Hasil deposisi dengan nilai parameter N D = 1000, F ext = 0, jumlah situs N S : (a) 100, (b) 200, (c) 300, dan (d) 800. Gaya luar yang diberikan dapat dimodifikasi sebagai fungsi spasial seperti berikut ini + Fext , x < 12 L, . (13) Fext = − Fext , x > 12 L. Implementasi Persamaan (13) diberikan dalam Gambar 6, yang untuk sementara belum memberikan pengaruh yang berarti. Penerapan Persamaan (13) dengan melakukan variasi nilai N S dapat memberikan bentuk-bentuk lain yang menarik sebagai mana diberikan dalam Gambar 7. Untuk itu digunakan |Fext| = 20×10-5 dengan berbagai nilai N S . Teramati bahwa struktu mirip kristal yang cukup rapih terbentuk saat N S < N B seperti dalam Gambar 7(a) dan 7(b).
73 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
(a)
(b)
(c) (d) Gambar 6. Hasil deposisi dengan nilai parameter N D = 1000, N S = 200, F ext > 0 untuk x < L/2, F ext < 0 untuk x > L/2, dengan |Fext|: (a) 4×10-5, (b) 6×10-5, (c) 8×10-5, dan (d) 10×10-5.
(a)
(b)
(c)
(d)
74 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
Gambar 7. Hasil deposisi dengan nilai parameter N D = 1000, F ext > 0 untuk x < L/2, F ext < 0 untuk x > L/2, dengan |Fext| = 20×10-5 dengan N S : (a) 25, (b) 50, (c) 100, dan (d) 200.
Terdapat pula bentuk gaya luar lain yang diberikan, yaitu 1 Fext = cL − x, . (14) 2 di mana c akan divariasikan. Parameter lainnya bernilai N D = 1000. Variabel x adalah komponen dari vektor posisi butiran yang terdeposisi.
(a)
(b)
(c)
(d)
Gambar 8. Hasil deposisi dengan nilai parameter N D = 1000, menggunakan Persamaan (14), c = 0.001, dan dengan N S : (a) 12, (b) 25, (c) 50, dan (d) 100.
75 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
Bentuk seperti dalam Gambar 8 terlihat bahwa proporsi bagian yang mirip kristal berkurang dengan semakin besarnya nilai N S /N B . 92 90 88 θ (°)
86 84 82 80 78 -25
-20
-15
-10
-5
0
F ext (10-5 )
Gambar 9. Sudut yang dibentuk oleh tumpukan butiran berbentuk garis seperti dalam Gambar 4 sebagai fungsi dari F ext . Sudut yang dibentuk oleh tumpukan butiran dalam Gambar 4 memiliki kebergantungan hampir linier dengan besar gaya luar F ext yang diberikan. Hubungan keduanya diberikan dalam Gambar 9. 0.18 0.16 h/l 0.14 0.12 0.1 0
200
400
600
800
1000
NS
Gambar 10. Nilai h/l untuk berbagai nilai N S sebagaimana hasil deposisi diberikan dalam Gambar 5. Tinggi struktu di tengah dalam Gambar 6 memiliki kebergantungan dari F ext yang diberikan, hubungan antara l/h dengan gaya luar yang diberikan disajikan dalam Gambar 11 berikut ini.
76 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
0.24 0.22 0.2 0.18 l/h 0.16 0.14 0.12 0.1 0.08 4
6
8
10
F ext (10-5 )
Gambar 11. Nilai l/h untuk berbagai nilai F ext sebagaimana hasil deposisi diberikan dalam Gambar 6 hanya untuk struktur bagian tengah. Sedangkan untuk Gambar 7 apabila dihitung kotak l × h seperti dalam Gambar 1, maka akan diperoleh hasil seperti dalam Gambar 12 berikut ini. 0.81 0.8 0.79 0.78 l/h 0.77 0.76 0.75 0.74 0
50
100
150
200
250
NS
Gambar 12. Nilai l/h untuk berbagai nilai N s sebagaimana hasil deposisi diberikan dalam Gambar 7.
77 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
80
1.25
70
1.2
60 1.15 l/h
50
1.1
k 40 30
1.05
20 1
10
0.95
0 0
20
40
60
80
NS
100
0
20
40
60
80
100
NS
Gambar 13. Nilai l/h (kiri) dan keteraturan k (kanan) untuk berbagai nilai N s dari Gambar 8. Struktur dalam Gambar 8 memberikan nilai l/h dan k seperti dalam Gambar 13, di mana nilai terendah dalam l/h bertepatan dengan nilai tertinggi dalam k. KESIMPULAN Variasi gaya luar F ext dan jarak situs deposisi N S mempengaruhi bagaimana butiran di tumbuhkan. Terdapat bentuk menyerupai sedikit kristal, struktur periodik, dan struktur dengan banyak rongga. Perbedaan bentuk ini dikendalikan dengan nilai F ext yang berkompetisi dengan nilai N S relatif terhadap N B . UCAPAN TERIMA KASIH Sosialisasi karya ini didukung oleh National Research Center for Nanotechnology dan bagian perhitungannya didukung oleh RIK ITB dengan Nomor 914/AL-J/DIPA/PN/SPK/2014. REFERENSI 1] T. K. Nagy, D. Swailes, “Random Walk Approach for Simulation of Particle Deposition from Turbulent Flows”, Periodica Polytechnica Mechanical Engineering 40 (2), 143-156 (1996). 2] L. Haris, S. N. Khotimah, F. Haryanto, S. Viridi, "Molecular Dynamics Simulation of Soft Grains: Malaria-Infected Red Blood Cells Motion within Obstructed 2-D Capillary Vessel", in Symposium on Biomathematics-2013, edited by H. Arimura et al., AIP Conference Proceedings 1587, American Institute of Physics, Melville, NY, 2014, pp. 43-46. 3] Z.-K. Hsiau, E. C. Kan, J. P. McVittie, R. W. Dutton, “Robust, Stable, and Accurate Boundary Movement for Physical Etching and
78 | CAS – Center for Advanced Sciences
Viridi et al., RDNI, Vol. 2, No.2, 2015, pp. 68-79
4]
5] 6]
7]
8]
Deposition Simulation”, IEEEE Transaction on Electron Devices 44 (9), 1375-1385 (1997). J. S. Sagu, T. A. N. Peiris, K. G. U. Wijayantha, “Rapid and Simple Potentiostatic Deposition of Copper (II) Oxide Thin Films”, Electrochemistry Communications 42 (), 68-71 (2014). J. Schäfer, S. Dippel, and D. E. Wolf, “Force Schemes in Simulations of Granular Materials”, Journal de Physique I 6 (1), 5-20 (1996). E. Sustini, S. N. Khotimah, F. Iskandar, S. Viridi, "Molecular Dynamics Simulation of Smaller Granular Particles Deposition on a Larger One Due to Velocity Sequence Dependent Electrical Charge Distribution", in The 4th Nanoscience and Nanotechnology Symposium-2011, edited by F. Iskandar et al., AIP Conference Proceedings 1415, American Institute of Physics, Melville, NY, 2011, pp. 209-213. Y. Wu, H. Yan, M. Huang, B. Messer, J. H. Song, P. Yang, “Inorganic Semiconductor Nanowires: Rational Growth, Assembly, and Novel Properties”, Chemistry – A European Journal 8 (6), 1260-1268 (2002). Y. Yu, Q. Zhang, J. Xie, J. Y. Lee, “Engineering the Architectural Diversity of Heterogeneous Metallic Nanocrystals”, Natura Communications 4 (), 1454 (2013).
79 | CAS – Center for Advanced Sciences