# Solving differential equations using implicit methods

• Sagekilla

#### Sagekilla

Hi all, I'm writing myself a ordinary differential equation solver and I've already implemented several explicit integrators, which were pretty easy for me to do. Now I've decided to work on some implicit methods (for any stiff equations) and I've run into some issues.

The most basic one is the Backward Euler:
$$y_{n+1} = y_{n} + hf(t_{n+1}, y_{n+1})$$

And as I understand, to solve this you need to use some iteration to get the value of $$y_{n+1}$$ to solve that equation for the next time step.

The problem I have is: How can I solve this? I thought I could use something like Newton's method, but I don't know how to properly apply it.

Any help, even a push in the right direction would be appreciated. Thanks!

Last edited:

## Answers and Replies

The backward method
$$y_{n+1} = y_{n} + hf(t_{n+1}, y_{n+1})$$

need to solve for yn+1 via NR method.

$$y_{n+1} - y_{n} - hf(t_{n+1}, y_{n+1})=0$$

Let $$g(z) = z - y_{n} - hf(t_{n+1}, z)$$

The iteration
$$z_{i+1} = z_i - \frac{g(z_i)}{g'(z_i)}$$
should converges to yn+1 using the initial estimation z0=yn.

Thanks for your help! I never understood how to properly apply Newton's method before. I implemented it as you had shown, with:

$$z_{i+1} = z_{i} - \frac{z_{i} - y_{n} - hf(t_{n+1}, z_{i})}{1 - hf'(t_{n+1}, z_{i})}$$

Is it safe to assume that I can say that: $$f'(t_{n+1}, z_{i}) = f(t_{n+1}, f(t_{n+1}, z_{i}))$$

I thought it out for a few minutes, and it seemed reasonable in the context of a few examples I put together:

If we let $$\frac{dy}{dt} = -ky$$ then $$f(t, y) = -ky$$ and $$f'(t, y) = -kf(t, y) = f(t, f(t, y))$$

I thought that we have this relation: tn+1=t n + h
where h is the step size increment.
In the case of $\frac{dy}{dt} = -ky$ , f(tn+1, z) = -kz

Yes, there is. But what I was asking was if I could replace the derivative of f inside of Newton's method by f(t, f(t, y)), since for first order ODEs it would be equivalent to differentiating (at least I think it is) all the y components of differential equation.

Because otherwise if I couldn't do that, then I know I'd have to numerically approximate the derivative as well.

On a second thought, I think it is easier to use the fixed-point interation
zn+1=g(zn)
to approximates yn+1.

This eliminate the requirement to compute the derivative of f.