Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

I have tried to solve a system of equations in the form:

$$\eqalign{ {dy_1\over dt }&= y_1 +dt*f(y_1,y_2)\cr {dy_2\over dt} &= y_2 +dt*g(y_1,y_2) } $$

using different schemes such as forward Euler and backward (implicit) Euler and Runge Kutta order 4. My results seem to be OK for the forward Euler scheme and Runge Kutta above a certain number of time steps (A). But for implicit Euler using A time steps, I get an error after a certain number of time steps (say A-100). I was wondering if anyone had any ideas of what the problem might be. I can't think of anything at all.

Thanks!


Sorry, the correct equations are:

$$\eqalign{ {dy_1\over dt }&= f(y_1,y_2)\cr {dy_2\over dt} &= g(y_1,y_2) } $$

where $y_1$ and $y_2$ have a size $n$ by $1$ ($n$ an integer greater than 1). I also just realized that if I use much smaller than that of explicit euler my reuslts become OK. I am now thinking it could be because I used forward euler to work out the initial RHS's at the new time step, to calculate the $y_1$ and $y_2$ at the new time steps.

share|improve this question
    
Could you be a bit specific and tell us the system that you were trying to solve with backward Euler? –  J. M. Jan 1 '12 at 14:04
    
Do you really want the $dt$'s on the left hand side? –  David Mitra Jan 1 '12 at 14:08
    
Sorry,The correct equations are: dy1/dt = f(y1,y2) dy2/dt = g(y1,y2) where y1 and y2 have a size n by 1 (n an integer greater than 1). I also just realized that if I use much smaller than that of explicit euler my reuslts become OK. I am now thinking it could be because I used forward euler to work out the initial RHS's at the new time step, to calculate the y1 and y2 at the new time steps. –  Hooman Jan 3 '12 at 13:29

1 Answer 1

up vote 2 down vote accepted

So basically what you get by applying Backward Euler is the following two equations, \begin{equation} y_{1}^{n+1} = y_{1}^n + h * f(y_1^{n+1}, y_2^{n+1}), \end{equation} \begin{equation} y_{2}^{n+1} = y_{2}^n + h * g(y_1^{n+1}, y_2^{n+1}). \end{equation} Theoretically, using this scheme your method should be more stable. So maybe you can start testing your solution by using larger time steps and checking if the method stays stable. Considering what you explained above I suspect either you're having some bug in your code or you are introducing some errors in the initial values by using Forward Euler to start the computations which is not allowed. $f$ and $g$ are both known functions I assume.

share|improve this answer
    
Thanks! f and g are known function. What other methods can be used to obtain the initial values? My Matlab code (just the time loop) at the moment is in this form (n is the number of time steps and h is the number of rows of y1 and y2): y1 = zeros(n,h); y2 = zeros(n,h); for i = 1:n y1old = y1; y2old = y2; f = A1*y1old + B1*y2old; g = A2*y1old + B2*y2old; y1new = y1old+f*dt; y2new = y2old+g*dt; fnew = A1*y1new + B1*y2new; gnew = A2*y1new + B2*y2new; y1 = y1new+fnew*dt; y2 = y2new+gnew*dt; y11(:,i) = y1new; y22(:,i) = y2new; end –  Hooman Jan 3 '12 at 13:28
    
your current comment is a bit hard to read, you should use code wrappers to indicate your code. You can find the syntax in the help section of the comments here. If by initial values you mean $y_1^0, y_2^0$ then from what you said, you are using zero's as the initial values. –  bamdadhosseini Jan 3 '12 at 19:53
1  
Since you are using Backwards Euler, you don't need the initial values of $f(y_1^0, y_2^0, t)$ and $g(y_1^0, y_2^0, t)$ anymore. Since $f$ and $g$ are known you get a L.H.S that is a function of $y_1^{n+1}, y_2^{n+1}$ and a R.H.S which is a function of $y_1^{n}, y_2^{n}$ (you need to move $f$ and $g$ to the left). This leads to a system of equation which should be solved to give $y_1^{n+1}, y_2^{n+1}$. Notice that the initial values are only used once. And you are already using zeros for that matter. –  bamdadhosseini Jan 3 '12 at 19:53
    
Thanks. So that means that the function would need to be solved inside a for loop. So for instance if I'm using Matlab I could use the built in "solve" function to solve a quadratic. But this could become complicated depending on the function and the No of variables.Am I missing a point here again? I would also be interested to know if there is a general approach that can be used with any function. Thanks again. –  Hooman Jan 15 '12 at 19:47
1  
The initial guess for the Newton iteration does not have to be equal to the initial conditions of the Euler method. The Newton iterations are sub iterations that you want to use in each step of Euler's method. –  bamdadhosseini Jan 17 '12 at 2:07

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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