Thread Subject: vectorizing a for loop in Euler's method (ODE solver)

Subject: vectorizing a for loop in Euler's method (ODE solver)

From: Phil Isaac

Date: 11 Aug, 2009 02:28:05

Message: 1 of 2

As an educational exercise, I'd like to know if the for loop in the following M-file code can be vectorized:

%%%
t0=0;
y0=4;
tMax = 5;
steps= 20;
rhsFunc = @(t,y)(t.* y.^2 .* sin(2*y.*t));
y=zeros(1,steps+1);
y(1) = y0;
deltaT = (tMax-t0)/steps;
t = t0:deltaT:tMax;
for i = 1:steps
  y(i+1) = y(i) + rhsFunc(t(i),y(i))*deltaT;
end
plot(t,y);
%%%

This is an implementation of Euler's method, a simple first order ODE solver.

I understand how to vectorize a loop using "filter" (or in some cases "cumsum" or "cumprod") that creates a vector whose elements depend on the previous element, in the case where those elements depend *linearly* on the previous elements. That is explained clearly in Support - Code Vectorization Guide:

http://www.mathworks.com/support/tech-notes/1100/1109.shtml

In the example code that I have provided, I have purposely made "rhsFunc" a non-linear function of y. I am having trouble working out how to vectorize that loop in this case, or indeed in the general case where "rhsFunc" is any (smooth) non-linear function of two variables, in particular the variable y.

Some additional comments:
- I know about the ODE solvers ode23, ode45, etc. I'm more interested in how to vectorize, rather than the context.
- The code I have provided already works fine. I have used a small number of steps, but if you change to "steps = 20000;", the plot produced is not bad (compare with ode45) and still runs in good time (of course not as efficient as Runge-Kutta!), and so there may be no point in vectorizing at all. Again, I'm interested in the question of how to vectorize, rather than the context.

Cheers,

Phil

Subject: vectorizing a for loop in Euler's method (ODE solver)

From: Bruno Luong

Date: 11 Aug, 2009 07:44:03

Message: 2 of 2

I have my doubt non-linear iterative can be vectorized.

Bruno

Tags for this Thread

Everyone's Tags:

Add a New Tag:

Separated by commas
Ex.: root locus, bode

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Tag Activity for This Thread
Tag Applied By Date/Time
vectorization Phil Isaac 10 Aug, 2009 22:29:05
rssFeed for this Thread
 

MATLAB Central Terms of Use

NOTICE: Any content you submit to MATLAB Central, including personal information, is not subject to the protections which may be afforded information collected under other sections of The MathWorks, Inc. Web site. You are entirely responsible for all content that you upload, post, e-mail, transmit or otherwise make available via MATLAB Central. The MathWorks does not control the content posted by visitors to MATLAB Central and, does not guarantee the accuracy, integrity, or quality of such content. Under no circumstances will The MathWorks be liable in any way for any content not authored by The MathWorks, or any loss or damage of any kind incurred as a result of the use of any content posted, e-mailed, transmitted or otherwise made available via MATLAB Central. Read the complete Terms prior to use.

Contact us at files@mathworks.com