Tugas Minggu 6_Laskar Brahimanov Riadi_5003251048

Exercise A

f <- function(x){
  1/sqrt(2*pi) * exp(-x^2/2)
}

a <- 0
b <- 1
n <- 4

if(n %% 2 != 0){
  stop("Simpson's rule requires n to be even")
}

h <- (b-a)/n
cat("Step size h =", h, "\n\n")

## Step size h = 0.25

x <- seq(a,b,length.out = n+1)
y <- f(x)

print(data.frame(i=0,x=x,y=y))

##   i    x         y
## 1 0 0.00 0.3989423
## 2 0 0.25 0.3866681
## 3 0 0.50 0.3520653
## 4 0 0.75 0.3011374
## 5 0 1.00 0.2419707

total <- 0
for(i in 1:(n+1)){
  idx <- i-1
  if(idx == 0 || idx == n){
    weight <- 1
  }
  else if(idx %% 2 == 1){
    weight <- 4
  }
  else{
    weight <- 2
  }
  total <- total + weight * y[i]
}
sn <- h/3 * total

cat("Simpson's Approximation =", sn, "\n")

## Simpson's Approximation = 0.3413555

exact <- pnorm(1)-pnorm(0)
cat("Exact value =", exact,"\n")

## Exact value = 0.3413447

cat("Error = ", sn - exact,"\n")

## Error =  1.074179e-05

Exercise B

f <- function(x){
  return(exp(-x))
}

a <- 0
b <- 2
n <- 10
h <- (b - a) / n
cat("Step size h =", h, "\n\n")

## Step size h = 0.2

x <- seq(a, b, length.out = n + 1)

# Trapezoidal Rule
total <- 0
for(i in 1:(n + 1)){
  if(i == 1 || i == (n + 1)){
    total <- total + f(x[i])
  } else {
    total <- total + 2 * f(x[i])
  }
}
tn <- (h / 2) * total

# Simpson's Rule
if(n %% 2 != 0){
  stop("Simpson's rule requires n to be even")
}

totals <- 0
for(i in 1:(n + 1)){
  idx <- i - 1
  if(idx == 0 || idx == n){
    weight <- 1
  } else if(idx %% 2 == 1){
    weight <- 4
  } else {
    weight <- 2
  }
  totals <- totals + weight * f(x[i])  
}
sn <- h / 3 * totals                    

exact <- 1 - exp(-2)
cat("Trapezoidal Approximation =", tn, "\n")

## Trapezoidal Approximation = 0.867545

cat("Simpson's Approximation   =", sn, "\n")

## Simpson's Approximation   = 0.8646724

cat("Exact value (1 - e^-2)    =", exact, "\n")

## Exact value (1 - e^-2)    = 0.8646647

cat("Error Trapezoidal         =", abs(tn - exact), "\n")

## Error Trapezoidal         = 0.002880296

cat("Error Simpson's           =", abs(sn - exact), "\n")

## Error Simpson's           = 7.649462e-06