Analisa Forward Modelling Dengan Berbagai Variasi Parameter.docx

ANALISA FORWARD MODELLING DENGAN BERBAGAI VARIASI PARAMETER DAN ANALISA INVERSI GEOLISTRIK

YOLANDA MUSTIKA BOHAL 03411640000012 KELAS A

TAHUN AJARAN 2018/2019 DEPARTEMEN TEKNIK GEOFISIKA FAKULTAS TEKNIK SIPIL LINGKUNGAN DAN KEBUMIAN

A. Horizontal Cylinder Script yang digunakan untuk melihat grafik model Horizontal Cylinder adalah clear all close all clc K=-125; %current dipole z=[2.5 5.5 8.5]; %kedalaman teta1=pi/5; %Sudut Polarisasi thd horinzontal 1 teta2=pi/3.6; %Sudut Polarisasi thd horinzontal 2 teta3=pi/2.5; %Sudut Polarisasi thd horinzontal 3 % Untuk Model Horizontal Cylinder, q=1 q=1; x=[-15:15]; %Posisi atau distance for i=1:length(x) for j=1:length(z) v1(i,j)=K*(x(i)*cos(teta1)+z(j)*sin(teta1))/(x(i).^2+z(j).^2).^q v2(i,j)=K*(x(i)*cos(teta2)+z(j)*sin(teta2))/(x(i).^2+z(j).^2).^q v3(i,j)=K*(x(i)*cos(teta3)+z(j)*sin(teta3))/(x(i).^2+z(j).^2).^q end end figure(1) plot(x,v1,'-') title('Horizontal Cylinder Anomaly Forward Modelling, Theta =36') xlabel('Jarak (m)') ylabel('Nilai SP (mV)') legend('z1=2.5','z2=5.5','z3=8.5') figure(2) plot(x,v2,'-') title('Horizontal Cylinder Anomaly Forward Modelling, Theta =50') xlabel('Jarak (m)') ylabel('Nilai SP (mV)') legend('z1=2.5','z2=5.5','z3=8.5') figure(3) plot(x,v3,'-') legend('z1=2.5','z2=5.5','z3=8.5') title('Horizontal Cylinder Forward Modelling, Theta =72') xlabel('Jarak (m)') ylabel('Nilai SP (mV)')

B. Vertical Cylinder Script yang digunakan untuk melihat grafik model Vertical Cylinder adalah clear all close all clc K=-125; %current dipole z=[2.5 5.5 8.5]; %kedalaman teta1=pi/5; %Sudut Polarisasi teta2=pi/3.6; %Sudut Polarisasi teta3=pi/2.5; %Sudut Polarisasi % Untuk Model Vertical Cylinder, q=0.5 q=0.5; x=[-15:15]; %Posisi atau distance for i=1:length(x) for j=1:length(z) v1(i,j)=K*(x(i)*cos(teta1)+z(j)*sin(teta1))/(x(i).^2+z(j).^2).^q v2(i,j)=K*(x(i)*cos(teta2)+z(j)*sin(teta2))/(x(i).^2+z(j).^2).^q v3(i,j)=K*(x(i)*cos(teta3)+z(j)*sin(teta3))/(x(i).^2+z(j).^2).^q end end figure(1) plot(x,v1,'-') title('Vertical Cylinder Anomaly Forward Modelling, Theta =36') xlabel('Jarak (m)') ylabel('Nilai SP (mV)') legend('z1=2.5','z2=5.5','z3=8.5') figure(2) plot(x,v2,'-') title('Vertical Cylinder Anomaly Forward Modelling, Theta =50') xlabel('Jarak (m)') ylabel('Nilai SP (mV)') legend('z1=2.5','z2=5.5','z3=8.5') figure(3) plot(x,v3,'-') legend('z1=2.5','z2=5.5','z3=8.5') title('Vertical Cylinder Forward Modelling, Theta =72') xlabel('Jarak (m)') ylabel('Nilai SP (mV)')

C. Sphere Script yang digunakan untuk melihat grafik model Sphere adalah clear all close all clc K=-125; %current dipole z=[2.5 5.5 8.5]; %kedalaman teta1=pi/5; %Sudut Polarisasi teta2=pi/3.6; %Sudut Polarisasi teta3=pi/2.5; %Sudut Polarisasi % Untuk Model Sphere, q=1.5 q=1.5; x=[-15:15]; %Posisi atau distance for i=1:length(x) for j=1:length(z) v1(i,j)=K*(x(i)*cos(teta1)+z(j)*sin(teta1))/(x(i).^2+z(j).^2).^q v2(i,j)=K*(x(i)*cos(teta2)+z(j)*sin(teta2))/(x(i).^2+z(j).^2).^q v3(i,j)=K*(x(i)*cos(teta3)+z(j)*sin(teta3))/(x(i).^2+z(j).^2).^q end end figure(1) plot(x,v1,'-') title('Sphere Anomaly Forward Modelling, Theta =36') xlabel('Jarak (m)') ylabel('Nilai SP (mV)') legend('z1=2.5','z2=5.5','z3=8.5') figure(2) plot(x,v2,'-') title('Sphere Anomaly Forward Modelling, Theta =50') xlabel('Jarak (m)') ylabel('Nilai SP (mV)') legend('z1=2.5','z2=5.5','z3=8.5') figure(3) plot(x,v3,'-') legend('z1=2.5','z2=5.5','z3=8.5') title('Sphere Forward Modelling, Theta =72') xlabel('Jarak (m)') ylabel('Nilai SP (mV)')

ANALISA A. Perbandingan Teta Terhadap Variasi q

Dari ketiga grafik diatas kita dapat melihat perbadingan dari ketiga model forward dengan variasi bentuk yaitu vertical cylinder, horizontal cylinder, dan sphere. Semakin besar teta atau sudut polarisasinya maka nilai SP nya semakin kecil dilihat ketika teta 36 dimulai dari range 50-100mV sedangkan model Sphere hanya dimulai di sekitara 0 mV lebih tidak sampai 1 mV. Dan kecuraman anomali yang ditimbulkan juga terlihat berbeda dimana pada model sphere setelah jarak nol dapat meningkat lagi secara tajam dan semakin kecilnya nilai q(bentuk) maka peningkatannya pun rendah seperti model horizontal cylinder dan vertical cylinder yang naik tidak terlalu significant. Selain itu pada bentuk vertical cylinder grafik dari masing-masing kedalamannya cenderung memiliki nilai yang tidak terlalu jauh dibandingkan dengan horizonal cylinder dan shepre. Karena kita dapat melihat bahwa melalui persamaan dalam membuat grafik q sebagai pangkat, sehingga semakin tinggi pangkatnya yakni q maka setiap grafik yang dihasilkan akan memiliki signifikansi yang cukup jauh satu sama lain. Kita dapat melihat juga bahwa range Vertical Cylinder cukup besar di range -150 hingga 100 mV dibandingkan Sphere yang hanya -15 hingga 2 mV jadi dapat ditarik kesimpulan semakin besar q maka range nilai SPnya semakin sempit.

B. Perbandingan Kedalaman Terhadap Variasi Teta

Ketiga gambar diatas secara umum dapat dilihat bahwa semakin besar tetanya maka semakin besar anomalinya diketahui melalui penunjaman teta 72 yang lebih tajam dan dalam dibandingkan dengan teta yang lebih kecil lainnya. Namun bisa dilihat perbedaan setiap kedalaman atau z akan semakin besar jika tetanya semakin kecil. Kita dapat melihatnya di gambar teta 36 dimana setiap grafik kedalamannya memliki nilai yang masing-masingnya cukup jauh. Sebaliknya pada teta yang besar seperti 72 antar grafik kedalamannya tidak terlalu memiliki nilai SP yang jauh karena grafiknya hampir berhimpit satu sama lainnya.

INVERSI DATA Script yang digunakan pada MATLAB untuk menginversi data % Inverse Modeling of Self Potential (SP) Method % Model Parameters and Variables: % q = Shape factor (dimensionless) % theta = The polarization angle between the axis of polarization % and the horizontal (in degree) % k = Electric current dipole moment (mV) % z = Depth (unit) % x = A discrete point along x-axis (unit) % v = Anomaly value at the origin, x0 (mV) % s = Observation point % vsp = SP value at observation point s, with xi=s (mV) % vsn = SP value at observation point s, with xi=-s (mV) % Y = Observed SP anomaly values at xi % x_0 = Distance at position xi=0 clear,clc q=1.5; %(sphere=1.5; horizontal cylinder=1; or vertical cylinder=0.5) s=6; %(changeable) xi=-12:12; %(depend on SP data) Y=[0.386802516 0.433437702 0.486124616 0.544187527 0.604956601 0.66152077 0.698347241 6.83E-01 0.556279455 2.20E-01 -4.34E-01 -1.412829378 -2.469135802 3.180992082 -3.334564962 -3.064648234 -2.618098265 -2.16175465 -1.765013908 1.442398158 -1.186869052 -0.985879967 -0.827384223 -0.701496354 -0.600576301]; %(depend on SP data) x_0=13; %(depend on SP data) v=Y(x_0); %(depend on SP data) vsp=Y(x_0+s); %(depend on SP data) vsn=Y(x_0-s); %(depend on SP data) %F and P are known numerical values from measured SP anomaly at x=0, x=s, %x=-s ---> The solution for q instead of the unknowns z and theta F=(vsp+vsn)/(2*v); P=(vsp-vsn)/(2*v); %Formula for W W=(s^(2*q-1))*(((xi.*P)+(s*F))./((xi.^2+(F^(1/q))*((s^2)-(xi.^2))).^q)); %The first derivative of W AA=(xi.*P+s*F)-s^(2*q-1); BB=(xi.^2+(s^2-xi.^2)-F^(1/q)).^q; CC=s^2; DD=xi.^2+(s^2-xi.^2)*F^(1/q); EE=(s^2-xi.^2)*(F^(1/q))*(log(F)); FF=q*(xi.^2+(s^2-xi.^2)*F^(1/q)); AB=AA./BB; CD=CC./DD; EF=EE./FF; Wd=(AB).*(log(CD)+(EF)); %Find q value GG=v*(W.^2); HH=rdivide((((s^2-xi.^2).*F^(1/q))*log(F)),(xi.^2+((F^(1/q))*(s^2-xi.^2)))); II=Y.*Wd; JJ=v*((W.^2).*log(rdivide(s^2,(xi.^2)+(F^(1/q))*(s^2-xi.^2)))); qc=rdivide(sum((GG).*(HH)),sum(II)-sum(JJ)) %Find z, theta, k, and e z2=(sqrt(((s^2)*(F^(1/qc)))/(1-(F^(1/qc))))) theta2=acotd((P/(s*F))*(sqrt(((s^2)*(F^(1/qc)))/(1-(F^(1/qc)))))) k2=v*((z2^(2*qc-1))/sind(theta2)) e=abs(qc-q) figure(1) plot(xi,Y,'*b')

Gambar dari script diatas dapat dilihat di gambar 1.

Gambar 1. Grafik Pemodelan Inversi Setlah dilakukan inversi pada data yang bisa kita lihat di script matlab, selanjutnya kita bisa menambahkan random noise 10% dan 20% kemudian diinversikan lagi. Dan script yang di gunakan ialah % Inverse Modeling of Self Potential (SP) Method % Model Parameters and Variables: % q = Shape factor (dimensionless) % theta = The polarization angle between the axis of polarization % and the horizontal (in degree) % k = Electric current dipole moment (mV) % z = Depth (unit) % x = A discrete point along x-axis (unit) % v = Anomaly value at the origin, x0 (mV) % s = Observation point % vsp = SP value at observation point s, with xi=s (mV) % vsn = SP value at observation point s, with xi=-s (mV) % Y = Observed SP anomaly values at xi % x_0 = Distance at position xi=0 clear,clc q=1.5; %(sphere=1.5; horizontal cylinder=1; or vertical cylinder=0.5) s=6; %(changeable) xi=-12:12; %(depend on SP data) Y=[0.386802516 0.433437702 0.486124616 0.544187527 0.604956601 0.66152077 0.698347241 6.83E-01 0.556279455 2.20E-01 -4.34E-01 -1.412829378 -2.469135802 3.180992082 -3.334564962 -3.064648234 -2.618098265 -2.16175465 -1.765013908 1.442398158 -1.186869052 -0.985879967 -0.827384223 -0.701496354 -0.600576301]; %(depend on SP data) u = rng r = 10/100*Y; Added Random Noise 10% rn = r.*(-1 + (1+1).*rand(1,25)); rng(u); Yn = Y+rn; x_0=19; %(depend on SP data) v=Yn(x_0); %(depend on SP data) vsp=Yn(x_0+s); %(depend on SP data) vsn=Yn(x_0-s); %(depend on SP data) %F and P are known numerical values from measured SP anomaly at x=0, x=s, %x=-s ---> The solution for q instead of the unknowns z and theta F=(vsp+vsn)/(2*v); P=(vsp-vsn)/(2*v); %Formula for W W=(s^(2*q-1))*(((xi.*P)+(s*F))./((xi.^2+(F^(1/q))*((s^2)-(xi.^2))).^q));

%The first derivative of W AA=(xi.*P+s*F)-s^(2*q-1); BB=(xi.^2+(s^2-xi.^2)-F^(1/q)).^q; CC=s^2; DD=xi.^2+(s^2-xi.^2)*F^(1/q); EE=(s^2-xi.^2)*(F^(1/q))*(log(F)); FF=q*(xi.^2+(s^2-xi.^2)*F^(1/q)); AB=AA./BB; CD=CC./DD; EF=EE./FF; Wd=(AB).*(log(CD)+(EF)); %Find q value GG=v*(W.^2); HH=rdivide((((s^2-xi.^2).*F^(1/q))*log(F)),(xi.^2+((F^(1/q))*(s^2-xi.^2)))); II=Y.*Wd; JJ=v*((W.^2).*log(rdivide(s^2,(xi.^2)+(F^(1/q))*(s^2-xi.^2)))); qc=rdivide(sum((GG).*(HH)),sum(II)-sum(JJ)) %Find z, theta, k, and e z2=(sqrt(((s^2)*(F^(1/qc)))/(1-(F^(1/qc))))) theta2=acotd((P/(s*F))*(sqrt(((s^2)*(F^(1/qc)))/(1-(F^(1/qc)))))) k2=v*((z2^(2*qc-1))/sind(theta2)) e=abs(qc-q) figure(1) plot(xi,Y,'-b');hold on plot(xi,Yn,'-r')

Gambar 2 (a) noise added 10% (grafik merah) (b) noise added 20% (grafik merah)

ANALISA Dari ketiga perbandingan model tersebut hasil command windownya adalah Model Normal qc = 1.4222 z2 =4.3106 theta2 = 31.0780 k2 = -70.8080 e = 0.0778

Model dengan Noise Added 10%

Model dengan Noise Added 20%

qc = 1.1057

qc = 1.0447

z2 = 15.6281

z2 = 14.4902

theta2 = -32.1576

theta2 = -34.0470

k2 = 89.7741

k2 = 54.4003

e = 0.3943

e = 0.4553

Dari gambar 2 dan tabel command window kita bisa melihat bahwa semakin besar noise yang ada maka model yang dihasilkan semakin tidak fit. Kita bisa membuktikannya melalui noise 20% memiliki bentuk grafik yang sangat berbeda dengan model normal (tidak berimpit) sedangkan noise dengan 10% tidak terlalu berbeda dengan grafik model yang asli. Selain itu nilai error yang dihasilkan cenderung semakin tinggi ketika noise yang ada semakin besar hal ini dibuktikan dari noise 10% besar e= 0.3943 sedangkan noise 20% besar e=0.4553.