# Nama : Atmariani Ardiningrum
# NRP : 5003251197
# Kelas : Komputasi Statistika D
# Trapezoidal Rule
trapezoidal <- function(f, a, b, n) {
h <- (b - a) / n
x <- seq(a, b, by = h)
y <- f(x)
integral <- (h / 2) * (y[1] + 2 * sum(y[2:n]) + y[n + 1])
return(integral)
}
# Simpson's 1/3 Rule
simpsons <- function(f, a, b, n) {
if (n %% 2 != 0) stop("n harus genap untuk aturan Simpson's 1/3")
h <- (b - a) / n
x <- seq(a, b, by = h)
y <- f(x)
sum_odd <- sum(y[seq(2, n, by = 2)])
sum_even <- sum(y[seq(3, n-1, by = 2)])
integral <- (h / 3) * (y[1] + 4 * sum_odd + 2 * sum_even + y[n + 1])
return(integral)
}
# --- A: NORMAL DISTRIBUTION ---
f_normal <- function(x) {
(1 / sqrt(2 * pi)) * exp(-x^2 / 2)
}
a_norm <- 0
b_norm <- 1
n_norm <- 4
res_simpson_A <- simpsons(f_normal, a_norm, b_norm, n_norm)
exact_norm <- pnorm(1) - pnorm(0)
cat("--- Exercise A (Normal Distribution) ---\n")
## --- Exercise A (Normal Distribution) ---
cat("Simpson's Rule Result:", res_simpson_A, "\n")
## Simpson's Rule Result: 0.3413555
cat("R pnorm() Result :", exact_norm, "\n")
## R pnorm() Result : 0.3413447
cat("Error :", abs(exact_norm - res_simpson_A), "\n\n")
## Error : 1.074179e-05
# --- B: EXPONENTIAL DISTRIBUTION ---
f_exp <- function(x) {
exp(-x)
}
a_exp <- 0
b_exp <- 2
n_exp <- 10
res_trap_B <- trapezoidal(f_exp, a_exp, b_exp, n_exp)
res_simpson_B <- simpsons(f_exp, a_exp, b_exp, n_exp)
exact_exp <- pexp(2, rate = 1) - pexp(0, rate = 1)
cat("--- Exercise B (Exponential Distribution) ---\n")
## --- Exercise B (Exponential Distribution) ---
cat("Trapezoidal Result :", res_trap_B, "\n")
## Trapezoidal Result : 0.867545
cat("Simpson's Result :", res_simpson_B, "\n")
## Simpson's Result : 0.8646724
cat("R pexp() Result :", exact_exp, "\n")
## R pexp() Result : 0.8646647