CODEWITHSUNDEEP
rmatlabfor-loop
AW27
AW27

Reputation: 505

Matlab to R For Loop Conversion

I'm currently trying to convert MATLAB code to R and run a loop to find the derivative of an equation.

MATLAB:

for i=1:length(x)
    dx(i,:) = lorenz(0,x(i,:),sigma,beta,rho);
end

My attempted code in R:


parameters <- c(s = 10, r = 28, b = 8/3)
state <- c(X = -8, Y = 7, Z = 27)

lorenz <- function(t, state, parameters) {
    with(as.list(c(state, parameters)), {
        dX <- s * (Y - X)
        dY <- X * (r - Z) - Y
        dZ <- X * Y - b * Z
        list(c(dX, dY, dZ))
    })
}

times <- seq(0.001, 1, by = 0.001)

out <- ode(y = state, times = times, func = lorenz, parms = parameters)

dx <- NULL
for (i in 1:1000){
    dx[i,] <- lorenz(0,state[i,],parameters)
}

I keep receiving the error "Error in state[i, ] : incorrect number of dimensions", any idea why?

Upvotes: 0

Views: 114

Answers (1)

stefan
stefan

Reputation: 124128

Try this

library(deSolve)

parameters <- c(s = 10, r = 28, b = 8/3)
state <- c(X = -8, Y = 7, Z = 27)

lorenz <- function(t, state, parameters) {
  with(as.list(c(state, parameters)), {
    dX <- s * (Y - X)
    dY <- X * (r - Z) - Y
    dZ <- X * Y - b * Z
    list(c(dX, dY, dZ))
  })
}

times <- seq(0.001, 1, by = 0.001)

out <- ode(y = state, times = times, func = lorenz, parms = parameters)

# Set up matrix
dx <- matrix(NA, nrow = 1000, ncol = 3)
for (i in 1:1000){
  # lorenz returns a list with one element. To assign to dx you have extract the list element using [[1]]
  dx[i,] <- lorenz(0,out[i,],parameters)[[1]]
}
head(dx)
#>                 [,1]          [,2]         [,3]
#>    [1,] 150.00000000  -15.00000000 -128.0000000
#>    [2,] 148.35404862  -15.83440159 -126.4953894
#>    [3,] 146.71639504  -16.62054372 -125.0045611
#>    [4,] 145.08747674  -17.35986863 -123.5280889
#>    [5,] 143.46759235  -18.05374748 -122.0664153
#>    [6,] 141.85715318  -18.70357479 -120.6200236

Created on 2020-03-20 by the reprex package (v0.3.0)

Upvotes: 2

Related Questions