PEMODELAN ELEMEN BATAS UNTUK KASUS ELEKTROMAGNETIK 2D Imran Hilman Mohammad*
ABSTRAK PEMODELAN KASUS RAMBATAN GELOMBANG ELEKTROMAGNETIK 2D MENGGUNAKAN METODE ELEMEN BATAS. Metode elemen batas telah dikenal sebagai salah satu metode numerik dalam memecahkan persamaan diferensial parsial. Metode ini memiliki keunikan dibandingkan metode konvensional seperti metode elemen hingga dan metode beda hingga yang membutuhkan diskretisasi seluruh domain pemodelan, yaitu diskretisasi hanya perlu dilakukan pada bidang-bidang batas pemodelan. Keunikan ini menjadikan metode elemen batas sangat efisien dalam pemodelan, karena hanya membutuhkan memori komputasi yang sedikit dibandingkan metode konvensional. Hal ini juga menyebabkan waktu komputasi yang dibutuhkan sangat sedikit dibandingkan metode konvensional. Makalah ini menjelaskan tentang pemodelan gelombang elektromagnetik 2D pada aplikasi kebumian dengan menggunakan metode elemen batas. Respon pemodelan dengan metode elemen batas dibandingkan dengan respon hasil pemodelan elemen hingga untuk melihat performa metode elemen batas. Kata kunci: Metode elemen batas, metode elemen hingga, pemodelan elektromagnetik, gelombang bidang ,resistivitas semu
ABSTRACT THE BOUNDARY ELEMENT MODELING FOR 2D ELECTROMAGNETIC CASE. The Boundary Element Method (BEM) has widely known as the numerical method for solving partial differential equation. The method is unique compared to konventional methods such as finite element or finite difference methods,since BEM needs the discretization only at the domain boundaries. The uniqueness leads the method very efficient, since only low computational memory is needed, and thus yields low time consuming for computational process. This paper presents the modeling of 2D electromagnetic wave in geophysics application using BEM. The modeling responses obtained by BEM are compared to those obtained by finite element method (FEM) to show the performance of BEM. Keywords: Boundary element method, finite element method, electromagnetic modeling, plane wave, apparent resistivity
PENDAHULUAN Teknik elemen batas (boundary element method) merupakan teknik memecahkan persamaan diferensial yang relatif baru dibandingkan teknik pemecahan menggunakan metode beda hingga maupun elemen hingga. Pemodelan elemen batas untuk kasus elektromagnetik pada bidang geofisika telah dilakukan semenjak dekade 80an. Xu dan Zhao (1987) memodelkan kasus MT 2D dengan permukaan flat dan *
Program Studi Geofisika – Universitas Padjadjaran Jatinangor, e-mail:
[email protected]
223
Lokakarya Komputasi dalam Sains dan Teknologi Nuklir, 10 Oktober 2012 (223-230)
homogen. Pemodelan efek topografi pada bumi homogen dilakukan Xu dan Zhou (1997) dengan menggunakan elemen batas linear. Untuk kasus bumi non homogen, pemodelan elemen batas dilakukan oleh Mohammad dkk (2011) yang mensimulasikan respon 2D pada kasus bumi berlapis, kasus kontak vertikal dan kasus anomali konduktif. Respon gelombang bidang akibat efek topografi untuk kasus bumi homogen dan bumi berlapis dengan elemen batas konstan dilakukan oleh Mohammad (2011). Pada makalah ini akan dimodelkan respon elektromagnetik akibat rambatan gelombang bidang untuk kasus kontak vertikal dengan lapisan penutup yang relatif dalam dan kontras resistivitas yang relatif dekat. Respon yang dihasilkan berupa resistivitas semu yang dihitung berdasarkan persamaan Cagniard (1953) akan dibandingkan dengan respon serupa hasil pemodelan elemen hingga, untuk melihat performa metode elemen batas.
FORMULASI ELEMEN BATAS UNTUK KASUS ELEKTROMAGNETIK 2D Persamaan gelombang bidang pada kasus geofisika dapat dinyatakan sebagai berikut (Xu dan Zhou 1997, Mohammad dkk 2011):
(∇ (∇
2 2
) )
+ μεω 2 − iμωσ E = 0 ⎫⎪ ⎬ + μεω 2 − iμωσ H = 0⎪⎭
(1)
Untuk kasus elektromagnetik dalam prospeksi geofisika, dimana proses difusi mendominasi pada frekuensi rendah sehingga σ >> ωε , maka persamaan (1) dapat disederhanakan dalam bentuk persamaan Helmholtz
(∇
2
)
+ k2 φ = 0
(2)
dengan k 2 ≈ −iωμσ Persamaan (2) merupakan persamaan diferensial parsial yang dapat dipecahkan menggunakan metode elemen batas. Dengan menggunakan teknik pembobotan residual, persamaan (2) dikalikan dengan suatu fungsi uji dan diintegralkan untuk seluruh domain pemodelan
∫ φ (∇ *
2
)
+ k 2 φdΩ = 0
(3)
Ω
Integrasi parsial persamaan (3) menghasilkan:
∫ φ (∇ *
Ω
224
2
)
(
)
(
)
+ k 2 φ d Ω = ∫ ∇ ∇ φφ * d Ω − ∫ ∇ φ ∇ φ * d Ω Ω
Ω
(4)
Pemodelan Kasus Rambatan Gelombang Elektromagnetik 2D Menggunakan ... (Imran Hilman Mohammad)
Dengan menggunakan teorema Gauss untuk pertama ruas kiri persamaan (4), maka akan didapatkan ∂φ * * 2 2 * (5) ∫Ω φ ∇ + k φdΩ = ∫Γ ∂nφ dΓ − Ω∫ ∇φ ∇φ dΩ
(
)
(
)
Integral parsial untuk suku kedua ruas sebelah kanan menghasilkan: ⎛ * ∂φ ∂φ * ⎞ ∂ 2φ * * 2 2 ⎜ ⎟ ∇ + Ω = − φ φ φ φ φ k d d Γ + ∫ ∫Γ ⎜⎝ ∂n ∂n ⎟⎠ Ω∫ ∂n dΩ Ω
(
)
(6)
Persamaan (2) menyatakan bahwa suku sebelah kiri dari persamaan (6) bernilai nol. Sehingga suku sebelah kanan dapat diuraikan sebagai berikut: − ∫φ Ω
⎛ ∂φ ∂ 2φ * ∂φ * ⎞ ⎟ dΓ dΩ = ∫ ⎜⎜ φ * −φ ∂n ∂n ∂n ⎟⎠ Γ⎝
(7)
Fungsi uji yang diperlukan dalam pembobotan persamaan (2) di atas memenuhi hubungan berikut:
(∇
2
+ k 2 )φ * = −δ (x, ξ )
(8)
dengan menggunakan sifat distribusi Dirac:
∫ f ( x)δ ( x, ξ )dx =
f (ξ )
(9)
Γ
maka suku sebelah kiri persamaan (7) dapat kita tuliskan:
∂ 2φ * − ∫φ dΩ = u (ξ ) ∂n Ω
(10)
sehingga secara lengkap kita tuliskan persamaan (7):
⎛ ∂φ ∂φ * ⎞ ⎟⎟ dΓ u (ξ ) = ∫ ⎜⎜ φ * −φ n n ∂ ∂ ⎠ Γ⎝
(11)
dan
⎛ ∂φ * ⎛ ∂φ * ∂φ ⎞ ϖ * ∂φ ⎞ ⎟⎟dΓ − ∫ ⎜⎜ φ ⎟dΓ −φ −φ* φ p = − ∫ ⎜⎜ φ 2π ∂n ∂n ⎠ ∂n ∂n ⎟⎠ Γ ⎝ Γ⎝
(12)
S
225
Lokakarya Komputasi dalam Sains dan Teknologi Nuklir, 10 Oktober 2012 (223-230)
Pada persamaan (12),
ϖ φ p menyatakan nilai medan yang dihitung pada titik p 2π
dikalikan harga sudut yang dibentuk antara titik p dengan sisi batas tempat titik p tersebut berada. Persamaan (12) merupakan persamaan integral batas yang menyatakan solusi nilai medan pada suatu titik pada batas pemodelan. Persamaan ini selanjutnya didiskretisasi menjadi elemen-elemen yang merepresentasikan integral batas dan dipecahkan menggunakan permasalahan syarat batas yang ditinjau. Hasil pemodelan adalah berupa nilai medan listrik dan medan magnet pada bidang batas domain pemodelan. Dengan menggunakan nilai medan listrik dan medan magnet tersebut, dapat dicari nilai resistivitas semu pada setiap titik pada batas pemodelan menggunakan persamaan Cagniard (Xu dan Zhou 1997, Mohammad dkk 2011).
HASIL DAN ANALISIS Pada makalah ini, model yang digunakan adalah suatu kontak vertikal, yaitu kontras resistivitas yang bersinggungan secara vertikal dengan lapisan penutup yang cukup tebal (500 meter). Kontras resistivitas antara dua bidang yang bersinggungan relatif kecil (100 Ohm dan 10 Ohm). Model ini mendekati kasus suatu patahan dalam yang tertutup suatu lapisan tebal, dengan batuan pembentuk patahan di kedua sisi memiliki kontras resistivitas yang relatif dekat. Pemodelan dilakukan baik untuk modus TE dimana medan listrik diasumsikan sejajar bidang patahan (dalam hal ini arah bidang batas pada kontak vertikal, yaitu tegak lurus bidang gambar) dan untuk modus TM dimana medan listrik diasumsikan tegak lurus bidang patahan. Respon pemodelan dihitung menggunakan metode elemen batas dan elemen hingga untuk melihat performa metode elemen batas. Frekuensi yang digunakan pada pemodelan ini adalah 0.1 Hz, 1 Hz dan 10 Hz. Respon elemen hingga dihitung menggunakan kode yang dikembangkan Srigutomo dan Sutarno (1998). Model kontak vertikal beserta responnya disajikan pada gambar-gambar berikut.
Gambar 1. Model kontak vertikal 226
Pemodelan Kasus Rambatan Gelombang Elektromagnetik 2D Menggunakan ... (Imran Hilman Mohammad)
Gambar 2. Perbandingan respon resistivitas semu modus TM antara metode elemen batas (garis putus-putus) dengan metode elemen hingga (garis tegas).
Gambar 3. Perbandingan respon resistivitas semu modus TE antara metode elemen batas (garis putus-putus) dengan metode elemen hingga (garis tegas). Pada Gambar 1 diasumsikan suatu kontak vertikal terdapat pada kedalaman 500 meter di bawah lapisan penutup. Bidang kontak vertikal berada pada titik 0 meter, 227
Lokakarya Komputasi dalam Sains dan Teknologi Nuklir, 10 Oktober 2012 (223-230)
dengan batuan dengan resistivitas 100 Ohm meter berada pada sisi sebelah kiri dan batuan dengan resistivitas 10 Ohm meter berada pada sisi sebelah kanan. Pada modus TM, arah medan listrik tegak lurus kontak vertikal menyebabkan timbulnya arus permukaan pada bidang patahan, sehingga menimbulkan respon yang menunjukkan kontras resistivitas. Gambar 2 menunjukkan respon yang dihasilkan oleh pemodelan elemen batas dan elemen hingga untuk modus TM. Pada frekuensi 0.1 Hz, terlihat bagaimana pemodelan dengan elemen batas dan elemen hingga sama-sama menunjukkan penurunan nilai resistivitas semu ketika mendekati bidang kontak vertikal. Namun nilai resistivitas semu respon elemen hingga menunjukkan hasil yang mendekati kontras resistivitas yang dimodelkan. Sementara respon elemen batas terlihat menurun dengan lambat untuk mendekati harga resistivitas batuan sebelah kanan. Hal serupa juga terlihat pada frekuensi 1 Hz. Untuk frekuensi 10 Hz, respon elemen batas dan elemen hingga menunjukkan kecenderungan yang mirip. Hal ini dikarenakan pada frekuensi 10 Hz, skindepth gelombang EM berada pada kedalaman dangkal sehingga lebih mencirikan resistivitas lapisan penutup ketimbang keberadaan kontak vertikal. Pemodelan modus TE menggunakan asumsi medan listrik yang digunakan sejajar arah bidang kontak vertikal. Berdasarkan hukum Ohm, keberadaan medan listrik yang sejajar bidang kontak akan mengakibatkan respon yang dihasilkan bersifat smooth dalam menunjukkan keberadaan kontak vertikal. Gambar 3 menunjukkan respon elemen batas yang dihasilkan memiliki kontras yang relatif lebih baik daripada respon elemen hingga untuk frekuensi 0.1 Hz. Pada frekuensi 1 Hz, kecenderungan respon elemen hingga menunjukkan adanya perbedaan resistivitas secara lateral, sementara respon elemen batas cenderung konstan. Sementara pada frekuensi 10 Hz, respon elemen batas dan elemen hingga menunjukkan kecenderungan yang sama, mencirikan kondisi pada lapisan penutup. Dari segi waktu komputasi, waktu rata-rata yang dibutuhkan pada elemen batas untuk menghitung respon pada satu harga frekuensi adalah 58 detik untuk modus TM dan 100 detik untuk modus TE. Sementara respon elemen hingga untuk satu harga frekuensi rata-rata membutuhkan waktu sekitar 300 detik untuk modus TM dan 500 detik untuk modus TE. Secara umum untuk kasus kontak vertikal ini, walaupun memiliki waktu komputasi yang lebih singkat, namun hasil pemodelan elemen batas belum menunjukkan konsistensi sebagaimana yang diharapkan. Pemodelan elemen batas masih memerlukan penyempurnaan untuk menghasilkan respon dengan presisi tinggi dan waktu komputasi yang rendah.
KESIMPULAN Pada makalah ini telah dijelaskan pemodelan respon elektromagnetik pada kasus kontak vertikal dengan lapisan penutup relatif tebal dan kontras resistivitas relatif dekat. Pemodelan elemen batas memiliki waktu komputasi yang relatif lebih singkat dibandingkan metode konvensional elemen hingga. Namun demikian respon
228
Pemodelan Kasus Rambatan Gelombang Elektromagnetik 2D Menggunakan ... (Imran Hilman Mohammad)
yang dihasilkan masih belum memenuhi hasil yang diharapkan. Dengan demikian masih dibutuhkan pengembangan lebih lanjut untuk pemodelan elemen batas agar dihasilkan skema pemodelan yang memiliki presisi tinggi, akurat dengan waktu komputasi yang singkat.
DAFTAR PUSTAKA
1. MOHAMMAD, I.H, SRIGUTOMO, W., dan SUTARNO, D., “Pemodelan Elektromagnetik 2D Menggunakan Metode Elemen batas”, Jurnal Matematika dan Sains, 16,(2011), 82-94. 2. MOHAMMAD, I.H, SRIGUTOMO, W., dan SUTARNO, D., “Topographic Effect Modeling of 2D MT Responses using Boundary Element Method,” Publikasi dalam Proceeding, Bandung, International Conferences of Physics and Its Application (ICPAP) (2011) 3. SRIGUTOMO, W., AND SUTARNO, D., “2D Electromagnetic Modeling using Finite Element Method: Application in Magnetotelluric Case,” Kontribusi Fisika Indonesia, 9 (2) (1998) 55-65. 4. XU, S, Z.,AND ZHAO, S.K., “Two Dimensional Magnetotelluric Modeling by the Boundary Element Method,” Journal of Geomagnetism and Geoelectricity 39,(1987), 677-698. 5. XU, S, Z., AND ZHOU, H., “Modeling of 2D Terrain Effect on MT by the Boundary Element Method”, Geophysical Prospecting (1997) 45, 931-943. DISKUSI SETIANTO Speed-up metode elemen batas terhadap metode konvensional (data dari computer) IMRAN HILMAN Speed up metode elemen batas terhadap metode konvensional relative berbeda untuk setiap kasus yang ditinjau. Pada kasus pemodelan elektromagnetik ini, sebagai perbandingan, untuk kasus bumi homogen dengan frekuensi 16 Hz, metode elemen batas membutuhkan waktu komputasi 90 detik, sedangkan metode elemen hingga membutuhkan waktu sekitar 2.000 detik. Pada kasus-kasus yang lain waktu komputasi meningkat seiring peningkatan kompleksitas sistem geologi yang ditinjau.
229
Lokakarya Komputasi dalam Sains dan Teknologi Nuklir, 10 Oktober 2012 (223-230)
SAHRUL HIDAYAT 1. Apa perbedaan modus TE/TM, apakah terdapat perbedaan interaksi medan E dan M dengan batuan? 2. Apakah yang dimaksud dengan syarat batas itu, apakah batas antar lapisan atau batas domain? IMRAN HILMAN 1. Modus TE dan modus TM dalam penjalaran gelombang EM muncul akibat adanya srtuktur bumi 2D, terutama pada kasus hukum Ohm, terjadi perbedaan antara pelajaran gelombang EM dengan arah medan listrik sejajar strike dan tegak lurus strike. Perbedaan ini terkait munculnya akumulasi arus permukaan pada bidang batas resistivitas pada kasus tegak lurus strike. 2. Syarat batas berlaku pada batas antar lapisan. Namun untuk memudahkan pemodelan, masing-masing lapisan dimodelkan sebagai satu domain tertentu dengan spesifikasi resistivitas tertentu.
DAFTAR RIWAYAT HIDUP
1. 2. 3. 4.
Nama Instansi / Unit Kerja Pekerjaan / Jabatan Riwayat Pendidikan
: Imran Hilman Mohammad : FMIPA Unpad : Dosen Fisika Unpad : - S1 Fisika ITB - S2 Fisika ITB - S3 Fisika ITB 5. Pengalaman Kerja :6. Organisasi Profesional : Himpunan Fisikawan Indonesia 7. Publikasi Ilmiah yang pernah disajikan/diterbitkan: 1. Mohammad, I.H, Srigutomo, W., dan Sutarno, D. (2011): Pemodelan Elektromagnetik 2D Menggunakan Metode Elemen batas, Jurnal Matematika dan Sains, 16, hal. 82-94. 2. Mohammad, I.H., Srigutomo, W., Sutarno, D dan Sumintadireja, P. (2012): The Modeling of 2D Controlled Source Audio Magnetotelluric (CSAMT) Responses Using Finite Element Method, Journal of Electromagnetic Analysis and Application (JEMAA), 4, hal. 293-304. 3. Mohammad, I.H., Srigutomo, W., Sutarno, D. (2011): Topographic Effect Modeling of 2D MT Responses using Boundary Element Method, Publikasi dalam Proceeding, Bandung, International Conferences of Physics and Its Application (ICPAP) 2011. 4. Mohammad I.H., Srigutomo, W., Sutarno, D., Sumintadireja D. (2012) 2D Modeling of Controlled Source Audio Magnetotelluric Responses, Bandung, Asian Physics Symposium (APS) 2012.
230
