Modelling an SIRD model to data with Death timestep totals

Question

hmod = (assign_df)

sir_model = function (t, y, parms)
{
  # label state variables from within the y vector of initial values
  S = y[1]        # susceptibles
  I = y[2]        # infectious
  R = y[3]        # recovered
  D = y[4]        # dead
  N = S + I + R + D
  
  # Pull parameter values from the parms vector
  beta = parms["beta"]
  gamma = parms["gamma"]
  alpha = parms["alpha"]
  
  # Define equations (dS etc. are simplified notation for dS/dt)
  dS = -(beta*S*I)/N
  dI = (beta*S*I)/N - gamma*I
  dR = alpha*gamma*I
  dD = (1-alpha)*gamma*I
  
  # return list of gradients   
  res = c(dS, dI, dR, dD)
  list(res)
}

par_init = c(beta=1, gamma=1, alpha =0.5)
init_sav = c(S=1079, I=1, R=0, D=0)
times_sav = hmod$time

n_obs <- length(hmod$D_obs)

opt_fun_sav = function(par){
  
  ## The following lines prevent 'optim_nm' from wandering into negative parameter space
  if (par[1]<0){
    par[1] = -par[1]
  } 
  if(par[2]<0){
    par[2] = -par[2]
  }
  ## Run the SIR model with the current parameter values 
  output_sav_try = ode(
    y=init_sav, 
    times = times_sav, 
    func = sir_model, 
parms= c(
      beta = par_init[1], 
      gamma = par_init[2],
      alpha = par_init[3]
    )
  )
  n_obs = nrow(hmod)
  
  output_savtry_df = as.data.frame(output_sav_try) # coerce to a data.frame object

  
  ## Create a data frame with observed and model prediction data for Savoonga
  sav_rmse_df = data.frame(
    time = hmod$time,
    obs = hmod$D_obs,
    # the following lines find the cumulative I prediction for each time in the mummps data
    pred = sapply(1:n_obs, function(i){output_savtry_df$D[
      output_savtry_df$time == hmod$time[i]
    ]})
  )
  # Calculate the RMSE for this prediction
  rmse_try = sqrt((1/n_obs)*(sum((sav_rmse_df$pred - sav_rmse_df$obs)^2)))  
  
  return(rmse_try)
  
}

optnm_result=optim_nm(opt_fun_sav,start=c(0.6,0.6),k=2, trace=TRUE)

Error in while (tol < End) { : missing value where TRUE/FALSE needed

I can't get around this error, I'm also not sure if other parts of my code is wrong The aim of this code is to optomise the prediction of ra low RMSE when compared with death counts in the dataset hmod

I've tried alot of things, removing pieces which seems to be causing issues, checking where False or true values could be, none of it has worked.

Use the code below to replicate the dataset of real death totals

time <- 1:40
D_obs <- round(seq(1, 200, length.out = 40))
hmod <- data.frame(time, D_obs)
print(hmod)