Задержка значений параметров / скорости притока и оттока в deSolve

#r #differential-equations

#r #дифференциальные уравнения

Вопрос:

Я использую R пакет deSolve для решения систем обыкновенных дифференциальных уравнений. В литературе по «системной динамике» задержки скоростей притока и оттока могут быть смоделированы с использованием среднего времени задержки. Например, скорость изменения запаса Y в момент времени t может быть:

 dy(t)/dt = inflow(t) - ( outflow(t) / D )
  

где время задержки D равно, например, 4 временным шагам. Предполагается, что задержка равна среднему времени задержки.

Однако другим способом моделирования задержек было бы предположить более дискретный случай события, когда отток равен количеству, поступившему на склад D в единицах времени ранее, таким образом:

 dy(t)/dt = inflow(t) - inflow(t - D)
  

В deSolve мы можем использовать функции lagvalue и lagderiv с dede функцией solver для задания дифференциальных уравнений задержки, которые используют запаздывающие значения переменных состояния, но, похоже, я не могу найти способ попросить deSolve использовать запаздывающие значения скоростей притока / оттока.

Для примера возьмем простую модель:

 m<- function(t,y,p){
    with(as.list(c(y,p)),{
        inflow <- 100
        outflow <- y*.5
        dy <- inflow - outflow
        return(list(c(dy), inflow=inflow, outflow=outflow))
    })}

fit <- ode(func=m, y=c(100),t=seq(0,10,1),p=c(), method="euler") 

   time        1 inflow  outflow
1     0 100.0000    100 50.00000
2     1 150.0000    100 75.00000
3     2 175.0000    100 87.50000
4     3 187.5000    100 93.75000
5     4 193.7500    100 96.87500
6     5 196.8750    100 98.43750
7     6 198.4375    100 99.21875
8     7 199.2188    100 99.60938
9     8 199.6094    100 99.80469
10    9 199.8047    100 99.90234
11   10 199.9023    100 99.95117

  

Используя dede , я могу сделать отток запаздывающим значением переменной состояния на D = 2 предыдущих временных шагах:

 m2<- function(t,y,p){
    with(as.list(c(y,p)),{
        inflow <- 100
        if(t < D) outflow <- y*.5
        if(t >= D) outflow <- lagvalue(t-D,1)*.5
        dy <- inflow - outflow
        return(list(c(dy), inflow=inflow, outflow=outflow))
    })}

fit2 <- dede(func=m, y=c(100),t=seq(0,10,1),p=c(D=2))

   time        1 inflow  outflow
1     0 100.0000    100 50.00000
2     1 139.3469    100 69.67344
3     2 163.2120    100 81.60602
4     3 177.6870    100 88.84349
5     4 186.4665    100 93.23323
6     5 191.7915    100 95.89575
7     6 195.0213    100 97.51064
8     7 196.9803    100 98.49013
9     8 198.1684    100 99.08422
10    9 198.8891    100 99.44456
11   10 199.3262    100 99.66312
  

Но теперь представьте, что я хочу, чтобы отток фактически соответствовал предыдущим временным шагам притока D=2 . Я хочу что-то вроде:

 **** Code will not run ****
m3<- function(t,y,p){
    with(as.list(c(y,p)),{
        inflow <- 100
        if(t < D) outflow <- 0
        if(t >= D) outflow <- lagvalue(t-D,inflow)
        dy <- inflow - outflow
        return(list(c(dy), inflow=inflow, outflow=outflow))
    })}
...
  

Насколько я могу видеть, deSolve этого не допускает. Есть ли простой способ разрешить это?

Причина, по которой я заинтересован в смешивании модели непрерывного и дискретного типов событий, заключается в моделировании цепочек поставок, где средняя временная задержка может быть неточной для определенных продуктов.

Комментарии:

1. Я могу придумать два сценария. 1) Приток не зависит от y , но является функцией времени. В этом случае приток может быть определен как функция вне вашего ode. Если отток является функцией от притока, то вы просто просматриваете приток с соответствующим запаздыванием и вычисляете отток, используя это значение. 2) Приток является функцией y . В этом случае, чтобы рассчитать отток, вы должны выполнить поиск y с соответствующим запаздыванием, рассчитать приток на основе этого значения, а затем рассчитать отток из этого притока. Или я упустил суть?

2. Не могли бы вы привести пример, имеющий отношение к моему коду? Как бы мне «просто посмотреть приток с соответствующей задержкой», используя пакет deSolve?

3. Определите функцию для притока, скажем, calc_inflow <- function(t)10 0.1 * t , затем вызовите ее из вашей модели inflow <- calc_inflow(t) . Если отток равен притоку на две единицы времени раньше, то outflow <- calc_inflow(t-2) .

4. Спасибо! По большей части это работает.