I am trying to implement a stiff ODE solver that uses step doubling with the backward Euler method but I want to make a modification and I am unsure how to incorporate it.

I am trying to solve $y'=f(y), y(a)=d$ and I am using an approximation to $f$ such that $y'\approx F(y),\;t\in (T,T+dt), y(T)=Y$ where

$$F(y)=f(y)+Df(Y)(y(t) - Y).$$

I understand how to implement the backward Euler method when we simply have $y'=f(t,y)$ but I am having trouble writing a program to solve the problem here because $F(y)$ is dependent on $y(t)$, which is obviously unknown so the usual formulation of $\Delta x=Df^{-1}(x)f(x)$ doesn't seem to work. I was wondering how I could go about writing a method for this?

  • 1
    $\begingroup$ Don't you mean $F(y) = f(Y) + Df(Y)(y - Y)$ since the Taylor Series is about $Y$? Assuming that, your problem becomes $y' = (f(Y) - AY) + Ay$, where $A = Df(Y)$. The first term is independent of $y$, making your equations linear with respect to $y$ since the latter term is linear. Should be pretty straight forward to do implicit integration then! $\endgroup$ – spektr Dec 5 '16 at 23:12
  • $\begingroup$ Thanks! That helped but honestly I am still a bit confused on how to do the implicit integration? The way I calculated I got that $Df(Y)\Delta y=-f(y)$ so I solve that using Newton-Raphson and then do $\Delta y + y_t=y_{t+1}$ but that seems to be giving me the wrong answer. Is there anything I am missing? Sorry, I am very new to this haha. $\endgroup$ – Twis7ed Dec 7 '16 at 3:30
  • $\begingroup$ Sorry, meant to say that $\Delta y(I-Df(Y)\Delta t)=\Delta t f(Y)$ $\endgroup$ – Twis7ed Dec 7 '16 at 3:37

Let's assume you have some differential equation that can be written as the following:

$$ \frac{dy}{dt} = f(t,y)$$

for some initial condition $y(a) = d$. We then want to integrate this quantity using Backward Euler, with the fundamental formula being derived in the following way:

$$ \begin{align} \frac{dy}{dt} &= f(t_{n+1},y_{n+1})\\ \frac{y_{n+1} - y_{n}}{\Delta t} &= f(t_{n+1},y_{n+1})\\ y_{n+1} &= y_{n} + \Delta t f(t_{n+1},y_{n+1}) \end{align}$$

Now we might want to make things a little easier than solving nonlinear equations for $y_{n+1}$, if we're willing to add in a little bit of error, by linearizing $f(t_{n+1},y_{n+1})$ around $y_{n}$. We can show the linearization for $f(t_{n+1},y_{n+1})$ about $y_{n}$ is:

$$ f(t_{n+1},y_{n+1}) = f(t_{n+1},y_{n} + (y_{n+1} - y_{n})) \approx f(t_{n+1},y_{n}) + \left.\frac{\partial f}{\partial y}\right|_{(t_{n+1},y_{n})}\left(y_{n+1} - y_{n}\right)$$

If we substitute this into the Backward Euler equation, and set $A_n = \left.\frac{\partial f}{\partial y}\right|_{(t_{n+1},y_{n})}$, we can work out the following: $$ \begin{align} y_{n+1} &= y_{n} + \Delta t f(t_{n+1},y_{n}) + \Delta tA_n\left(y_{n+1} - y_{n}\right) \\ (I + \Delta t A_n)y_{n+1} &= (I - \Delta t A_n)y_{n} + \Delta t f(t_{n+1},y_{n}) \\ y_{n+1} &= y_{n} + \Delta t (I - \Delta t A_n)^{-1}f(t_{n+1},y_{n}) \end{align}$$

The final formula above should be all you need, assuming you can compute $A_n$.

| cite | improve this answer | |
  • $\begingroup$ Thank you so much! This was extremely helpful and helped me learn why the motivation behind the modification. $\endgroup$ – Twis7ed Dec 9 '16 at 4:07
  • $\begingroup$ @Twis7ed no problem! Glad it helped! $\endgroup$ – spektr Dec 9 '16 at 4:11

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.