gistfile1.txt

On 11/20/2013 1:59 AM, Brandon Degner wrote:
Radu

For conducting the backward euler method in matlab, we have to run something that will give us an estimate as to what the slope of the next point is.  That something for us I assume is newton raphson.  I’m getting confused as to how we are supposed to setup the equations for newton raphson?  Any information on this would be much appreciated.

Thanks

Brandon Degner


Brandon,

This would have been a great question to ask on the forum!

In any case, here's the story: indeed, since BE is an implicit method, obtaining the solution at each time step requires the solution of a (typically) nonlinear problem. 

Suppose you want to solve using Backward Euler (BE) the following IVP:
y' = f(t,y)
y(t_0) = y0
and suppose you have obtained the (numerical, approximate) solution at time t_n. Now you want to find the solution at time t_{n+1}.  Using the BE integration formula, you get
y_{n+1} = y_n + h * f(t_{n+1}, y_{n+1})
or, rearranging:
y_{n+1} - h * f(t_{n+1},y_{n+1}) - y_n = 0
Since y_n is known at this point (recall, we assumed you have already solved the problem at the previous time step) and since t_{n+1} is also known,  the above represents a nonlinear problem in the unknown y_{n+1}.  

To make things even more clear, denote your unknown by X.  Then, you find X by solving the nonlinear problem:  
g(X) = 0
where
g(X) = X - h*f(t_{n+1},X) - y_n

You can use Newton-Raphson for this or use the Matlab function 'fsolve'.  The pieces that are missing in order to solve this nonlinear problem are:
The partial derivative of g with respect to X.  Assuming for simplicity that X has dimension 1 (like in the assignment), this is:
g'(X) = 1 - h * df/dX
For the assignment problem, where f(t,y) = sin(y), this becomes:
g'(X) = 1 - h * cos(X)

A good initial guess for X at the current time is the solution at the previous time: X_{guess} = y_n
I hope this helps.  Good luck,
--Radu