Rでオイラー法を解く

#one reaction

step_num <- 10000

dt <- 0.001

S_initial <- 100

P_initial <- 0

k1 <- 1

S <- numeric(step_num+1)

P <- numeric(step_num+1)

t <- numeric(step_num+1)

S[1] <- S_initial

P[1] <- P_initial

t[1] <- 0

for (i in 1:step_num){

  S[i+1] <- S[i] - k1 * S[i] *dt

  P[i+1] <- P[i] + k1 * S[i] *dt

  t[i+1] <- t[i] + dt

}

plot(t, S, ylim =c(0,120))

points(t, P, col="red")