Search code examples
rmatlabfunctionode

Solving Linear System in R


I'm currently trying to solve a 2D linear system in R. I have the below function:

times <- seq(0,25,.01)
state.linear <- c(X = 2, Y = 0) # Initial State

A <- matrix(c(-.1, 2, -2, -.1), 2, 2)
b <- c(2,0)

  linear2D <- function(t, A, b) {
    with(as.list(c(A, b)), {
      dX <- A*b
      list(c(dX))
    })
  }
out.linear2D <- ode(y = state.linear, times = times, func = linear2D, parms = A)

I'd like to put it through the ODE solver from the deSolve package but am running into trouble with the initial function. I'm trying to translate the code from MATLAB:

n = 2;          % 2D system
A = [-.1 2; -2 -.1];  % dynamics
rhs = @(x)A*x;   % right hand side
tspan=[0:.01:25];   % time span
x0 = [2; 0];        % initial conditions
options = odeset('RelTol',1e-10,'AbsTol',1e-10*ones(1,n));
[t,x]=ode45(@(t,x)rhs(x),tspan,x0,options);

But am pretty unclear how exactly this should be translated. If anyone could help translating the MATLAB code to R in order to put it through the ODE function in R that'd be a big help.

Thanks in advance.


Solution

  • Does the following answer your question? Note the following:

    • %*%is the matrix product
    • the default method = "lsoda" would be better than method = "ode45"
    library("deSolve")
    
    times <- seq(0, 25, .01)
    state.linear <- c(X = 2, Y = 0) # Initial State
    
    A <- matrix(c(-.1, 2, -2, -.1), 2, 2)
    b <- c(2, 0)
    
    linear2D <- function(t, A, b) {
      with(as.list(c(A, b)), {
        dX <- A %*% b
        list(c(dX))
      })
    }
    
    out.linear2D <- ode(y = state.linear, times = times, 
                        func = linear2D, parms = A, method="ode45",
                        atol = 1e-10, rtol = 1e-10)
    
    ## time series
    plot(out.linear2D)
    
    ## phase plot
    plot(out.linear2D[, -1], type="l")