RunModel_Lag returns 2 values for a one time step run

The example below tries to run a daily SD model for the time step 1990-01-01 (adapted from tests\testthat\test-RunModel_Lag.R):

library(airGR)
data(L0123001)
BasinAreas <- c(BasinInfo$BasinArea, BasinInfo$BasinArea)
# Qupstream = sinusoid synchronised on hydrological year from 0 mm to mean value of Qobs
Qupstream <- floor((sin((seq_along(BasinObs$Qmm)/365*2*3.14))+1) * mean(BasinObs$Qmm, na.rm = TRUE))
InputsModel <- CreateInputsModel(
  FUN_MOD = RunModel_GR4J,
  DatesR = BasinObs$DatesR,
  Precip = BasinObs$P,
  PotEvap = BasinObs$E,
  Qupstream = matrix(Qupstream, ncol = 1),
  LengthHydro = 1000,
  BasinAreas = BasinAreas
)
Ind_Run <- seq(which(format(BasinObs$DatesR, format = "%Y-%m-%d") == "1990-01-01"),
               which(format(BasinObs$DatesR, format = "%Y-%m-%d") == "1990-01-01"))
RunOptions <- suppressWarnings(CreateRunOptions(FUN_MOD = RunModel_GR4J,
                                                InputsModel = InputsModel,
                                                IndPeriod_Run = Ind_Run))
Param <- c(1., 257.237556, 1.012237, 88.234673, 2.207958)
OutputsModel <- RunModel(InputsModel = InputsModel, RunOptions = RunOptions, Param = Param, FUN_MOD = RunModel_GR4J)
OutputsModel$Qsim

The output of this code is:

> OutputsModel$Qsim
[1] 1.709998 1.715785

The error comes from the lines in R\RunModel_Lag.R:

OutputsModel$Qsim <- OutputsModel$Qsim +
      c(IniStates[[upstream_basin]][-length(IniStates[[upstream_basin]])],
        Qupstream[1:(LengthTs - floor(PT[upstream_basin]))]) *
      HUTRANS[1, upstream_basin] +
      c(IniStates[[upstream_basin]],
        Qupstream[1:(LengthTs - floor(PT[upstream_basin]) - 1)]) *
      HUTRANS[2, upstream_basin]

Especially from the expression (LengthTs - floor(PT[upstream_basin]) - 1) which is equal to zero when LengthTs = 1. Initial states and Qupstream should be merge only if we need to address them.

I tried locally to solve the problem with this solution:

Qupstream1 <- c(IniStates[[upstream_basin]][-length(IniStates[[upstream_basin]])], Qupstream[1:(LengthTs - floor(PT[upstream_basin]))])
Qupstream2 <- IniStates[[upstream_basin]]
if(LengthTs - floor(PT[upstream_basin]) - 1 > 0) Qupstream2 <- c(Qupstream2, Qupstream[1:(LengthTs - floor(PT[upstream_basin]) - 1)])
OutputsModel$Qsim <- OutputsModel$Qsim +
                     Qupstream1 * HUTRANS[1, upstream_basin] +
                     Qupstream2 * HUTRANS[2, upstream_basin]