0
$\begingroup$

I am applying a 4th order Runge-Kutta code, using the method of lines, to solve the following:

\begin{equation} \frac {\partial y_1}{\partial t} = y_2 y_3 - C_1 y_1 \end{equation}

\begin{equation} \frac {\partial y_2}{\partial t} = y_3y_1 - C_2 y_2 \end{equation}

\begin{equation} \frac {\partial y_3}{\partial z} = y_4 \end{equation}

\begin{equation} \frac {\partial y_4}{\partial z} = \frac {\partial ^2 y_2}{\partial t^2} + \frac {\partial ^2 y_3}{\partial t^2} + y_1y_3 - y_2y_3 + y_3^2y_2 - y_4\end{equation}

(where $ y_n$ are the solutions I am looking for, $ C_n $ are real constants)

My code appears to work and solve the system when I just remove the $ y_4 $ term from equation (4), but otherwise, the solution goes to infinity. I have tried rescaling parameters/variables, decreasing the spatial intervals, and changing the boundary/initial conditions (except trivial ones), but these seem not to work. For testing purposes, I also changed equation (4) to:

$$ \frac {\partial y_4}{\partial z} = \frac {\partial ^2 y_2}{\partial t^2} + \frac {\partial ^2 y_3}{\partial t^2} + y_1y_3 - y_2y_3 + y_3^2y_2 - \text {SCALE}*y_4 $$

where SCALE = $10^{-5}$. With this change, the output did not go to infinity, although it had an unexpected output (namely there were frequent large oscillations). I could not change SCALE to values larger than $10^{-5}$, otherwise the output would go to infinity.

RK4 seems to work well for almost everything except this one term in the 4th equation, and I was just wondering if you had any ideas? I was considering possibly applying an implicit method on only the 3rd and 4th equations, although due to the nonlinearity, I haven't been able to isolate for the $y_3$ when using higher order BDF methods.

Improve this question
$\endgroup$
1
  • 2
    $\begingroup$ Runge-Kutta methods solve systems of first-order ODEs. That is not what you have here. You need to rewrite it as a first-order system. $\endgroup$ – David Ketcheson Dec 10 '16 at 18:25
-2
$\begingroup$

This is definitely a stiffness issue. Instead of using RK4, what about using a reasonably well-optimized adaptive timestepping method? In many cases that will work well enough. If not, the places to go are Rosenbrock methods are BDF methods. You don't need to isolate anything to use BDF methods: you know the second time derivatives (since you know the first derivatives), so you can write down the vector function for how everything changes and just plop that into an ODE solver (with all 4 equations). The ODE solver will solve the root finding problem for you.

Improve this answer
$\endgroup$
6
  • $\begingroup$ Thank you for the response. I tried an adaptive tilmestep using the Fehlberg method but it was going too slowly as there were small oscillations present and the time steps were on the order of nanoseconds while my domain is days. I will look into implicit methods, but I don't have an explicit expression for $\frac {\partial ^2 y_2}{\partial t^2}$ (nor the first derivative in time), so will it work? $\endgroup$ – Mathews24 Oct 9 '16 at 16:13
  • $\begingroup$ The same expression that you're plugging into the ODE solver for the RK method should just plug into the BDF method. Also, I don't know what ODE solver you're using but you shouldn't use a Fehlberg method, those are really awful in comparison to modern methods like Dormand-Prince and are never suggested. $\endgroup$ – Chris Rackauckas Oct 9 '16 at 16:17
  • $\begingroup$ Hmmm, interesting. I looked back and it was actually Cash-Karp that I used, but isn't the difference just 4th order vs 5th order? Also, I was under the impression that although lower order methods were less accurate, they had a greater radius of convergence? Maybe I am simply lacking the understanding, but my equation for $y_4$ is a partial differential equation with respect to $z$ (and my equation for $y_3$ is a PDE in $z$), where as BDF methods would require them all to be differential equations of one variable (e.g. $t$), no? $\endgroup$ – Mathews24 Oct 9 '16 at 20:24
  • $\begingroup$ That can be true. Sometimes a higher order method can have the same or a larger radius of convergence. Cash-Karp and RKF are both order 4/5, but they have much different properties. See this post. Any ODE solver you use will propogate it forward in terms of the derivatives you specify it, i.e. if you natively throw this into an ODE solver it will go forward in time for $y_1$ and $y_2$, and in $z$ for $y_3$ and $y_4$. $\endgroup$ – Chris Rackauckas Oct 9 '16 at 20:30
  • $\begingroup$ But what I am saying is, if you made a function f to put into one ODE solver, there's no reason why that same f cannot be used with another method. $\endgroup$ – Chris Rackauckas Oct 9 '16 at 20:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.