METODE ELEKTROMAGNETIK Mencari Grafik hasil metode Magnetotelluric (MT) 1 Dimensi dari variasi data dengan MATLAB
Nama
: Dhika Rosari Purba
NIM
: 102120008
1. PENDAHULUAN Survey geofisika dilakukan untuk memperoleh informasi mengenai distribusi parameterparameter fisik bawah permukaan seperti kecepatan gelombang elastik, rapat massa, kemagnetan, kelistrikan dan lain lain. Dalam survey geofisika dengan menggunakan metoda elektromagnetik (EM) sifat fisik yang dapat ditemukan adalah konduktivitas atau resistivitas (tahanan-jenis) batuan. Dari hasil penelitian, kita dapat melihat adanya hubungan antara tahanan-jenis dengan porositas, kandungan fluida (air atau gas) dan temperatur formasi batuan. Pengaruh masing-masing faktor tersebut terhadap tahanan-jenis formasi batuan sangat kompleks. Oleh sebab itu, salah satu metoda untuk menentukan parameter fisis tersebut, dibutuhkan metoda EM yang dapat digunakan untuk keperluan eksplorasi sumber daya alam seperti mineral, minyak dan gas bumi, geotermal dan keperluan studi lingkuan yang lain. Metoda magnetotellurik (MT) merupakan salah satu metoda eksplorasi geofisika yang bersifat pasif memanfaatkan medan elektromagnetik alam sebagai sumber gelombang atau energy untuk mengetahui struktur tahanan jenis bawah permukaan. Medan EM di alam ditimbulkan oleh berbagai proses fisik yang cukup kompleks sehingga spektrum frekuensinya sangat lebar (10-5Hz – 104Hz). Metode Magnetotellurik bergantung pada fenomena listrik-magnet terutama terhadap konduktivitas bumi, yang berasal dari aktivitas petir(>1Hz), resonansi lapisan ionosfer bumi (<1Hz) dan bintik hitam matahari (<<1Hz). Prinsip kerja metode magnetotellurik adalah proses induksi elektromagnetik yang terjadi pada anomaly bawah permukaan. Medan EM yang menembus bawah permukaan akan menginduksi anomaly konduktif bawah permukaan bumi yang menghasilkan medan listrik E dan medan magnetic sekunder H (arus eddy) yang nanti direkam oleh alat magnetotellurik. Informasi mengenai konduktivitas medium yang terkandung dalam data MT dapat diperoleh dari penyelesaian persamaan Maxwell.
2. PENURUNAN PERSAMAAN MAXWELL ∇xE=−
!𝐁 !! !𝐃
(1a)
∇ x H = j + !!
(1b)
∇⋅D = q
(1c)
∇⋅B = 0
(1d)
dimana : E : medan listrik (Volt/m) B : fluks atau induksi magnetik(Weber/m2atauTesla) H : medan magnet (Ampere/m) j : rapat arus (Ampere/m2) D : perpindahan listrik (Coulomb/m2) q : rapat muatan listrik (Coulomb/m3) Hubungan antara intensitas medan dengan fluks yang terjadi pada medium dinyatakan oleh persamaan : B = 𝜇H
(2a)
D = 𝜀E
(2b) 𝐄
j = 𝜎E = !
(2c)
dimana : 𝜇: permeabilitas magnetik (Henry/m) 𝜀: permitivitas listrik (Farad/m) : 𝜎: konduktivitas (Ohm-1/m atau Siemens/m) 𝜌: tahanan-jenis (Ohm.m) Untuk menyederhanakan masalah,kita asumsikan pendekatan yang dilakukan sifat fisik medium diasumsikan tidak bervariasi terhadap waktu dan posisi (homogen isotropik) sehingga akumulasi muatan seperti pada persamaan (1c) tidak akan terjadi, sehingga : ∇xE=-
!!𝐇
(3a )
!! !𝐄
∇ x H = 𝜎E + 𝜀 !!
(3b)
∇⋅E = 0
(3c)
∇⋅H = 0
(3d)
Dalam persamaan (3) tampak pada seluruh persamaan hanya menggunakan variable E dan H. Dengan operasi curl pada persamaan (3a) dan (3b) didapat : !𝑬
!!𝑬
∇ × ∇ × 𝑬 = −𝜇𝜎 !" − 𝜇𝜀 !" ! ∇ × ∇ × 𝑯 = −𝜇𝜎
!𝑯 !"
− 𝜇𝜀
!!𝑯 !" !
(4a) (4b)
Dengan menggunakan identitas vector : ∇ × ∇ × 𝑥 = ∇∇. 𝑥 − ∇! 𝑥 ; didapat Persamaan gelombang Helmholtz untuk E dan H : 𝜕𝑬
∇! 𝑬 = 𝜇𝜎 𝜕𝑡 + 𝜇𝜀
2
𝜕 𝑬 𝜕𝑡
𝜕𝑯
2
(5a)
2
𝜕 𝑯
∇! 𝑯 = 𝜇𝜎 𝜕𝑡 + 𝜇𝜀 2 𝜕𝑡
(5b)
Dimana E dan H adalah suatu fungsi dalam ruang dan waktu :
E(r,t) = Eo (r) eiwt
(6a)
H(r,t) = Ho (r) eiwt
(6b)
Dimana : Eo = amplitude medan listrik Ho = amplitudo medan magnet w= frekuensi gelombang elektromagnetik
Sehingga persamaan Helmholtz (Persamaan 5 sebelumnya) menjadi :
E = ( iw𝜇𝜎 – w2𝜇𝜀 ) E
(7a)
H = (iw𝜇𝜎 - w2𝜇𝜀) H
(7b)
Pada kondisi umum dalam eksplorasi geofisika suku yang mengandung 𝜀 (perpindahan listrik)
dapat diapaikan terhadap suku yang mengandung 𝜎 (konduksi listrik). Sebab nilai w𝜇𝜎 >> w2𝜀. Sehingga : E = ( iw𝜇𝜎 )E = k2 E
(8a)
H = ( iw𝜇𝜎 )H = k2 H
(8b)
Dimana k = ± iw𝜇𝑜𝜎 = ± 𝛼 + 𝑖𝛽 Dengan 𝛼 = 𝛽 =
!!𝑜! !
Dengan melakukan pendekatan model bumi yang halfspace homogenous isotropic, dimana diskontinuitas tahanan jenis hanya terjadi pada batas udara dengan bumi. Sehingga komponen horizontal medan listrik dan medan magnet hanya bervariasi terhadap kedalaman, didapat persamaan : ! ! !" !! !
= 𝑘 ! 𝐸𝑥
(9)
Solusi elementernya adalah : Ex= A e-kz + Bekz
(11a)
Ex = A 𝑒 !!"# 𝑒 !!" + B 𝑒 !"# 𝑒 !"
(11b)
Dimana x dan y adalah sumbu koordinat kartesian dengan z adalah kedalaman positif (vertical ke bawah) Eksponensial sendiri yang mengandung : -
komponen imajiner : menyatakan variasi sinusoidal gelombang EM terhadap kedalaman komponen real : menyatakan factor atenuasi menurut kedalaman
Dengan dekomposisi persamaan 3a, 6b dan 11a dari medan magnet berikut : !
Hy = − !"#$
!"# !"
!
Hy = − !"#! A e-kz + Bekz
(12)
Sehingga persamaan 12 merupakan solusi persamaan difusi untuk medan magnet pd persamaan 8b. Untuk bumi homogen, nilai B adalah = 0. Sebab sumber medan EM bersifat ekstem dimana nilainya akan =0 di kedalaman tak hingga. Sedangkan A nilainya bernilai, sebab menyatakan atenuasi gelombang EM terhadap kedalaman (z+ arahnya ke bawah). -
Impedansi : Adalah perbandingan anatara medan listrik dan medan magnet yang saling tegak lurus. Dari persamaan 11 dan 12 : !" Zxy = !" = 𝑖𝑤𝜇𝑜𝜌 (13a) !"
Zyx = !" = − 𝑖𝑤𝜇𝑜𝜌
(13b)
Untuk impedansi bumi Homogen,impedansinya adalah impedansi intrinstik (Zo) dimana : !"
Zo = Zxy = !" = 𝑖𝑤𝜇𝑜𝜌 Sehingga didapat resisitivitas dan fasanya : !
1. Resisitivitas (𝜌) = !"𝑜 |𝑍𝑜|!
(14a)
!" !"
2. Fasa = tan-1 [ !" !" ] -
(14b)
Skin Effect dan Penetration Depth 1. Skin effect adalah eksponensial dari gelombang EM yang teratenuasi oleh kedalaman lapisan tanah 2. Skin depth(δ ) didefinisikan sebagai kedalaman medan pada suatu medium homogen dimana amplitude gelombang EM telah tereduksi menjadi 1/e. A(δ) = Ae-αδ
δ=
!
= !"𝑜 2𝜌
!! !"#
(15)
TEKNIK PEMODELAN Pada kali ini akan dicari kurva resistivitas semu dan fasa terhadap frekuensi berdasarkan pemodelan pada data berikut : 𝝆1 𝝆1 30000 ohm-‐m h = 65 m 𝝆2 0.9 ohm-‐m h = 40 m
𝝆𝟑 4200 ohm-‐m h = 120 m
Gambar 1. Model 1D lapisan bawah permukaan yang terdiri 3 lapisan horizontal homogen, dimana lapisan terakhir adalah half-space dengan ketebalan tak berhinga
Teknik yang digunakan dalam memodelkan kurva resistivitas semu terhadap frekuensi dan kurva fasa terhadap frekuensi adalah dengan pemodelan forward (forward modeling). Teknik ini dilakukan dengan menghitung respon dari suatu model untuk dibandingkan dengan data impedansi (resistivitas semu dan fasa) pengamatan. Dengan cara coba-coba, dapat diperoleh suatu model yang responnya paling cocok dengan data, sehingga model tersebut dapat dianggap mewakili kondisi bawah permukaan. Model 1D berupa model berlapis horizontal, dimana model ini terdiri dari beberapa lapisan. Setiap lapisan memiliki tahanan jenis yang homogen. Dalam hal ini, parameter model 1D adalah tahanan jenis dan ketebalan tiap lapisan. Kode Program Menggunakan Matlab 1. Koding Pemrograman untuk Perumusan Fungsi Resistivitas Semu dan Fasa function [resistivitassemu, fasa] = fungsi(rho, tebal, frek) mu = 4*pi*1E-7; %Deklarasi Permeabilitas ruang vakum(H/m) w = 2 * pi * frek; % Deklarasi Frekuensi Sudut (Radian) n=length(rho); %Banyaknya Layer yaitu lapisan j, lapisan 1, lapisan 2 dan lapisan 3. impedansi = zeros(n,1); % mendeklarasikan impedansi Zn = sqrt(sqrt(-1)*w*mu*rho(n)); % Zn - Impedansi Lapisan Dasar untuk memulai rekursi impedansi(n) = Zn; % lakukan iterasi untuk menghitung impedansi dan parameter fisis lainnya untuk tiap-tiap lapisan. for j = n-1:-1:1 rho = rho(j); ketebal = tebal(j); induksij = sqrt(sqrt(-1)* (w * mu * (1/rho))); % induksijParameterInduksi impij = induksij * rho; % impijj - Impedansi Intrinsik ej = exp(-2*ketebalan*induksij); % ej - Faktor Eksponensial x = impedansi(j + 1);
rj = (impijj - x)/(impijj + x); % rj - Koefisien Refleksi re = rj*ej; % re - Earth R.C. implj = impij * ((1 - re)/(1 + re)); % implj - Impedansi Layar impedansi(j) = Zj; end Z = impedansi(1); absZ = abs(Z); resistivitassemu = (absZ * absZ)/(mu * w); % deklarasi persamaan untuk mendapat nilai resistivitas semu fasa = atan2(imag(Z),real(Z)); % deklarasi persamaan untuk mendapat nilai fasa
2. Lalu buat koding untuk memplot hasil pada kodingan sebelumnya dengan terlebih dahulu menginput nilai resistivitas pada tiap lapisan, ketebalan dan frekuensi gelombang sehingga dapat membentuk grafik resistivitas semu terhadap frekuensi clear all; clc; figure(1); close; rho = [30000 0.9 4200]; %nilai resistivitas yang diinput tebal = [65 40 120]; %nilai resistivitas yang diinput log = -5:0.1:4; frek = 10.^log; resistivitassemu = zeros(length(frek),1); fasa = zeros(length(frek),1); for i = 1:length(frek); [resistivitassemu(i) fasa(i)] = fungsi(rho, tebal, frek(i)); end %koding untuk plot grafik untuk Resistivitas Semu dan Fasa figure(2); subplot(2,1,1); loglog(frek,resistivitassemu, 'bs','LineWidth',1,'MarkerSize',2,'MarkerFaceColor','b'); grid on ylabel('Resistivitas Semu (Ohm m)'); xlabel('Frekuensi (Hz)'); title('Grafik Resistivitas Semu(Ohm m) Terhadap Frekuensi(Hz)'); subplot(2,1,2); loglog(frekuensi,fasa, 'bs','LineWidth',1,'MarkerSize',2,'MarkerFaceColor','b'); grid on ylabel('Fasa (rad)'); xlabel('Frekuensi (Hz)'); title('Grafik Fasa(Ohm m) Terhadap Frekuensi(Hz)');
3. Setelah di-run tampak hasil plot grafik resistivitas semu terhadap frekuensi yaitu :
Dalam skala lebih besar :