Построение графика решения системы дифференциальных уравнений, меняющей параметр времени относительно R - PullRequest
Новая
       87

Построение графика решения системы дифференциальных уравнений, меняющей параметр времени относительно R

0 голосов
/ 07 апреля 2020

я делал этот код для построения решения модели Курамото, моя цель состояла в том, чтобы построить решение, меняющее коэффициент связи K, а затем собрать все графики вместе, чтобы увидеть, как эволюция связи и моя попытка были здесь: 

x11()
plot(times1, s1,type="l",xlim=c(0,180),col="white")
value <- seq(0, 3, by=0.1)
pe=1
for(numero in value){
  parameters <- c(N = 2, w1=0.4838171, w2=0.8842176, K=numero)
  state <- c(theta1=0, theta2=0)

  FN <- function (t, state, parameters) {
    with(as.list(c(state,parameters)),{
      dtheta1 <- 2*pi*w1 + (K/N)*(sin(theta1-theta1)+sin(theta2-theta1))
      dtheta2 <- 2*pi*w2 + (K/N)*(sin(theta1-theta2)+sin(theta2-theta2))


      list(c(dtheta1,dtheta2))
    })
  }
  times1 <- seq(1800*(value[pe]*10),1800*(value[1+pe]*10), by=1)*T
  if (pe<=29){
    pe=pe+1
  }
  library(deSolve)
  out <- ode(y = state,times = times1,func = FN,parms=parameters)
  head(out)


  s1 = (1/length(state))*(cos(out[,2])+cos(out[,3]))
  lines(times1, s1,type="l")
}

Но теперь я хочу иметь K, которое зависит от времени, и оно меняется, то есть взять K как некоторую функцию, например K = K (t), и построить график решения, которое в то же время, что время вектор меняет изменения параметра K и график, который меняется в зависимости от решения ... но я не знаю, как решить систему с параметром, зависящим от времени, я читал здесь некоторые вопросы, но ни один из них не помог, можете ли вы дать мне несколько помогите или как сделать параметр зависимым от этого ?, большое спасибо!

1 Ответ

0 голосов
/ 09 апреля 2020

Общий подход состоит в том, чтобы определить «форсирующую функцию», которая интерполируется с R approxfun(). Поскольку ваш код содержит несколько проблем, которые мешают воспроизводимости, я покажу это на другом примере. По сути, это модель роста логистики c с параметром несущей способности K, зависящим от времени: 

library(deSolve)

FN <- function (t, state, parms, K_fun) {
  with(as.list(c(state, parms)), {
    K <- K_fun(t)
    dN <- r * N * (1 - N / K)
    list(c(dN))
  })
}

parms <- c(r = 1)
state <- c(N = 1)

times <- seq(0, 50)
value <- seq(10, 20, length.out = length(times))
K_fun <- approxfun(times, value, rule = 2)
out <- ode(y = state, times = times, func = FN, parms = parms, K_fun = K_fun)
plot(out)

Подробнее об этом можно узнать на странице справки ?events пакета deSolve и в следующем учебное пособие: https://tpetzoldt.github.io/deSolve-forcing/deSolve-forcing.html

...