Solve variable coefficients second order linear ODE?

Question

For the variable coefficients second order linear ODE

$x''(t)+\beta_1(t)x'(t)+\beta_0 x(t)=0$

I have the numerical values (in terms of vectors) for $\beta_1(t)$ and $\beta_0(t)$, does anyone know some R package to do that? And some simple examples to illustrate would be great as well.

I googled to find 'bvpSolve' can solve constant coefficients value.

you should be able to do this with the standard `ode()` function in the `deSolve` package, I think ... the gradient function you need to define takes the arguments `t`, `y`, `parms`, where `t` is the time. The tricky part, I guess, is deciding whether/how you want to interpolate your forcing functions -- are they piecewise constant? piecewise linear? smooth? — Ben Bolker, Apr 02 '13 at 07:12
@BenBolker Could you please give some details? And there is no forcing function in my current example. And for my coefficients functions, it is an expansion from the b-splines. — Lerong, Apr 02 '13 at 14:49
By "forcing functions", I (sloppily) meant your time-dependent parameters, not a time-dependent constant term ... could you give a reproducible example (see http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example ) of a simple set of equations and parameters? — Ben Bolker, Apr 02 '13 at 14:54
it's still not clear to me what "I have the numerical values (in terms of vectors) for $\beta_1(t)$ and $\beta_0(t)$" means. I *may* have answered your question below. — Ben Bolker, Apr 02 '13 at 19:30

Ben Bolker · Answer 1 · 2013-04-02T19:35:50.273

In order to use deSolve, you have to make your second-order ODE

x''(t) + \beta_1(t) x'(t) + \beta_0 x(t)=0

into a pair of coupled first-order ODEs:

x'(t) = y(t)
y'(t) = - \beta_1(t) y(t) - \beta_0 x(t)

Then it's straightforward:

gfun <- function(t,z,params) {
      g <- with(as.list(c(z,params)),
       {
         beta0 <- sin(2*pi*t)
         beta1 <- cos(2*pi*t)
       c(x=y,
         y= -beta1*y - beta0*x))
     list(g,NULL)
}
library("deSolve")
run1 <- ode(c(x=1,y=1),times=0:40,func=gfun,parms=numeric(0))

I picked some initial conditions (x(0)=1, x'(0)=1) arbitrarily; you might also want to add parameters to the model (i.e. make parms something other than numeric(0))

PS if you're not happy doing the conversion to coupled first-order ODEs by hand, and want a package that will seamlessly handle second-order ODEs, then I don't know the answer ...

Solve variable coefficients second order linear ODE?

1 Answers1

Linked