[R-sig-dyn-mod] Which parameter identifies the 'time' to reach steady state with runsteady?

[Vince Pile]

Hello,
New to working with deSolve and rootSolve and ran the example under
runsteady and tried to determine the time to reach steady state but despite
all my attempts to manipulate the 'istate' or 'rstate' values to match the
plot output could not work it out.
Thanks in advance for your help!
Vince
## =============================================
## A simple sediment biogeochemical model
## =============================================

model<-function(t, y, pars) {

  with (as.list(c(y, pars)),{

    Min       = r*OM
    oxicmin   = Min*(O2/(O2+ks))
    anoxicmin = Min*(1-O2/(O2+ks))* SO4/(SO4+ks2)

    dOM  = Flux - oxicmin - anoxicmin
    dO2  = -oxicmin      -2*rox*HS*(O2/(O2+ks)) + D*(BO2-O2)
    dSO4 = -0.5*anoxicmin  +rox*HS*(O2/(O2+ks)) + D*(BSO4-SO4)
    dHS  = 0.5*anoxicmin   -rox*HS*(O2/(O2+ks)) + D*(BHS-HS)

    list(c(dOM, dO2, dSO4, dHS), SumS = SO4+HS)
  })
}

# parameter values
pars <- c(D = 1, Flux = 100, r = 0.1, rox = 1,
          ks = 1, ks2 = 1, BO2 = 100, BSO4 = 10000, BHS = 0)
# initial conditions
y <- c(OM = 1, O2 = 1, SO4 = 1, HS = 1)

# time
times = seq(0, 40, 0.05)

# dynamic soln
dynOut <- lsoda(y = y, func = model, time = times, parms = pars, verbose =
TRUE)

# steady state soln
ssOut <- runsteady(y = y, func = model, parms = pars, times = c(0, 1000),
verbose=TRUE)

# plot
head(dynOut)
matplot(dynOut[,1], dynOut[,c("OM", "O2", "HS")], type="l", lty=1:4,
col=1:4, lwd=3)
legend("topleft", legend=c("OM", "O2", "HS"), lty=1:4, col=1:4, lwd=3)
abline(h=seq(0,900,100), v=seq(0,40,5), lty=2, col="grey")

	[[alternative HTML version deleted]]