55 lines
1.5 KiB
R
55 lines
1.5 KiB
R
# Project 3 for the University of Tulsa's CS-7863 Sci-Stat Course
|
|
# Numerical Ordinary Differential Equations
|
|
# Professor: Dr. McKinney, Spring 2023
|
|
# Noah L. Schrick - 1492657
|
|
|
|
## 1. Approximate deriv. of sin(x) at x=pi from h=10^-1 -> 10^-20
|
|
forward.approx <- function(func, x, h){
|
|
approx <- (func(x+h) - func(x))/h
|
|
}
|
|
|
|
forward.approx.table <- matrix(nrow = 0, ncol = 3)
|
|
colnames(forward.approx.table) <- c("h", "approximation", "error")
|
|
for(h in 10^(seq(-1,-20,-1))){
|
|
approx <- forward.approx(sin, pi, h)
|
|
error <- abs(cos(pi)-approx)
|
|
forward.approx.table <- rbind(forward.approx.table, c(h, approx, error))
|
|
}
|
|
|
|
plot(abs(log(forward.approx.table[,"h"],10)), forward.approx.table[,"error"],
|
|
xlab="h [10^-x]", ylab="error (logscale)", type="o", log="y")
|
|
|
|
# Repeat with central difference approx
|
|
|
|
## 2.
|
|
# a) Runge-Kutta
|
|
|
|
# b) Numerically solve the decay ode and compare to Euler and RK error
|
|
|
|
## 3. Use library function ode45 to solve the decay numerically
|
|
|
|
## 4. Solve the predator-prey model numerically
|
|
|
|
# a) k1=0.01, k2=0.1, k3=0.001, k4=0.05, prey(0)=50, pred(0)=15. t=0 -> 200
|
|
|
|
# plot
|
|
|
|
# b) Use Euler and RK to solve with h=10
|
|
|
|
# plot comparing Prey solutions to ode45
|
|
|
|
# c) Use k3=0.02
|
|
|
|
## 5. Solve SIR model numerically from t=0 -> 20
|
|
|
|
# a) a=0.5, b=1, S(0)=0.9, I(0)=0.1, R(0)=0
|
|
|
|
# plot
|
|
|
|
# b) a=3
|
|
|
|
## 6. Decomp of dinitrogen pentoxygen into nitrogen dioxide and molecular oxygen
|
|
|
|
# a) k1=1.0, k2=0.5, k3=0.2,k4=1.5, [N2O5]o=1, all other IC's=0, t=0 -> 10
|
|
|
|
# b) Increase k4 to make the intermediate species -> 0 |