6
BAB 2 LANDASAN TEORI
2.1.
Persamaan Integral Pengintegrasian
(integration),
bersama
dengan
pendiferensiasian
(differentiation), merupakan konsep matematika yang menjadi jantung dari kalkulus. Berdasarkan Microsoft Encarta Premium 2006 Dictionary Tools, mengintegrasikan (to integrate) berarti “membawa seluruhnya, bagian per bagian, menjadi sebuah bentuk utuh; menyatukan; untuk mengindikasikan jumlah total ...” Secara matematik, pengintegrasian dilambangkan dengan: b
I = ∫ f ( x) dx
Pers. (2.1)
a
yang berarti integral dari fungsi f(x) dengan variabel bebas x, dievaluasi berdasarkan batas x = a dan x = b. Fungsi f(x) dalam pers. (2.1) disebut integran (integrand). Berdasarkan definisi yang diberikan oleh kamus bahasa, “arti” dari pers. (2.1) adalah nilai total atau penjumlahan (summation) dari f(x) dx dalam jangkauan x = a hingga b. Berdasarkan faktanya, simbol
∫
sebenarnya merupakan sebuah huruf kapital
S yang diberi sentuhan gaya yang ditujukan untuk menandakan hubungan tertutup antara pengintegrasian dan penjumlahan. Gambar 2.1 berikut merupakan sebuah penggambaran secara grafik dari konsep
tersebut. Untuk fungsi yang terletak di atas sumbu x, integral yang dinyatakan oleh pers. (2.1) mengacu kepada area di bawah kurva f(x) antara x= a dan b.
7
Gambar 2.1 Penggambaran secara grafik dari integral f(x) antara batas x = a dan b. Integral tersebut ekivalen dengan daerah di bawah kurva
Fungsi yang diintegrasikan akan memiliki satu dari tiga bentuk khusus berikut ini: •
Sebuah fungsi kontinu sederhana, seperti sebuah polinomial, sebuah eksponensial, atau sebuah fungsi trigonometri.
•
Sebuah fungsi kontinu yang rumit, yang sulit atau tidak mungkin diintegralkan secara langsung..
•
Sebuah fungsi yang ditabulasikan, di mana nilai-nilai dari x dan f(x) diberikan pada sejumlah titik tertentu (diskrit), sebagaimana seringkali terjadi pada kasus dengan data lapangan atau data hasil eksperimen. Pada kasus pertama, integral dari sebuah fungsi sederhana mungkin dapat
dievaluasi secara analitik menggunakan kalkulus. Untuk kasus kedua, solusi analitik seringkali tidak praktis, dan terkadang mustahil, untuk diperoleh. Dalam hal ini, demikian juga dengan kasus ketiga dari data diskrit, metode pendugaan/perkiraan harus
8 digunakan. Pendekatan berorientasi visual digunakan untuk mengintegrasikan data tertabulasi dan fungsi yang rumit pada masa sebelum komputer ditemukan. Sebuah pendekatan berdasarkan intuisi yang sederhana adalah dengan menempatkan fungsi tersebut pada sebuah kisi-kisi (grid) (Gambar 2.2) dan menghitung jumlah kotak-kotak yang digunakan sebagai penduga dari daerah tersebut. Jumlah tersebut dikalikan dengan daerah di mana setiap kotak memberikan sebuah perkiraan kasar mengenai daerah total di bawah kurva. Perkiraan ini dapat diperhalus, dengan usaha tambahan, menggunakan sebuah kisi-kisi yang lebih baik.
Gambar 2.2 Penggunaan dari sebuah kisi-kisi (grid) untuk menaksir sebuah integral
9 Pendekatan yang masuk akal lainnya adalah dengan membagi daerah tersebut menjadi bagian-bagian berbentuk vertikal, atau strips, dengan tinggi yang setara dengan nilai fungsi pada titik tengah setiap bagian (Gambar 2.3). Daerah dari persegi panjang tersebut dapat kemudian dihitung dan dijumlahkan untuk memperkirakan.jumlah total dari daerah tersebut. Pada pendekatan ini, diasumsikan bahwa nilai pada titik tengah setiap bagian memberikan sebuah pendugaan yang sah dari rata-rata tinggi dari fungsi untuk setiap strip. Sebagaimana dengan metode kisi-kisi, perkiraan yang diperhalus adalah mungkin dengan menggunakan lebih banyak (dan lebih tipis) strips untuk menaksir integral tersebut.
Gambar 2.3 Penggunaan dari persegi panjang atau strips untuk menaksir sebuah integral
2.2.
Metode dan Integrasi Numerik Walaupun pendekatan sederhana yang demikian memiliki sarana untuk
pendugaan yang cepat, teknik numerik alternatif juga tersedia untuk tujuan yang sama. Tidak mengejutkan, cara yang paling sederhana dari metode ini mirip dalam hal jiwa dari teknik nonkomputer.
10 Integrasi numerik atau metode kuadratur (quadrature) tersedia untuk memperoleh integral-integral. Metode-metode ini, yang mana lebih mudah untuk diimplementasikan dibandingkan dengan pendekatan kisi-kisi, sangat mirip jiwanya dengan metode strip. Yaitu, tinggi fungsi dikalikan dengan lebar strip dan dijumlahkan untuk memperkirakan integral tersebut. Bagaimanapun, melalui pemilihan yang cerdas dari faktor pembobot (weighting factor), perkiraan yang dihasilkan dapat dibuat lebih akurat dibandingkan dengan metode strip yang sederhana. Sebagaimana dalam metode strip sederhana, integrasi numerik menggunakan data pada titik-titik tertentu. Karena informasi berbentuk tabel sudah berbentuk demikian, maka hal itu dapat cocok dengan banyak dari pendekatan-pendekatan numerik. Sebagaimana ditunjukkan oleh gambar 2.4, tabel
tersebut
dievaluasi
dengan
dapat sebuah
kemudian metode
numerik.
Gambar 2.4 Aplikasi dari sebuah metode integrasi numerik: (a) Fungsi kontinu yang rumit. (b). Tabel nilai diskrit dari f(x) yang dihasilkan dari sebuah fungsi. (c). Penggunaan salah satu metode numerik (metode strip) untuk menaksir integral yang berdasarkan kepada titik-titik diskrit. Untuk sebuah fungsi berbentuk tabel, data sudah dalam bentuk tabel (b); untuk itu, langkah (a) tidak diperlukan
11 2.3.
Integrasi Romberg Pada subbab 2.1., diketahui bahwa fungsi-fungsi yang akan diintegrasikan secara
numerik memiliki dua bentuk yang khas: sebuah tabel berisi nilai-nilai atau sebuah fungsi. Bentuk dari data memiliki pengaruh terhadap pendekatan yang dapat digunakan untuk mengevaluasi integralnya. Untuk informasi berbentuk tabel (yang tidak dibahas dan digunakan dalam skripsi ini), kita dibatasi oleh jumlah dari titik-titik yang diberikan. Berlawanan dengan hal itu, Jika sebuah fungsi tersedia, maka kita dapat menghasilkan sebanyak mungkin nilai dari f(x) yang dibutuhkan untuk mencapai tingkat akurasi yang dapat diterima (seperti gambar 2.4). Integrasi Romberg merupakan teknik yang digunakan dalam integrasi numerik untuk menganalisis kasus di mana fungsi yang akan diintegrasikan tersedia. Teknik ini memiliki keunggulan untuk menghasilkan nilai-nilai dari fungsi untuk mengembangkan skema yang efisien bagi pengintegrasian secara numerik. Integrasi Romberg didasarkan pada ekstrapolasi Richardson (Richardson’s extrapolation), yang mana merupakan metode untuk mengkombinasikan dua perkiraan integral secara numerik untuk memperoleh nilai ketiga, yang lebih akurat. Integrasi Romberg merupakan algoritma komputasi untuk mengimplementasikan ekstrapolasi Richardson melalui sebuah cara yang sangat efisien. Teknik ini bersifat rekursif dan dapat digunakan untuk menghasilkan sebuah perkiraan integral dalam batas toleransi kesalahan (error tolerance) yang sudah ditentukan terlebih dahulu. Integrasi Romberg didasarkan pada aplikasi beturut-turut (rekursif) dari aturan trapesium (trapezoidal rule). Walau bagaimanapun, melalui manipulasi secara
12 matematik, hasil yang superior dapat dicapai dengan sedikit usaha. 2.3.1. Ekstrapolasi Richardson Ekstrapolasi Richardson merupakan metode yang menggunakan dua perkiraan dari sebuah integral untuk mengkomputasi pendugaan ketiga, yang lebih akurat. Perkiraan dan kesalahan (error) yang diasosiasikan dengan aturan trapesium
multi-aplikasi
(multiple-application
trapezoidal
rule)
dapat
digambarkan secara umum sebagai
I = I ( h) + E ( h)
Pers. (2.2)
di mana I adalah nilai yang sebenarnya dari integral tersebut, I(h) adalah pendugaan dari sebuah aturan trapesium dengan aplikasi bersegmen n dengan lebar langkahnya h = (b – a)/n, dan E(h) adalah kesalahan pemotongan (truncation error). Jika kita membuat dua perkiraan yang berbeda menggunakan lebar langkah h1 dan h2 dan memiliki nilai yang sebenarnya untuk error-nya, I (h1 ) + E (h1 ) = I (h2 ) + E (h2 )
Pers. (2.3)
Sekarang ingat bahwa error dari aturan trapesium multi-aplikasi dapat diperkirakan sebagai
Ea = −
(b − a) 3 f" 12n 2
Dengan n = (b - a)/h menjadi
Pers. (2.4)
13
E≅−
b−a 2 h f" 12
Pers. (2.5)
Jika diasumsikan bahwa f " adalah konstan, tidak dipengaruhi oleh lebar langkah. Pers. (2.5) dapat digunakan untuk menentukan bahwa rasio dari kedua error adalah E (h1 ) h12 ≅ E (h2 ) h22
Pers. (2.6)
Perhitungan ini memiliki efek penting dalam penghilangan bagian f " dari perhitungan. Dengan demikian, hal itu memungkinkan kita untuk menggunakan informasi yang dinyatakan oleh pers. (2.5) tanpa terlebih dahulu mengetahui turunan kedua dari fungsi tersebut. Untuk melakukan hal itu, kita ubah pers. (2.6) menjadi ⎛h ⎞ E (h1 ) ≅ E (h2 )⎜⎜ 1 ⎟⎟ ⎝ h2 ⎠
2
Pers. (2.7)
yang mana dapat disubstitusikan ke dalam pers (2.3): 2
⎛h ⎞ I (h1 ) + E (h2 )⎜⎜ 1 ⎟⎟ ≅ I (h2 ) + E (h2 ) ⎝ h2 ⎠
Pers. (2.8)
yang dapat disederhanakan menjadi
E (h2 ) ≅
I (h1 ) − I (h2 ) 1− (h1 / h2 ) 2
Pers. (2.9)
Hasilnya, kita telah mengembangkan sebuah perkiraan dari error pemotongan
14 dalam hal pendugaan integral dan lebar langkahnya. Perkiraan ini dapat disubstitusikan ke dalam I = I (h2 ) + E (h2 )
Pers. (2.10)
untuk menghasilkan sebuah pendugaan integral yang telah diperbaiki:
I ≅ I (h2 ) +
1 [ I (h2 ) − I (h1 )] (h1 / h2 ) 2 −1
Pers. (2.11)
Hal ini menunjukkan (Ralston dan Rabinowitz, 1978) bahwa error dari pendugaan ini adalah O(h4). Hasilnya, kita telah mengkombinasikan dua perkiraan aturan trapesium dari O(h2) untuk menghasilkan sebuah penduga yang baru dari O(h4). Untuk kasus khusus di mana intervalnya dibagi dua (h2 = h1/2), persamaan ini menjadi
I ≅ I (h2 ) +
1 [ I (h2 ) − I (h1 )] 2 −1 2
Pers. (2.12)
atau, dibentuk menjadi, 4 1 I ≅ I (h2 ) − I (h1 ) 3 3
Pers. (2.13)
Pers. (2.4) memberikan sebuah cara untuk mengkombinasikan dua aplikasi
dari aturan trapesium dengan error O(h2) untuk menghitung penduga ketiga dengan error O(h4). Pendekatan ini adalah sebuah subhimpunan dari sebuah metode yang lebih umum untuk mengkombinasikan integral-integral dari O(h4) pada basis dari tiga penduga aturan trapesium. Kedua penduga tersebut dapat,
15 dalam perubahannya, dikombinasikan untuk menghasilkan sebuah nilai yang lebih baik dengan O(h6). Untuk kasus khusus di mana penduga trapesium yang asli didasarkan kepada pembagian terus-menerus dari lebar langkah, persamaan yang digunakan untuk keakuratan O(h6) adalah
I≅
16 1 Im − Il 15 15
Pers. (2.14)
di mana Im dan Il masing-masing adalah penduga yang lebih (more) dan kurang (less) akurat. Dengan cara yang sama, dua hasil O(h6) dapat dikombinasikan untuk menghitung sebuah integral O(h8) menggunakan
I≅
64 1 Im − Il 63 63
Pers. (2.15)
2.3.2. Algoritma Integrasi Romberg
Perhatikan bahwa koefisien dari masing-masing persamaan ekstrapolasi (Pers. (2.13), Pers. (2.14), dan Pers. (2.15)) bertambah hingga mendekati satu. Hasilnya, mereka merepresentasikan faktor-faktor pembobot yang, sejalan dengan meningkatnya akurasi, menempatkan bobot yang secara relatif lebih besar pada penduga integral superior. Formulasi-formulasi tersebut dapat diekspresikan dalam bentuk yang umum yang cocok untuk diimplementasikan pada komputer:
I j ,k ≅
4 k −1 I j +1,k −1 − I j ,k −1 4 k −1 − 1
Pers. (2.16)
di mana Ij + 1,k – 1 dan Ij,k – 1, masing-masing adalah integral yang lebih dan kurang
16 akurat, dan Ij,k adalah integral yang diperbaiki. Indeks k menandakan tingkat dari pengintegrasian, di mana k = 1 mengacu kepada penduga aturan trapesium yang asli, k = 2 mengacu kepada O(h4), k = 3 kepada O(h6), dan seterusnya. Indeks j digunakan untuk membedakan antara penduga yang lebih akurat (j + 1) dan penduga yang kurang akurat (j). Sebagai contoh, untuk k = 2 dan j = 1, pers. 2.16) menjadi
I1, 2 ≅
4 I 2,1 − I1,1 3
Pers. (2.17)
yang ekivalen dengan pers. (2.13).
Gambar 2.5 Penggambaran secara grafik dari barisan penduga integral f(x) = 0,2 + 25x – 200x2 + 675x3 – 900x4 + 400x5 dari a = 0 hingga b = 0,8, yang dihasilkan menggunakan integrasi Romberg. (a) Iterasi pertama. (b) Iterasi kedua. (c) Iterasi ketiga
Bentuk umum yang direpresentasikan oleh pers. (2.16) merupakan hasil karya Romberg, dan aplikasi sistematik untuk mengevaluasi integralnya dikenal dengan nama Integrasi Romberg. Gambar 2.5 merupakan penggambaran secara grafik dari barisan penduga integral yang dihasilkan menggunakan pendekatan ini. Setiap matriks mengacu kepada sebuah iterasi tunggal. Kolom pertama mengandung evaluasi terhadap aturan trapesium yang menggambarkan Ij,i, di
17 mana j = 1 adalah untuk sebuah aplikasi bersegmen tunggal (lebar langkahnya adalah b – a), j = 2 adalah untuk sebuah aplikasi dengan segmen dua (lebar langkahnya adalah (b – a)/2), j = 3 adalah untuk sebuah aplikasi bersegmen empat (lebar langkahnya adalah (b – a)/4), dan seterusnya. Kolom matriks lainnya dihasilkan oleh pengaplikasian pers. (2.16) secara sistematis untuk memperoleh penduga yang lebih baik dari suatu integral secara berturut-turut. Sebagai contoh, iterasi pertama (Gambar 2.5a) melibatkan perhitungan penduga aturan trapesium bersegmen satu dan dua (I1,1 dan I2,1). Pers. (2.16) kemudian digunakan untuk mencari elemen I1,2 = 1.367467, yang merupakan sebuah error dari O(h4). Sekarang, kita harus memeriksa untuk menentukan apakah hasil ini cukup untuk keperluan kita. Sebagaimana dengan metode pendugaan lainnya, sebuah
kriteria
pengakhiran
(termination)
atau
penghentian
(stopping)
dibutuhkan untuk menentukan keakuratan dari suatu hasil. Salah satu metode yang dapat digunakan untuk tujuan ini adalah sebagai berikut
εa =
perkiraan saat ini − perkiraan sebelumnya 100% perkiraan saat ini
Pers. (2.18)
atau
εa =
I1,k − I 2,k −1 I1,k
100%
Pers. (2.19)
di mana ε a adalah sebuah penduga dari persen error relatif. Hasilnya,
18 sebagaimana telah dilakukan sebelumnya pada proses iterasi lainnya, kita membandingkan penduga yang baru dengan nilai sebelumnya. Ketika pergantian antara nilai yang lama dan baru sebagaimana direpresentasikan oleh ε a , berada di bawah kriteria error ε s yang telah ditentukan sebelumnya, maka perhitungan dihentikan. Untuk gambar 2.5a, evaluasi ini menandakan perubahan sebanyak 21,8% terhadap hasil iterasi pertama. Tugas dari iterasi kedua (Gambar 2.5b) adalah untuk memperoleh penduga O(h6)—I1,3. Untuk itu, sebuah penduga aturan trapesium tambahan, I3,1 = 1,4848, ditentukan.
Kemudian
penduga
tersebut
dikombinasikan
dengan
I2,1
menggunakan pers. (2.16) untuk menghasilkan I2.2 = 1,623467. Hasilnya, pada gilirannya, dikombinasikan dengan I1,2 untuk menghasilkan I1,3 = 1,640533. Pers. (2.17)
dapat
diaplikasikan
untuk
menentukan
bahwa
hasil
tersebut
merepresentasikan sebuah perubahan sebanyak 1.0% ketika dibandingkan dengan hasil I1,2 sebelumnya. Iterasi ketiga (Gambar 2.5c) melanjutkan proses tersebut dengan cara yang sama. Pada kasus ini, sebuah penduga trapesium ditambahkan ke dalam kolom pertama, dan kemudian pers. (2.16) diaplikasikan untuk menghitung secara berturut-turut integral yang lebih akurat sepanjang diagonal yang lebih rendah. Setelah hanya tiga kali iterasi, karena mengevaluasi sebuah polinomial berordo lima, hasilnya (I1,4 = 1,640533) adalah tepat. Integrasi Romberg lebih efisien dibandingkan dengan aturan trapesium dan aturan Simpson (tidak dibahas dalam skripsi ini). Sebagai contoh, untuk
19 menentukan integral yang ditunjukkan oleh f(x) = 0,2 + 25x – 200x2 + 675x3 – 900x4 + 400x5 dari a = 0 hingga b = 0,8, aturan Simpson 1/3 membutuhkan sebuah aplikasi dengan segmen sebanyak 256 untuk menghasilkan sebuah penduga dari 1,640533. Penaksiran yang lebih baik mustahil dilakukan karena error yang disebabkan oleh pembulatan angka (di belakang koma). Berlawanan dengan hal itu, integrasi Romberg mencapai hasil yang tepat (hingga tujuh angka yang signifikan) dengan menggunakan aturan trapesium bersegmen satu, dua, empat, dan delapan; yaitu, dengan pengevaluasian fungsi sebanyak 15 kali saja. Perlu diingat bahwa integrasi Romberg dirancang untuk kasus di mana fungsi yang akan diintegrasikan sudah diketahui. Hal ini disebabkan karena pengetahuan fungsi tersebut mengijinkan pengevaluasian yang dibutuhkan untuk implementasi pertama dari aturan trapesium. 2.4.
Kuadratur Gauss dan Gauss-Legendre
Kuadratur Gauss merupakan pendekatan metode numerik yang digunakan untuk mengatasi kelemahan yang dimiliki oleh persamaan Newton-Cotes, yaitu di mana persamaan Newton-Cotes ini (contoh: aturan trapesium) melakukan pendugaan terhadap suatu integral dengan berdasarkan kepada fungsi bernilai genap. Akibatnya, nilai dasar yang digunakan dalam persamaan ini sudah ditentukan sebelumnya atau bersifat tetap. Sebagai contoh, seperti pada gambar 2.6, aturan trapesium didasarkan kepada pengambilan daerah pendugaan di bawah garis lurus yang menghubungkan nilai fungsi pada akhir interval integrasi. Rumus yang digunakan untuk menghitung area ini adalah
20 I ≅ (b − a )
f (a ) + f (b ) 2
Pers. (2.20)
di mana a dan b adalah batas integrasi dan b – a adalah jarak/lebar dari interval integrasi. Karena aturan trapesium harus melewati titik akhir, maka akan muncul kasus-kasus seperti gambar 2.6 di mana rumus tersebut menghasilkan error yang tinggi.
Gambar 2.6 (a) Penggambaran secara grafik dari aturan trapesium sebagai daerah di bawah garis lurus yang menghubungkan titik-titik akhir. (b) Sebuah pendugaan integral yang lebih baik yang diperoleh dengan cara mengambil area di bawah garis lurus yang melewati dua titik yang ada di antaranya. Dengan memposisikan titik-titik tersebut dengan bijak, maka error positif dan negatif menjadi berimbang, dan sebuah pendugaan integral yang lebih baik akan dihasilkan.
Sekarang, seandainya halangan titik dasar yang tetap tersebut dihilangkan dan kita dengan bebas menempatkan titik-titik tersebut dengan bijak, maka kita dapat menghasilkan sebuah garis lurus yang akan menyimbangkan error positif dan negatif. Seperti gambar 2.6, kita akan melihat sebuah pendugaan integral yang telah diperbaiki. Kuadratur Gauss adalah nama untuk sebuah teknik untuk mengimplementasikan strategi tersebut. Rumus kuadratur Gauss ini disebut rumus Gauss-Legendre. Rumus
21 Gauss-Legendre dikembangkan dengan menggunakan teknik yang dapat menurunkan rumus integrasi numerik seperti aturan trapesium yang menggunakan metode koefisien tak tertentu. 2.4.1. Metode Koefisien Tak Tertentu
Metode tertentu
koefisien
memberikan
pendekatan
yang
tak sebuah
memiliki
perangkat untuk menurunkan teknik integrasi lainnya seperti kuadratur Gauss.
Gambar 2.7 Dua integral yang dievaluasi dengan tepat oleh aturan trapesium: (a) sebuah konstan dan (b) sebuah garis lurus.
Sebagai contoh, pers. (2.20) dinyatakan sebagai I ≅ c 0 f (a ) + c1 f (b )
Pers. (2.21)
di mana c adalah konstan. Ingat bahwa aturan trapesium akan menghasilkan hasil yang tepat ketika fungsi yang diintegrasikan adalah sebuah konstan atau garis lurus. Dua persamaan sederhana yang menggambarkan kasus tersebut adalah y = 1 dan y = x. Keduanya diilustrasikan pada gambar 2.7. Karena itu, persamaan
22 berikut akan menjadi
c0 + c1 = ∫
(b − a ) / 2
(b − a ) / 2
1 dx
dan
− c0
b−a b − a (b − a ) / 2 + c1 =∫ x dx (b − a ) / 2 2 2
atau, dengan mengevaluasikan integral, c 0 + c1 = b − a
dan
− c0
b−a b−a + c1 =0 2 2
Berikut adalah dua persamaan dengan dua variabel yang tidak diketahui yang dapat diselesaikan dengan
c 0 = c1 =
b−a 2
yang mana, ketika disubstitusikan kembali ke dalam pers. (2.21), sehingga memberikan
I=
b−a b−a f (a ) + f (b ) 2 2
yang sama dengan aturan trapesium.
23 2.4.2. Penurunan Rumus Gauss-Legendre Dua Titik
Seperti halnya penurunan aturan trapesium di atas, tujuan kuadratur Gauss adalah untuk menentukan koefisien dari sebuah persamaan dari bentuk I ≅ c 0 f ( x 0 ) + c1 f (x1 )
Pers. (2.22)
di mana c adalah koefisien yang tidak diketahui. Walau bagaimanapun, berlawanan dengan aturan trapesium yang menggunakan titik-titik akhir a dan b yang tetap, argumen fungsi x0 dan x1 yang tidak tetap pada titik akhir, namun tidak diketahui (Gambar 2.8). Dengan demikian, kita sekarang memiliki empat variabel yang tidak diketahui yang harus dievaluasi, dan akibatnya, kita membutuhkan empat kondisi untuk menentukan variabel yang tidak diketahui tersebut secara tepat.
Gambar 2.8 Penggambaran secara grafik dari variabel x0 dan x1 untuk pengintegrasian dengan menggunakan kuadratur Gauss
Seperti halnya aturan trapesium, kita dapat memperoleh dua dari kondisi tersebut dengan mengasumsikan bahwa pers. (2.22) cocok dengan integral dari sebuah konstan dan sebuah fungsi linier secara tepat. Kemudian, untuk sampai
24 pada kedua kondisi lainnya, kita hanya melanjutkan argumen ini dengan mengasumsikan bahwa kondisi ini juga cocok dengan integral dari sebuah fungsi parabolik (y = x2) dan fungsi kubik (y = x3). Dengan melakukan hal itu, kita menentukan keempat variabel yang tidak diketahui dan sebagai pertukarannya menurunkan sebuah rumus integrasi linier dua titik yang tepat untuk kubikkubik. Keempat persamaan yang akan diselesaikan adalah c 0 f ( x 0 ) + c1 f ( x1 ) = ∫ 1 dx = 2
Pers. (2.23)
c 0 f ( x 0 ) + c1 f ( x1 ) = ∫ x dx = 0
Pers. (2.24)
1
−1
1
−1
c 0 f ( x 0 ) + c1 f ( x1 ) = ∫ x 2 dx =
2 3
Pers. (2.25)
c 0 f ( x 0 ) + c1 f ( x1 ) = ∫ x 3 dx = 0
Pers. (2.26)
1
−1
1
−1
Pers. (2.24) hingga pers. (2.26)
dapat diselesaikan secara simultan dengan
c0 = c1 = 1
x0 = −
x1 =
1 3
Pers. (2.27)
= − 0.5773503...
1 = 0.5773503... 3
Pers. (2.28)
Pers. (2.29)
yang mana bisa disubstitusikan ke dalam pers. (2.22) untuk menghasilkan rumus Gauss-Legendre dua titik
25 ⎛ −1 ⎞ I ≅ f ⎜⎜ ⎟⎟ + ⎝ 3⎠
⎛ 1 ⎞ f ⎜⎜ ⎟⎟ ⎝ 3⎠
Pers. (2.30)
Dengan demikian, kita sampai pada hasil yang menarik bahwa penambahan yang sederhana nilai fungsi pada
x =1 / 3
dan
− 1/ 3
menghasilkan sebuah pendugaan integral yang memiliki keakuratan ordo ketiga. Ingat bahwa batas integrasi pada pers. (2.23) hingga (2.26) adalah dari -1 hingga 1. Hal ini dilakukan untuk menyederhanakan perhitungan dan untuk rumus tersebut seumum mungkin. Sebuah penggantian sederhana terhadap variabel dapat digunakan untuk menerjemahkan batas lainnya dari integrasi ke dalam bentuk ini. Hal ini dilakukan dengan mengasumsikan bahwa sebuah variabel baru xd dihubungkan dengan variabel asli x dengan cara yang linier, sebagaimana dalam x = a0 + a1 xd
Pers (2.31)
Jika batas bawah, x = a, sama dengan xd = -1, nilai ini dapat disubstitusikan ke dalam pers. (2.31) untuk menghasilkan a = a0 + a1 (− 1)
Pers (2.32)
Dengan cara yang sama, batas atas, x = b, sama dengan xd = 1, untuk memberikan b = a0 + a1 (1) Pers. (2.32) dan (2.33) dapat diselesaikan secara simultan menjadi
Pers (2.33)
26 a0 =
b+a 2
Pers (2.34)
a1 =
b−a 2
Pers (2.35)
dan
yang dapat disubstitusikan ke dalam pers. (2.31) untuk menghasilkan
x=
(b + a )+ (b − a )xd 2
Pers. (2.36)
Persamaan ini dapat diturunkan untuk memberikan
dx =
b−a dxd 2
Pers. (2.36) dan (2.37)
Pers. (2.37)
dapat disubstitusikan untuk x dan dx, masing-masing, dalam
persamaan yang akan diintegrasikan. Penyubstitusian ini secara efektif mengubah interval integrasi tanpa mengubah nilai dari integrasi tersebut. 2.4.3. Rumus dengan Jumlah Titik yang Lebih Banyak
Lebih jauh dibandingkan rumus dua titik dideskripsikan dalam subbab sebelumnya, versi jumlah titik yang lebih banyak dapat dikembangkan dalam bentuk umum I ≅ c 0 f ( x 0 ) + c1 f ( x1 ) + ... + d n − 1 f ( x n − 1 )
Pers. (2.38)
di mana n adalah jumlah dari titik. Nilai untuk c dan x untuk (≤) enam titik dirangkum dalam tabel 2.1.
27 Karena kuadratur Gauss memerlukan evaluasi fungsi pada titik dengan jarak yang tidak seragam pada interval integrasi, hal ini tidak ditujukan untuk kasus fungsi yang tidak diketahui. Oleh karena itu, metode ini tidak cocok untuk permasalahan teknik yang berhubungan dengan data berbentuk tabel. Bagaimanapun, ketika fungsinya diketahui, keefisienannya dapat menjadi sebuah keuntungan. Hal ini sangat benar ketika sejumlah besar evaluasi integral harus dilakukan.
Tabel 2.1 Faktor pembobot c dan argumen fungsi x yang digunakan dalam Gauss-Legendre
2.4.4. Analisis Error untuk Gauss-Legendre
Error dari rumus Gauss-Legendre dispesifikasikan secara umum oleh (Carnahan et al., 1969)
28 2 2 n + 3 [(n + 1)!]
4
Et =
(2n + 3)[(2n + 2)!]
3
f
(2 n + 2 )
(ξ )
Pers. (2.39)
di mana n adalah jumlah titik dikurangi satu dan f(2n + 2)(ξ) adalah turunan ke (2n + 2) dari fungsi setelah perubahan terhadap variabel dengan ξ terletak di suatu tempat pada interval -1 hingga 1. Perbandingan terhadap pers. (2.38) dengan tabel 2.1 menandakan kelebihan dari kuadratur Gauss dibandingkan dengan rumus
Newton-Cotes, menyediakan turunan dengan ordo yang lebih tinggi yang tidak meningkat secara nyata dengan meningkatnya n. 2.5.
Metode Simulasi Monte Carlo
Metode simulasi Monte Carlo merupakan sebuah kelas dari algoritma komputer untuk
melakukan
simulasi
terhadap
sifat-sifat
dari
sistem
kefisikaan
dan
kematematikaan. Metode simulasi Monte Carlo sangat penting dalam perkomputasian fisika dan bidang terapan yang berhubungan. Metode ini sudah teruji dalam menyelesaikan persamaan integro-diferensial yang mendefinisikan bidang radiasi, dan metode ini telah digunakan dalam perhitungan global illumination yang menghasilkan citra yang terlihat nyata dari model buatan berdimensi tiga, untuk penerapan pada permainan video, arsitektur, desain, film yang dihasilkan dengan komputer, efek khusus dalam film, bidang bisnis, ekonomi, dan bidang lainnya. Metode ini terutama sangat berguna untuk sistem pembelajaran dengan sejumlah besar derajat bebas berpasangan, seperti cairan, material tidak teratur, dan benda padat yang berpasangan sangat kuat. Lebih jauh, metode Monte Carlo berguna untuk memodelkan fenomena dengan ketidakpastian yang signifikan pada input, seperti penghitungan dari resiko pada bisnis.
29 Yang menarik, metode Monte Carlo tidak memerlukan bilangan acak yang sejati agar berjalan dengan baik. Banyak dari teknik-teknik yang sangat berguna menggunakan barisan bilangan acak palsu (pseudo-random sequences) yang bersifat deterministik (tanpa keacakan), membuatnya mudah untuk diuji dan melakukan simulasi ulangan. Satu-satunya kualifikasi yang dibutuhkan untuk membuat simulasi yang baik adalah agar barisan bilangan pseudo-random itu terlihat “cukup acak” dalam beberapa hal.
Gambar 2.9 Hasil plotting menggunakan R Language terhadap 1000 bilangan pseudo-random
Proses membuat bilangan tersebut tampak acak disebut proses pseudo-random, sebuah proses yang terlihat acak tetapi tidak. Barisan pseudo-random secara khusus memperlihatkan keacakan secara statistik yang dihasilkan oleh sebuah proses komputasi yang seluruhnya bersifat deterministik/tanpa keacakan. Apa arti dari hal ini tergantung pada penerapannya, namun secara khusus hal tersebut harus melewati sebuah tes
30 statistik berseri. Menguji apakah bilangan-bilangan tersebut terdistribusi secara seragam atau mengikuti distribusi lainnya yang diinginkan ketika jumlah elemen yang cukup besar dari barisan tersebut dipertimbangkan, merupakan salah satu yang paling sederhana dan yang paling umum. Sebuah algoritma Monte Carlo merupakan metode Monte Carlo yang digunakan untuk mencari solusi dari permasalahan matematik (yang mungkin memiliki banyak variabel) yang tidak mudah untuk dihitung dengan integral kalkulus atau metode numerik lainnya. Keefisienannya secara relatif meningkat dibandingkan dengan metode numerik lainnya bila dimensi dari permasalahan tersebut meningkat. Karena perulangan dari algoritma dan sejumlah besar perhitungan dilibatkan, Monte Carlo merupakan sebuah metode yang cocok untuk melakukan perhitungan dengan komputer, dengan menggunakan banyak teknik dari simulasi komputer. Bagian ini akan menjelaskan tentang teori probabilitas di balik integrasi Monte Carlo. Oleh karena itu akan digunakan contoh hypercube C = [0,1)d dengan unit dimensi r d. Integrannya adalah sebuah fungsi f( x ) yang diasumsikan real dan tidak negatif dan
tentunya dapat diintegrasikan terhadap C. Dengan demikian akan didefinisikan sebagai berikut: r r I m = ∫ f ( x)m d d x
m = 1, 2, 3, ...
Pers (2.40)
C
sehingga I1 adalah integral yang dibutuhkan. Perlu diingat bahwa Im tidak selalu harus tertentu (finite) untuk m ≥ 2. Dalam Monte Carlo, diasumsikan bahwa titik-titik integrasi sebanyak N dipilih secara bebas, dari distribusi peluang seragam terhadap C. Hal itu
31 r r r berarti bahwa himpunan titik X = { x1 , x2 , ..., x N } di mana integrasi didasarkan,
diasumsikan menjadi sebuah anggota yang khas dari kelompok himpunan titik-titik tersebut, sedemikian rupa sehingga distribusi peluang N yang dikombinasikan menunjuk ke kelompok himpunan titik-titik tersebut bersifat seragam, bebas, dan bersitribusi identik. PN ( x1 , x 2 , ..., x N , ) = 1
Pers. (2.41)
Rata-rata dari kelompok himpunan di atas akan diambil, diasumsikan bahwa r sebuah himpunan titik dari X telah dihasilkan, dan nilai-nilai dari integran f( x ) sudah r dihitung. Hal ini dilambangkan dengan fj ≡ f( x j), j = 1, 2 , ... , N. Dari sini kita dapat
menghitung analogi diskrit dari integral Im yang dapat dihitung dalam waktu linier: N
Sm = ∑ ( f j )m
Pers. (2.42)
j =1
Penduga Monte Carlo dari integral tersebut menjadi
E1 =
1 S1 N
Pers. (2.43)
Nilai harapan E1 dari himpunan kelompok titik-titik di atas adalah
E1 =
1 N
∑ i
fi =
r
∫ f (x) d
d
r x = I1
Pers. (2.44)
Cd
yang mana merupakan integral yang dibutuhkan. Ini adalah basis dari metode Monte Carlo. Kegunaannya terlihat jika kita menghitung ragam dari E1
32
σ ( E1 ) 2 = E12 − E1
2
=
1 ( I 2 − I 12 ) N
Pers. (2.45)
Karena nilai di atas menurun sejalan dengan N-1, metode Monte Carlo sebenarnya konvergen untuk N yang besar. Ingat bahwa O(N0), hubungan E12 , dan E1
2
saling
menghilangkan: hal ini merupakan sebuah fenomena biasa dalam pendugaan ragam jenis ini.1 Ragam σ ( E1 ) 2 diduga dengan penduga error ordo pertama
E2 =
1 1 S 2 − 3 S12 2 N N
Pers. (2.46)
sehingga kita memiliki
σ (E 2 ) 2 =
(
1 I 4 − 4 I 3 I 1 − I 22 + 8 I 2 I 12 − 4 I 12 3 N
)
Pers. (2.47)
yang mana penduganya adalah
E4 =
(
1 N 3 S 4 − 4 N 2 S 3 S1 − N 2 S 22 + 8 NS 2 S12 − 4S14 7 N
)
Pers. (2.48)
yang juga dapat dihitung dalam waktu linier; maka kita memiliki
(
E 4 = σ (E 2 ) + O N − 4 2
)
Pers. (2.49)
Galat dari metode Monte Carlo ini adalah
1
Perlu diketahui bahwa apa yang diduga adalah rata-rata dari squared error, bukan error itu sendiri, dan
penguadratan dan perataan tidak menyebabkan perubahan. Karena itu, inilah alasan mengapa penduga berordo dua dikatakan relevan.
33
εn =
∫
C
f ( x) dx −
d
1 N
N
∑ f (x ) i =1
Pers. (2.50)
i
Jika nilai N besar, maka ⎡ ⎛ ⎞ ⎤ ε n = ⎢ ∫ ⎜⎜ f ( x ) − ∫ f ( x ) dx ⎟⎟ dx ⎥ Cd ⎠ ⎦⎥ ⎣⎢c s ⎝ 2.6.
1/ 2
Pers. (2.51)
Metode Simulasi Quasi-Monte Carlo
Dalam analisis numerik, sebuah metode quasi-Monte Carlo merupakan metode untuk perhitungan integral (atau permasalahan lainnya) yang didasarkan pada barisan bilangan dengan ketidakcocokan/ketidaksesuaian yang rendah (low-discrepancy sequences). Hal inilah yang membedakan metode quasi-Monte Carlo dengan metode
Monte Carlo yang menggunakan barisan bilangan acak palsu (pseudo-random).
Gambar 2.10 Hasil plotting menggunakan R Language terhadap 1000 bilangan quasirandom menurut barisan Sobol (Sobol sequence)
34 Monte Carlo dan quasi-Monte Carlo dinyatakan dengan cara yang sama. Masalahnya adalah untuk menduga fungsi f sebagai rataan dari fungsi yang dievaluasi pada sehimpunan titik-titik x1, ... , xN.
∫
f (u ) du ≈
Cd
1 N ∑ f ( x) N i =1
Pers. (2.52)
di mana C adalah kubus dengan dimensi d, Cd = [0,1] x ... x [0,1]. Sehingga masingmasing xi adalah sebuah vektor dari elemen sebanyak d. Dalam metode Monte Carlo, xi adalah bilangan-bilangan pseudo-random. Dalam metode quasi-Monte Carlo, himpunan tersebut adalah akibat dari sebuah barisan low-discrepancy atau disebut dengan bilangan quasi-random. Error pendugaan dari metode di atas dibatasi dengan sebuah syarat yang
proporsional terhadap ketidaksesuaian (discrepancy) dari himpunan x1, ... , xN, oleh pertidaksamaan Koksma-Hlawka.
1 N
N
∑ f (x ) → ∫ f (x ) dx i =1
N
≤ V HK ( f ) • D N∗ ( x N )
Pers. (2.53)
Cd
di mana DN∗ ( x N ) adalah discrepancy dan VHK ( f ) adalah variasi (Hardy Krause). Discrepancy adalah ukuran simpangan dari keseragaman suatu barisan titik-titik dalam D = [0,1] d . Discrepancy menyebabkan galat dalam metode quasi-Monte Carlo. Discrepancy suatu barisan secara khusus yang digunakan untuk metode quasi-Monte Carlo dibatasi oleh sebuah waktu konstan
35
(log N )d N
Pers. (2.54)
Sebagai perbandingan, dengan peluang bernilai satu, discrepancy yang diharapkan dari barisan acak seragam (seperti yang digunakan dalam metode Monte Carlo) memiliki ordo konvergensi log log N 2N
Pers. (2.55)
menurut hukum dari algoritma teriterasi. Hal ini memperlihatkan bahwa keakuratan dari metode quasi-Monte Carlo meningkat lebih cepat dibandingkan dengan metode Monte Carlo. Walau bagaimanapun, menurut Morokoff and Caflisch (1995, p218-230), keuntungan dari quasi-Monte Carlo masih di bawah harapan dibandingkan dengan teorinya. Namun, masih menurut penelitian Morokoff and Caflisch, metode quasi-Monte Carlo dapat menghasilkan hasil yang lebih akurat dibandingkan dengan metode Monte Carlo dengan jumlah titik yang sama. Morokoff and Caflisch menambahkan bahwa keuntungan dari quasi-Monte Carlo jauh lebih baik jika integrannya halus (smooth) dan jumlah dimensi d dari integral adalah kecil. Untuk memperjelas, katakanlah terdapat sebuah fungsi D(x) dari himpunan titik yang bertambah dengan ketidakseragamannya. D(x) = 0 terjadi jika himpunan titik tersebut seragam sempurna dalam semua pernyataan, situasi yang ideal yang tidak pernah dapat diperoleh untuk himpunan titik apapun. Metode quasi-Monte Carlo didasarkan pada penggunaan himpunan titik x yang mana D(x) memiliki beberapa nilai s
36 (sangat banyak) yang lebih kecil dari s , nilai yang mungkin diharapkan untuk titik yang bebas, terdistribusi identik, dan seragam yang sebenarnya. Bagaimana
menggunakan
himpunan
titik
quasi-random
setelah
kita
memperolehnya? Hal yang pasti adalah dengan terlebih dahulu menentukan di kelompok manakah himpunan quasi-random dari titik x menjadi anggota yang khas. Distribusi dengan banyak titik dari himpunan titik PN yang demikian tidak lagi merupakan kesatuan yang sederhana, yang akan berarti kebebasan dari titik-titik tersebut di dalam himpunan titiknya. Karena itu, kita tulis distribusi dengan banyak titik sebagai r r r r r r 1 PN (s; x1 , x2 , ... , xN ) = 1 − FN (s; x1 , x2 , ... , xN ) N
Pers. (2.56)
di mana kita telah mengantisipasi sebuah faktor 1/N sebelum korelasi dengan banyak titik FN. FN(s; ...) harus simetris total; selain itu kita harus memiliki r r r r r r r r Fk (s; x1 , x2 , ... , xK ) = ∫ Fk +1 (s; x1 , x2 , ... , xk , xk +1 ) d d xk +1
Pers. (2.57)
C
Untuk memenuhi syarat minimum bahwa integral quasi-Monte Carlo haruslah tidak bias, kita harus memiliki r P1 (s; x1 ) =1
Pers. (2.58)
sehingga r r
∫ F ( s; x , x ) d 2
1
2
d
r x2 = 0
Pers. (2.59)
C
Hal ini membangun sifat-sifat dari kelompok himpunan titik x yang mana
37 seharusnya pendugaan yang dilakukan didasarkan. Kita akan mengindikasikan sifat alamiah quasi-Monte Carlo mengenai penduga dengan menggunakan (q). Berikut adalah penduga pertama dari integral
E1(q ) =
1 N
N
∑f i =1
Pers. (2.60)
j
Penjumlahan akan berjalan dari 1 hingga N. Berdasarkan rataan (q) terhadap kelompok quasi-random yang telah didiskusikan di atas, kita kemudian memiliki E1(q )
(q )
r r r = ∫ f ( x ) P1 (s; x ) d d x = J1
Pers. (2.61)
C
sebagaimana sebelumnya: berdasarkan fakta bahwa distribusi satu titik adalah seragam, pendugaan quasi-Monte Carlo memang tidak bias sebagaimana yang dimiliki oleh metode Monte Carlo. Perbedaan yang jelas antara kedua metode tersebut tampak pada pendugaan error ordo satu. Kita definisikan
α (xi , x j )=1 + F2 (s; xi , x j ) r r
r r
Pers. (2.62)
kemudian kita memiliki
σ (E1(q ) )(q ) =
(
)
1 ⎛ 1 ⎞ I 2 − ∫ f 1 f 2α 12 + O⎜ 2 ⎟ N ⎝N ⎠
Pers. (2.63)
r r r r r r f α12 = ∫ f ( x1 ) f ( x2 )α ( x1 , x2 )d d x1d d x2
Pers. (2.64)
2
di mana
∫f
1 2
dan lainnya.
C
38 Keuntungan dari metode quasi-Monte Carlo sekarang menjadi jelas: jika kita r r dapat memastikan bahwa α12 > 1, yaitu saat x1 dan x2 “dekat” pada suatu pengertian,
maka error quasi-Monte Carlo akan lebih kecil dibandingkan metode Monte Carlo. Sebuah himpunan titik quasi-Monte Carlo yang baik adalah himpunan titik yang mana setiap titiknya saling menolak untuk beberapa perluasan. Pendugaan error ordo pertama secara sederhana
E2(q ) =
1 N2
N
∑ i =1
1 N ∑ f i f jα ij N 3 i =1
fi 2 −
Pers. (2.65)
Mudah untuk menunjukkan bahwa sesungguhnya E2(q )
(q )
( )( ) + O( N
= σ E1(q )
2
q
−2
)
Pers. (2.66)
ragam dari penduga E2(q ) dapat dievaluasi menjadi
σ (E2(q ) ) = 2
1 N3
(∫ f
4
i
− 4∫ f i 3 f jα ij − ∫ f i 2 f j2α ij + 4∫ f i 2 f k f lα ikα kl
) ( )
+ 4 ∫ f i 2 f k f lα ikα il − 4∫ f i f j f k f lα ijα jkα kl + O N −4 yang mana penduga yang bersesuaian adalah
E4(q ) =
N N 1 ⎛ 3 N 4 2 3 2 ⎜ N f 4 N f f α N f i 2 f j2α ij − − ∑ ∑ ∑ i i j ij 7 ⎜ N ⎝ i =1 i =1 i =1
N
N
i =1
i =1
+ 4 N ∑ f i 2 f k f lα ikα kl + 4 N ∑ f i 2 f k f lα ikα il
Pers. (2.67)
39 N ⎞ − 4∑ f i f j f k f lα jkα kl ⎟⎟ i =1 ⎠
Pers. (2.68)
2.6.1. Barisan Sobol (Sobol Sequence)
Barisan Sobol (Sobol sequence) dihasilkan dari sebuah himpunan bilangan pecahan biner khusus sepanjang w bit, vij dengan i = 1, 2, ... , w dan j = 1, 2, ... , d. Bilangan vij disebut direction numbers.
Untuk menghasilkan direction numbers untuk dimensi j, digunakan sebuah polinom primitif (tidak dapat disederhanakan lagi) terhadap bidang F2 dengan elemen {0, 1}. Polinom primitif tersebut dalam dimensi j akan menjadi p j ( x ) = x q + a1 x q −1 + K + aq −1 x + 1
Pers. (2.69)
Direction numbers tersebut dalam dimensi j dihasilkan menggunakan
hubungan q berulang berikut
(
vij = a1vij−1 ⊕ a2 vij− 2 ⊕ K ⊕ aq −1vij− q +1 ⊕ vij− q ⊕ vij− q 2 q
)
Pers. (2.70)
di mana i > q. ⊕ menandakan operasi XOR. Bilangan v1j ⋅ 2 w , v2j ⋅ 2 w , ... , vqj ⋅ 2 w masing-masing dapat berupa bilangan integer ganjil lebih kecil dari 2, 22, ... , dan 2q. Barisan Sobol xnj ( n = ∑i = 0 bi 2i , bi ∈ {0, 1}) dalam dimensi j dihasilkan oleh w
xnj = b1v1j ⊕ b2 v2j ⊕ K ⊕ bw vwj
Pers. (2.71)
Sebuah polinom primitif yang berbeda harus digunakan untuk
40 menghasilkan barisan Sobol di masing-masing dimensi.
