Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Getting rid of for loop

Subject: Getting rid of for loop

From: Peder Schmedling

Date: 27 Sep, 2010 14:12:04

Message: 1 of 12

Hi,
can someone please help me get rid of the for loop in this code?

I tried but got stuck on:
allZ = rand(2,N);
allZ(allZ<=0.025)=fernM1*z;


CODE -------------------

%
% Graphics demo
%
figure(1) % Establish a graphic window - use the same on each time
clf % Clear it and set the axis
axis([-5 5 0 10]);
hold on % Points accumulate in the window as an animation

% Initial values
z = [0;1];

%
% 4 Matrices are key to the IFS Fractal description
%
fernM1 = [0 0 ; 0 0.16];
fernM2 = [0.2 -0.26 ; 0.23 0.22];
fernM3 = [-0.15 0.28 ; 0.26 0.24];
fernM4 = [0.85 0.04 ; -0.04 0.85];

% How many points are we making this time? Allocate the array
N = 50000;
allZ = zeros(2, N);
for kk = 1:N
    dd = rand();
    if dd <= 0.025
        z = fernM1 * z;
    elseif (dd > 0.025) && (dd <= 0.125)
        z = fernM2 * z + [0 ; 1.6];
    elseif (dd > 0.125) && (dd <= 0.225)
        z = fernM3 * z + [0 ; 0.44];
    else
        z = fernM4 * z + [0 ; 1.6];
    end
    allZ(:, kk) = z;
end

% plot the rest
plot(allZ(1,:), allZ(2,:), '.g','MarkerSize', 4)
display('The fern is done')

END CODE -------------------

Subject: Getting rid of for loop

From: Oleg Komarov

Date: 27 Sep, 2010 16:42:07

Message: 2 of 12

"Peder Schmedling" <peder.schmedling@sensonor.no> wrote in message <i7q8nk$lg7$1@fred.mathworks.com>...
> Hi,
> can someone please help me get rid of the for loop in this code?
>
> I tried but got stuck on:
> allZ = rand(2,N);
> allZ(allZ<=0.025)=fernM1*z;
>
>
> CODE -------------------
>
> %
> % Graphics demo
> %
> figure(1) % Establish a graphic window - use the same on each time
> clf % Clear it and set the axis
> axis([-5 5 0 10]);
> hold on % Points accumulate in the window as an animation
>
> % Initial values
> z = [0;1];
>
> %
> % 4 Matrices are key to the IFS Fractal description
> %
> fernM1 = [0 0 ; 0 0.16];
> fernM2 = [0.2 -0.26 ; 0.23 0.22];
> fernM3 = [-0.15 0.28 ; 0.26 0.24];
> fernM4 = [0.85 0.04 ; -0.04 0.85];
>
> % How many points are we making this time? Allocate the array
> N = 50000;
> allZ = zeros(2, N);
> for kk = 1:N
> dd = rand();
> if dd <= 0.025
> z = fernM1 * z;
> elseif (dd > 0.025) && (dd <= 0.125)
> z = fernM2 * z + [0 ; 1.6];
> elseif (dd > 0.125) && (dd <= 0.225)
> z = fernM3 * z + [0 ; 0.44];
> else
> z = fernM4 * z + [0 ; 1.6];
> end
> allZ(:, kk) = z;
> end
>
> % plot the rest
> plot(allZ(1,:), allZ(2,:), '.g','MarkerSize', 4)
> display('The fern is done')
>
> END CODE -------------------

A hint:

N = 50000;
allZ = zeros(2, N);
dd = rand(N,1);

idx = dd <= 0.025;
allZ(:,idx) = repmat(fernM1 * z,1,sum(idx));

idx = (dd > 0.025 & dd <= 0.125);
allZ(:,idx) = repmat(fernM2 * z + [0 ; 1.6],1,sum(idx));

And so on

Oleg

Subject: Getting rid of for loop

From: Sean

Date: 27 Sep, 2010 16:55:23

Message: 3 of 12

"Peder Schmedling" <peder.schmedling@sensonor.no> wrote in message <i7q8nk$lg7$1@fred.mathworks.com>...
> Hi,
> can someone please help me get rid of the for loop in this code?
>
> I tried but got stuck on:
> allZ = rand(2,N);
> allZ(allZ<=0.025)=fernM1*z;
>
>
> CODE -------------------
>
> %
> % Graphics demo
> %
> figure(1) % Establish a graphic window - use the same on each time
> clf % Clear it and set the axis
> axis([-5 5 0 10]);
> hold on % Points accumulate in the window as an animation
>
> % Initial values
> z = [0;1];
>
> %
> % 4 Matrices are key to the IFS Fractal description
> %
> fernM1 = [0 0 ; 0 0.16];
> fernM2 = [0.2 -0.26 ; 0.23 0.22];
> fernM3 = [-0.15 0.28 ; 0.26 0.24];
> fernM4 = [0.85 0.04 ; -0.04 0.85];
>
> % How many points are we making this time? Allocate the array
> N = 50000;
> allZ = zeros(2, N);
> for kk = 1:N
> dd = rand();
> if dd <= 0.025
> z = fernM1 * z;
> elseif (dd > 0.025) && (dd <= 0.125)
> z = fernM2 * z + [0 ; 1.6];
> elseif (dd > 0.125) && (dd <= 0.225)
> z = fernM3 * z + [0 ; 0.44];
> else
> z = fernM4 * z + [0 ; 1.6];
> end
> allZ(:, kk) = z;
> end
>
> % plot the rest
> plot(allZ(1,:), allZ(2,:), '.g','MarkerSize', 4)
> display('The fern is done')
>
> END CODE -------------------

Like this?

z = [0;1];
N = 50000;
allZ = zeros(2,N);
dd = rand(1,N);
allZ(:,dd < 0.025) = fernM1*z;
allZ(:,(dd > 0.025) & (dd <= 0.125)) = fernM2*z+[0;1.6];
allZ(:,(dd > 0.125) & (dd <= 0.225)) = fernM3*z+[0;0.44];
allZ(:,(dd > 0.225)) = fernM4*z+[0;1.6];
%Not tested in ML

Subject: Getting rid of for loop

From: Matt Fig

Date: 27 Sep, 2010 16:55:23

Message: 4 of 12

"Oleg Komarov" <oleg.komarovRemove.this@hotmail.it> wrote in message
> A hint:
>
> N = 50000;
> allZ = zeros(2, N);
> dd = rand(N,1);
>
> idx = dd <= 0.025;
> allZ(:,idx) = repmat(fernM1 * z,1,sum(idx));
>
> idx = (dd > 0.025 & dd <= 0.125);
> allZ(:,idx) = repmat(fernM2 * z + [0 ; 1.6],1,sum(idx));
>
> And so on
>
> Oleg

I don't think this will work because z itself is changing in the loop.

Subject: Getting rid of for loop

From: Peder Schmedling

Date: 27 Sep, 2010 17:18:05

Message: 5 of 12

"Matt Fig" <spamanon@yahoo.com> wrote in message <i7qi9r$l4q$1@fred.mathworks.com>...
> "Oleg Komarov" <oleg.komarovRemove.this@hotmail.it> wrote in message
> > A hint:
> >
> > N = 50000;
> > allZ = zeros(2, N);
> > dd = rand(N,1);
> >
> > idx = dd <= 0.025;
> > allZ(:,idx) = repmat(fernM1 * z,1,sum(idx));
> >
> > idx = (dd > 0.025 & dd <= 0.125);
> > allZ(:,idx) = repmat(fernM2 * z + [0 ; 1.6],1,sum(idx));
> >
> > And so on
> >
> > Oleg
>
> I don't think this will work because z itself is changing in the loop.

Thats what I concluded with also (using my limited knowledge) but I was hoping that there were some gurus out there who could prove me wrong..

Oleg: I tried your code and it's kind of up the alley I was going. But it fails, probably due to what Matt is pointing out..

Subject: Getting rid of for loop

From: Sean

Date: 27 Sep, 2010 17:49:12

Message: 6 of 12

> > I don't think this will work because z itself is changing in the loop.
>
> Thats what I concluded with also (using my limited knowledge) but I was hoping that there were some gurus out there who could prove me wrong..
>
> Oleg: I tried your code and it's kind of up the alley I was going. But it fails, probably due to what Matt is pointing out..

That is true and my example will not work for the same reason.

Subject: Getting rid of for loop

From: Jan Simon

Date: 27 Sep, 2010 20:08:03

Message: 7 of 12

Dear Peder,

> can someone please help me get rid of the for loop in this code?

Do you want to get rid of the loop of accelerate the function?

For the later case:
1. Move the RAND call outside the loop:
  d = rand(1, N);
  ...
  dd = d(kk);

2. For the most likely interval, the most IF checks are processed. Re-order them:
  if dd > 0.225 % Catch 77% of cases at first
    z = fernM4 * z + [0 ; 1.6];
  elseif (dd > 0.125) % Not needed: && (dd <= 0.225)
     z = fernM3 * z + [0 ; 0.44];
  elseif (dd > 0.025) % Not needed: && (dd <= 0.125)
      z = fernM2 * z + [0 ; 1.6];
  else % if dd <= 0.025
     z = fernM1 * z;
  end

Think: In the real life, you'd never check if d<0.125, if you have already sure, that d is not >=0.125 !

But this saves just 10% computing time on my Matlab2009a for the computations - the drawing needs much more time at all.

Kind regards, Jan

Subject: Getting rid of for loop

From: Bruno Luong

Date: 27 Sep, 2010 21:05:23

Message: 8 of 12

I think Matlab does not have recursion built-in function beside FILTER, which cannot process varying coefficients like yours. The only enhancement I can think of if get rid of calling RAND within the loop:

% Initial values
z = [0;1];

%
% 4 Matrices are key to the IFS Fractal description
%
fernM1 = [0 0 ; 0 0.16];
fernM2 = [0.2 -0.26 ; 0.23 0.22];
fernM3 = [-0.15 0.28 ; 0.26 0.24];
fernM4 = [0.85 0.04 ; -0.04 0.85];

N = 50000;
allZ = zeros(2, N);
dd = rand(1,N);
[~, where] = histc(dd, [0 0/025 0.125 0.225 1]);
fern = cat(3,fernM1,fernM2,fernM3,fernM4);
a = cat(2, [0; 0], [0 ; 1.6], [0 ; 0.44], [0 ; 1.6]);
for kk = 1:N
    w = where(kk);
    z = fern(:,:,w)*z + a(:,w);
    allZ(:,kk) = z;
end

% Bruno

Subject: Getting rid of for loop

From: Bruno Luong

Date: 27 Sep, 2010 21:11:20

Message: 9 of 12

Sorry for that typo at the line:

[~, where] = histc(dd, [0 0.025 0.125 0.225 1]);

% Bruno

Subject: Getting rid of for loop

From: Jan Simon

Date: 27 Sep, 2010 23:23:06

Message: 10 of 12

Dear Bruno,

> for kk = 1:N
> w = where(kk);
> z = fern(:,:,w)*z + a(:,w);
> allZ(:,kk) = z;
> end

I've tried exactly this, but was slower than with the IF method on Matlab 2009a and a single core CPU. Perhaps it is a benefit to omit the "+ a(:, w)" if it is zero.
Is it faster on your machine?

BTW. I'm impressed by the resulting picture anytime I see it again.

Kind regards, Jan

Subject: Getting rid of for loop

From: Pedro80

Date: 28 Sep, 2010 05:14:52

Message: 11 of 12

On Sep 28, 1:23 am, "Jan Simon" <matlab.THIS_Y...@nMINUSsimon.de>
wrote:
> Dear Bruno,
>
> > for kk = 1:N
> >     w = where(kk);
> >     z = fern(:,:,w)*z + a(:,w);
> >     allZ(:,kk) = z;
> > end
>
> I've tried exactly this, but was slower than with the IF method on Matlab 2009a and a single core CPU. Perhaps it is a benefit to omit the "+ a(:, w)" if it is zero.
> Is it faster on your machine?
>
> BTW. I'm impressed by the resulting picture anytime I see it again.
>
> Kind regards, Jan

Yes, me too :-)
This should be the Barnsley fern: http://en.wikipedia.org/wiki/Barnsley's_fern

Maybe one of you guys can come up with a "Matlab way" of computing the
fern with N number of data points when yuo see the actual model on the
wiki page..

Anyways thank you all for your input so far, I'm learning by the hour!

Subject: Getting rid of for loop

From: Bruno Luong

Date: 28 Sep, 2010 05:56:09

Message: 12 of 12

"Jan Simon" <matlab.THIS_YEAR@nMINUSsimon.de> wrote in message <i7r90q$7sn$1@fred.mathworks.com>...
> Dear Bruno,
>
> > for kk = 1:N
> > w = where(kk);
> > z = fern(:,:,w)*z + a(:,w);
> > allZ(:,kk) = z;
> > end
>
> I've tried exactly this, but was slower than with the IF method on Matlab 2009a and a single core CPU. Perhaps it is a benefit to omit the "+ a(:, w)" if it is zero.
> Is it faster on your machine?
>

Yes it's 30% faster on my machine. I think the extra time is to copy the slice of array fern at each loop, not because of "+ a(:, w)". Here is an even faster way on my machine because no copy is necessary:

% Initial values
z = [0;1];

%
% 4 Matrices are key to the IFS Fractal description
%
fernM1 = [0 0 ; 0 0.16];
fernM2 = [0.2 -0.26 ; 0.23 0.22];
fernM3 = [-0.15 0.28 ; 0.26 0.24];
fernM4 = [0.85 0.04 ; -0.04 0.85];

tic
N = 50000;
allZ = zeros(2, N);
dd = rand(1,N);
[~, where] = histc(dd, [0 0.025 0.125 0.225 1]);
fern = { fernM1,fernM2,fernM3,fernM4 };
a = cat(2, [0; 0], [0 ; 1.6], [0 ; 0.44], [0 ; 1.6]);
for kk = 1:N
    w = where(kk);
    z = fern{w}*z + a(:,w);
    allZ(:,kk) = z;
end
toc

% Bruno

Tags for this Thread

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.

Contact us