# implementation of method of line and Runge-Kutta to the given equation

$$\frac{d^2 u}{dx^2} +A \frac{d^2}{dx^2}\left(\frac{du}{dt}\right)=B$$

I want to solve the equation given above. I need to first discretize it by the Method of Lines and then evolve the resulting ODE using a Runge-Kutta method. My question is, how can I reduce this special equation to two different ODEs?

You can convert your equation into the following system of two PDE

$$\frac{\partial u}{\partial t} = v$$ $$0= \frac{\partial^2 u}{\partial x^2} + A \frac{\partial^2 v}{\partial x^2} - B$$

There are now two dependent variables, $u$ and $v$.

The discretization in the spatial dimension (i.e. method of lines (MOL)) is straightforward; central difference will work. The resulting set of equations is a system of differential algebraic equations (DAE) rather than an ODE. There are many available methods for solving DAE but explicit methods, like classical Runge-Kutta, are not typically used.

A widely used DAE solver is IDA from the Sundials suite. Another well-regarded DAE solver is an implicit Runge-Kutta solver, RADAU5. If you are determined to code your own solution, the most basic implicit ODE solver, backward Euler, will also likely solve this DAE system but will, of course, require a relatively small time step for an accurate solution.

Finally, if you have access to Matlab or Octave, there are 1D PDE solvers that will solve these two equations and let you avoid the difficulties of coding your own MOL solution.

• Hi Bill, i have tried to solve it before just you have suggested but i have two different variables (u,v) and wrote it as below. (u_(i+1,k+1)-2u_(i,k+1)+u_(i-1,k+1))/h^2 +A (v_(i+1,k+1)-2v_(i,k+1)+v_(i-1,k+1))/h^2 =B, but i do not know how two code it in case of two variables. what is the idea? i have just started to work on such a systems. Jan 9, 2017 at 15:15

As an example, you should be able to build up a system of equations to solve for the time derivatives if you use something like the following scheme:

\begin{align} A \frac{\partial ^2}{\partial x^2}\left(\frac{\partial u_i}{\partial t}\right) &= B - \frac{\partial^2 u_i}{\partial x^2}\\ % \frac{A}{\Delta x^2} \frac{\partial u_{i-1}}{\partial t} - \frac{2A}{\Delta x^2}\frac{\partial u_i}{\partial t} + \frac{A}{\Delta x^2} \frac{\partial u_{i+1}}{\partial t} &\approx B - \frac{u_{i-1} - 2u_{i} + u_{i+1}}{\Delta x^2}\\ % A \frac{\partial u_{i-1}}{\partial t} - 2A\frac{\partial u_i}{\partial t} + A \frac{\partial u_{i+1}}{\partial t} &\approx B \Delta x^2 - \left(u_{i-1} - 2u_{i} + u_{i+1}\right)\\ \end{align}

where $u_i$ is the solution at the $i^{th}$ spatially discretized location at some time and where this equation holds for all discretized locations. With the final equation written, you should then be able to put together what should be a tridiagonal system of equations so you can obtain the values for the time derivatives.