Email : ronaldokefas@gmail.com
RPubs : https://rpubs.com/kefasronaldo/
Department : Business Statistics
Address : ARA Center, Matana University Tower
Jl. CBD Barat Kav, RT.1, Curug Sangereng, Kelapa Dua, Tangerang, Banten 15810.
Hitung integral fungsi \(f(x)=sin^{2}(x)\) pada domain \(x\in [0, \pi]\)
Disini perhitungan akan menggunakan metode rieman dengan m = 4 \[ \begin{align} \int_0^\pi \sin^{2}(x)dx&\approx \sum^4_{i=1}f(i\frac{|\pi|}{4}-\frac{|\pi|}{8}) ~\frac{|\pi|}{4}\\ &=\left(f(\frac{\pi}{8})+f(\frac{3\pi}{8})+f(\frac{5\pi}{8})+f(\frac{7\pi}{8})\right)\frac{|\pi|}{4} \end{align} \] \[ \begin{align} f(\frac{\pi}{8})&=\sin^{2}(\frac{\pi}{8})\\ &=\frac{2-\sqrt{2}}{4}\\ f(\frac{3\pi}{8})&=\sin^{2}(\frac{3\pi}{8})\\ &=\frac{2+\sqrt{2}}{4}\\ f(\frac{5\pi}{8})&=\sin^{2}(\frac{5\pi}{8})\\ &=\frac{2+\sqrt{2}}{4}\\ f(\frac{7\pi}{8})&=\sin^{2}(\frac{7\pi}{8})\\ &=\frac{2-\sqrt{2}}{4} \end{align} \] \[ \begin{align} \int_0^\pi \sin^{2}(x)dx&\approx\left(\frac{2-\sqrt{2}}{4}+\frac{2+\sqrt{2}}{4}+\frac{2+\sqrt{2}}{4}+\frac{2-\sqrt{2}}{4}\right)\frac{|\pi|}{4}\\ &=1.570796 \end{align} \]
Tuliskan fungsi R yang dapat melakukan integrasi Riemann dengan aturan titik kiri!
Algoritma untuk Fungsi integral Riemann adalah sebagai berikut:
riemann <- function(f, a, b, m = 100){
n_width <- (b-a)/m
x <- seq(a, b-n_width, length.out = m) + n_width/2
y <- f(x)
return(sum(y)*abs(b-a)/m)
}Sekarang saya akan menguji fungsi ini menggunakan soal dari nomor 1
riemann(function(x) (sin(x))^(2), a=0, b=pi, m=4)## [1] 1.570796
Hasil perhitungan dengan fungsi dan secara manual adalah sama yang membuktikan bahwa perhitungan manual dan fungsi sudah tepat. hehe
Buatlah sebuah fungsi R yang dapat melakukan integrasi adaptif menggunakan metode Simpson 1/3!
Algoritma untuk Fungsi integral Simpson 1/3 adalah sebagai berikut:
simpson <- function(fun, a, b, n=100) {
# numerical integral using Simpson's rule
# assume a < b and n is an even positive integer
h <- (b-a)/n
x <- seq(a, b, by=h)
if (n == 2) {
s <- fun(x[1]) + 4*fun(x[2]) +fun(x[3])
} else {
s <- fun(x[1]) + fun(x[n+1]) +
2*sum(fun(x[seq(2,n,by=2)])) +
4 *sum(fun(x[seq(3,n-1, by=2)]))
}
s <- s*h/3
return(s)
}
simpson_adapt <- function(fun, a, b, tol=1e-8) {
# numerical integration using adaptive quadrature
quadrature_internal <- function(S.old, fun, a, m, b, tol, level) {
level.max <- 100
if (level > level.max) {
cat ("recursion limit reached: singularity likely\n")
return (NULL)
}
S.left <- simpson(fun, a, m, n=2)
S.right <- simpson(fun, m, b, n=2)
S.new <- S.left + S.right
if (abs(S.new-S.old) > tol) {
S.left <- quadrature_internal(S.left, fun, a, (a+m)/2, m, tol/2, level+1)
S.right <- quadrature_internal(S.right, fun, m, (m+b)/2, b, tol/2, level+1)
S.new <- S.left + S.right
}
return(S.new)
}
level = 1
S.old <- (b-a) * (fun(a) + fun(b))/2
S.new <- quadrature_internal(S.old, fun, a, (a+b)/2, b, tol, level+1)
return(S.new)
}Sekarang saya akan menguji fungsi ini menggunakan soal dari nomor 1
simpson_adapt(function(x) (sin(x))^(2), a=0, b=pi)## [1] 1.570796
Hasil perhitungan dengan fungsi Adaptive Simpson 1/3 menghasilkan nilai yang sama sehingga fungsi Adaptive simpson 1/3 ini sudah bekerja dengan baik.
Kerjakan kembali soal 3 dengan menggunakan metode Simpson 3/8!
Algoritma untuk Fungsi integral Simpson 3/8 adalah sebagai berikut:
simpson38 <- function(f, a, b, m=90){
h <- (b-a)/m # jarak selang
x <- a # awal selang
I <- f(a)+f(b)
sigma <- 0
if(m%%3 != 0){
stop("jumlah panel harus kelipatan 3")
}else{
for(i in 1:(m-1)){
x <- x+h
if(i%%3==0){
sigma <- sigma + 2*f(x)
}else{
sigma <- sigma + 3*f(x)
}
}
}
return((3*h/8)*(I+sigma))
}Sekarang saya akan menguji fungsi ini menggunakan soal dari nomor 1
simpson38(function(x) (sin(x))^(2), a=-pi, b=0, m=6)## [1] 1.570796
Hasil perhitungan dengan fungsi Simpson 3/8 menghasilkan nilai yang sama sehingga fungsi simpson ini sudah bekerja dengan baik.
Sekarang kita akan melakukan perhitungan secara manual menggunakan m = 6
\[ \begin{align} h &=\frac{b-a}{m}=\frac{\pi}{6} \\ f_0 &=f(0)=\sin^2 (0)=0 \\ f_1 &= f\left(\frac{\pi}{6} \right) = \sin^2 \left(\frac{\pi}{6} \right)=0.25 \\ f_2 &= f\left(\frac{\pi}{3} \right) = \sin^2 \left(\frac{\pi}{3} \right)=0.75 \\ f_3 &= f\left(\frac{\pi}{2} \right) = \sin^2 \left(\frac{\pi}{2} \right)=1 \\ f_4 &= f\left(\frac{2\pi}{3} \right) = \sin^2 \left(\frac{2\pi}{3} \right)=0.75 \\ f_5 &= f\left(\frac{5\pi}{6} \right) = \sin^2 \left(\frac{5\pi}{6} \right)=0.25 \\ f_6 &= f\left(\pi\right) = \sin^2 \left(\pi\right)=0 \end{align} \]
| n | 0 | 1 | 2 | 3 | 4 | 5 | 6 |
|---|---|---|---|---|---|---|---|
| x | 0 | \(\frac{\pi}{6}\) | \(\frac{\pi}{3}\) | \(\frac{\pi}{2}\) | \(\frac{2\pi}{3}\) | \(\frac{5\pi}{6}\) | \(\pi\) |
| y | 0 | 0.25 | 0.75 | 1 | 0.75 | 0.25 | 0 |
\[ \begin{align} h&=\frac{b-a}{m}=\frac{\pi}{6} \\ \int_0^\pi \sin^{2}(x)~dx&\approx \frac{3}{8}\times \frac{\pi}{6}\left(f_0+3\sum^{5}_{i=1,i\ne3, 6, 9,\cdots} f_i+2\sum^{3}_{i=3, 6, 9,\cdots} f_i+f_6\right)\\ &=\frac{\pi}{16}\left(f_0+3(f_1+f_2+f_4+f_5)+2f_3+f_6\right)\\ &=\frac{\pi}{16}\left(0+3(0.25+0.75+0.75+0.25)+2*1+0\right)\\ &=\frac{\pi}{16}\left(3(2)+2\right) = \frac{\pi}{2}=1.570796 \end{align} \]
Hasil perhitungan manual sesuai dengan fungsi, mengartikan perhitungan sudah tepat.
Fungsi monte_int() hanya mampu melakukan integrasi pada domain positif. Buatlah algoritma baru sehingga metode ini dapat melakukan integrasi pada domain negatif!
monte_int <- function(f, a, b, m=1e6){
x <- runif(m, min=a, max=b)
return((b-a)*sum(f(x))/m)
}Kita akan menguji apakah fungsi ini bekerja pada domain negatif atau tidak dengan menerapkannya pada soal dari nomor 1 yang berada dalam domain negatif \([-\pi, 0]\). Jadi soal yang akan kita gunakan untuk pengujian adalah:
\[ \begin{align} \int^0_{-\pi} \sin^{2}(x)~dx=\cdots \end{align} \]
monte_int(function(x) (sin(x))^2, a=-pi, b=0)## [1] 1.571055
Sekarang kita akan menguji hasil perhitungan fungsi ini untuk memastikan apakah fungsi ini benar bekerja di domain negatif atau tidak. Disini akan diuji dengan metode Simpson 3/8
simpson38(function(x) (sin(x))^2, a=-pi, b=0)## [1] 1.570796
Hasil yang dihasilkan memang tidak sama identik, tetapi telah mendekati, karena perhitungan Monte Carlo memang memiliki kualitas yang rendah dibanding perhitungan yang lainnya. Karena hasil Monte telah mendekati Simpson 3/8, maka perhitungan di domain negatif bekerja.
Metode Monte Carlo memiliki kualitas yang lebih rendah dibandingkan dengan metode lainnya karena pendekatan metode ini adalah pendekatan bilangan acak yang mana mengakibatkan kemungkinan terdapat error dalam menjumpai nilai sejati cukup besar. Sekalipun begitu metode Monte Carlo memiliki keuntungan karena dapat menangani integrasi berganda (bivariat).