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:
Vectorization Help...I think

Subject: Vectorization Help...I think

From: Jonny O'Connell

Date: 17 Jul, 2010 18:06:03

Message: 1 of 9

Hi,

I'm really new to MALTAB, I only sat the module this semester, but somehow I've managed to get myself a summer project writing a parallel evolutionary algorithm in it!

Anyway I've done some basic profiling on my code and found that one line is accounting for about ~45% of my run time. It's a find function, but it's pretty much a copy and paste from the help documentation (obviously adapted for my needs), but I'm not sure if I can improve it any further.

From the literature I've read on MATLAB, vectorization is king, and should be used anywhere it can be for better performance. I've had some success on various parts of the program, but I'm not sure if what I have is a candidate for vectorization and I was hoping someone could cast an experienced eye over it.

% Perform swapping until cycle is complete
    while cycleIndex ~= seedNumber;
        cycleIndex = find(firstParent == secondParent(currentNumber));
        offspring(cycleIndex) = firstParent(cycleIndex);
        currentNumber = cycleIndex;
    end

The culprit is the 'find' line, the first line in the loop.

Any pointers would be good!

Cheers

Subject: Vectorization Help...I think

From: Bruno Luong

Date: 17 Jul, 2010 19:01:06

Message: 2 of 9

"Jonny O'Connell" <jonnydotoconnelldotx@gmail.com> wrote in message <i1sreb$6fs$1@fred.mathworks.com>...
> Hi,
>
> I'm really new to MALTAB, I only sat the module this semester, but somehow I've managed to get myself a summer project writing a parallel evolutionary algorithm in it!
>
> Anyway I've done some basic profiling on my code and found that one line is accounting for about ~45% of my run time. It's a find function, but it's pretty much a copy and paste from the help documentation (obviously adapted for my needs), but I'm not sure if I can improve it any further.
>
> From the literature I've read on MATLAB, vectorization is king, and should be used anywhere it can be for better performance. I've had some success on various parts of the program, but I'm not sure if what I have is a candidate for vectorization and I was hoping someone could cast an experienced eye over it.
>
> % Perform swapping until cycle is complete
> while cycleIndex ~= seedNumber;
> cycleIndex = find(firstParent == secondParent(currentNumber));
> offspring(cycleIndex) = firstParent(cycleIndex);
> currentNumber = cycleIndex;
> end
>

Try this:

[trash loc] = ismember(secondParent,firstParent);
while cycleIndex ~= seedNumber;
    cycleIndex = loc(cycleIndex);
    offspring(cycleIndex) = firstParent(cycleIndex);
end

Bruno

Subject: Vectorization Help...I think

From: Jonny O'Connell

Date: 17 Jul, 2010 22:14:04

Message: 3 of 9

> Try this:
>
> [trash loc] = ismember(secondParent,firstParent);
> while cycleIndex ~= seedNumber;
> cycleIndex = loc(cycleIndex);
> offspring(cycleIndex) = firstParent(cycleIndex);
> end
>
> Bruno

Thank-you! I think you accidentally missed the stopping condition, but after adding that back in I did a quick test. For the whole loop originally, 320 seconds, with your alteration, the loop and the ismember line are only 145....over twice as quick!

Many thanks again,

Jonny

Subject: Vectorization Help...I think

From: Roger Stafford

Date: 17 Jul, 2010 23:05:06

Message: 4 of 9

"Jonny O'Connell" <jonnydotoconnelldotx@gmail.com> wrote in message <i1t9vc$k9s$1@fred.mathworks.com>...
> > Try this:
> >
> > [trash loc] = ismember(secondParent,firstParent);
> > while cycleIndex ~= seedNumber;
> > cycleIndex = loc(cycleIndex);
> > offspring(cycleIndex) = firstParent(cycleIndex);
> > end
> >
> > Bruno
>
> Thank-you! I think you accidentally missed the stopping condition, but after adding that back in I did a quick test. For the whole loop originally, 320 seconds, with your alteration, the loop and the ismember line are only 145....over twice as quick!
>
> Many thanks again,
>
> Jonny
- - - - - - -
  It is difficult to give you a reliable answer without some additional information about the firstParent and secondParent arrays. For example, it would be easy to fall into an infinite loop if, say,

 firstParent(1) == secondParent(2)
 firstParent(2) == secondParent(1)

Also if more than one firstParent is equal to a given secondParent element, then the find function could yield multiple values, or if no firstParent is found, it could produce empty results. I think you need to tell us more about the nature of these two arrays.

  The ismember function presumably uses an order n*log(n) algorithm where n is the arrays' length, whereas the find function is an order n algorithm. For Bruno's suggestion to save time with large n, it would be necessary that either a single use of the while loop is made and the number of trips through it before exiting is likely to be at least order log n, or some equivalent multiple usage of the while loop is made using the results of a single call to ismember.

  I think there is much more that you can tell us about the statistics of your problem that would be useful in giving you reliable advice.

Roger Stafford

Subject: Vectorization Help...I think

From: Jonny O'Connell

Date: 18 Jul, 2010 02:36:03

Message: 5 of 9

> It is difficult to give you a reliable answer without some additional information about the firstParent and secondParent arrays. For example, it would be easy to fall into an infinite loop if, say,
>
> firstParent(1) == secondParent(2)
> firstParent(2) == secondParent(1)
>
> Also if more than one firstParent is equal to a given secondParent element, then the find function could yield multiple values, or if no firstParent is found, it could produce empty results. I think you need to tell us more about the nature of these two arrays.
>
> The ismember function presumably uses an order n*log(n) algorithm where n is the arrays' length, whereas the find function is an order n algorithm. For Bruno's suggestion to save time with large n, it would be necessary that either a single use of the while loop is made and the number of trips through it before exiting is likely to be at least order log n, or some equivalent multiple usage of the while loop is made using the results of a single call to ismember.
>
> I think there is much more that you can tell us about the statistics of your problem that would be useful in giving you reliable advice.
>
> Roger Stafford


Okay I'll try and put what I'm doing into a bit of context. This function takes two permutations of size n (size normally in the range 100-500) generated using the 'randperm' function, to be the firstParent and secondParent.

Lets say we have the two parents:
First: 8 7 6 4 1 2 5 3 9 10
Second: 2 5 1 7 3 8 4 6 10 9

Taking a randomly selected starting point, that element from the first parent is copied into an 'offspring'. So say this algorithm chooses to start at index 3, 6 is copied into location 3 of an offspring.

Offspring: - - 6 - - - - - - -

The algorithm then looks at index location 3 in the second parent, and finds 1. It locates 1 in the first parent and copies that into the offspring, at its first parent index of 5.

Offspring: - - 6 - 1 - - - - -

It repeats, looks at index 5 in the second parent and finds 3. 3 is located in the first parent at index 8 and copied too the offspring.

Offspring - - 6 - 1 - - 3 - -

Now when the algorithm comes onto the next loop of the cycle, and checks index 8 of the second parent it finds 6 which is the value it started with. At this point the cycle halts and any empty spaces are filled with values from the second parent.

Offspring: 2 5 6 7 1 8 4 3 10 9

My entire function to achieve this:

function offspring = cycle_crossover(firstParent, secondParent)
   
    % Initilisation of starting point and output matrix
    seedNumber = randi(length(firstParent), 1, 1);
    offspring = zeros(1, length(firstParent));
    offspring(seedNumber) = firstParent(seedNumber);
    currentNumber = seedNumber;
    cycleIndex = 0;

    % Perform swapping until cycle is complete
    [~, locations] = ismember(secondParent, firstParent);
    while cycleIndex ~= seedNumber;
        cycleIndex = locations(currentNumber);
        offspring(cycleIndex) = firstParent(cycleIndex);
        currentNumber = cycleIndex;
    end
    
    % Keep remaining items in their original location
    zeroLocations = ismember(offspring, 0);
    offspring(zeroLocations) = secondParent(zeroLocations);
end


I changed the 'copy remaining items' section code also based on what Bruno provided, and that too seemed to speed up significantly.

For the scenario laid out, is that how you believe the code should look?

Thanks again!

Jonny

Subject: Vectorization Help...I think

From: Roger Stafford

Date: 18 Jul, 2010 05:10:06

Message: 6 of 9

"Jonny O'Connell" <jonnydotoconnelldotx@gmail.com> wrote in message <i1tpaj$kl0$1@fred.mathworks.com>...
> > ......
> > I think there is much more that you can tell us about the statistics of your problem that would be useful in giving you reliable advice.
> > Roger Stafford
>
> Okay I'll try and put what I'm doing into a bit of context. This function takes two permutations of size n (size normally in the range 100-500) generated using the 'randperm' function, to be the firstParent and secondParent.
>
> Lets say we have the two parents:
> First: 8 7 6 4 1 2 5 3 9 10
> Second: 2 5 1 7 3 8 4 6 10 9
>
> Taking a randomly selected starting point, that element from the first parent is copied into an 'offspring'. So say this algorithm chooses to start at index 3, 6 is copied into location 3 of an offspring.
>
> Offspring: - - 6 - - - - - - -
>
> The algorithm then looks at index location 3 in the second parent, and finds 1. It locates 1 in the first parent and copies that into the offspring, at its first parent index of 5.
>
> Offspring: - - 6 - 1 - - - - -
>
> It repeats, looks at index 5 in the second parent and finds 3. 3 is located in the first parent at index 8 and copied too the offspring.
>
> Offspring - - 6 - 1 - - 3 - -
>
> Now when the algorithm comes onto the next loop of the cycle, and checks index 8 of the second parent it finds 6 which is the value it started with. At this point the cycle halts and any empty spaces are filled with values from the second parent.
>
> Offspring: 2 5 6 7 1 8 4 3 10 9
>
> My entire function to achieve this:
>
> function offspring = cycle_crossover(firstParent, secondParent)
>
> % Initilisation of starting point and output matrix
> seedNumber = randi(length(firstParent), 1, 1);
> offspring = zeros(1, length(firstParent));
> offspring(seedNumber) = firstParent(seedNumber);
> currentNumber = seedNumber;
> cycleIndex = 0;
>
> % Perform swapping until cycle is complete
> [~, locations] = ismember(secondParent, firstParent);
> while cycleIndex ~= seedNumber;
> cycleIndex = locations(currentNumber);
> offspring(cycleIndex) = firstParent(cycleIndex);
> currentNumber = cycleIndex;
> end
>
> % Keep remaining items in their original location
> zeroLocations = ismember(offspring, 0);
> offspring(zeroLocations) = secondParent(zeroLocations);
> end
>
>
> I changed the 'copy remaining items' section code also based on what Bruno provided, and that too seemed to speed up significantly.
>
> For the scenario laid out, is that how you believe the code should look?
>
> Thanks again!
>
> Jonny
- - - - - - - - - -
  If you are sure that both firstParent and secondParent will be permutations on the same integers 1 to n, that would take care of the two difficulties I listed earlier. However, much more importantly than that, it should not be necessary to use ismember with its order n*log(n) algorithm in that case. By taking the inverse of one of these permutations applied to the other one you can easily provide the equivalent of the 'location' array that ismember would otherwise furnish, and inverses of permutations require only an order n algorithm. If n were large that could make an important difference in computation time required.

  I haven't worked out the details yet. If in general you will need to use firstParent and secondParent arrays that are not restricted in the above respect, please let me know before I go to all the trouble of figuring it out. My aging brain cells object to doing unnecessary work these days.

Roger Stafford

Subject: Vectorization Help...I think

From: Roger Stafford

Date: 18 Jul, 2010 06:34:13

Message: 7 of 9

"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in message <i1u2be$sdc$1@fred.mathworks.com>...
> "Jonny O'Connell" <jonnydotoconnelldotx@gmail.com> wrote in message <i1tpaj$kl0$1@fred.mathworks.com>...
> > > ......
> > > I think there is much more that you can tell us about the statistics of your problem that would be useful in giving you reliable advice.
> > > Roger Stafford
> >
> > Okay I'll try and put what I'm doing into a bit of context. This function takes two permutations of size n (size normally in the range 100-500) generated using the 'randperm' function, to be the firstParent and secondParent.
> >
> > Lets say we have the two parents:
> > First: 8 7 6 4 1 2 5 3 9 10
> > Second: 2 5 1 7 3 8 4 6 10 9
> >
> > Taking a randomly selected starting point, that element from the first parent is copied into an 'offspring'. So say this algorithm chooses to start at index 3, 6 is copied into location 3 of an offspring.
> >
> > Offspring: - - 6 - - - - - - -
> >
> > The algorithm then looks at index location 3 in the second parent, and finds 1. It locates 1 in the first parent and copies that into the offspring, at its first parent index of 5.
> >
> > Offspring: - - 6 - 1 - - - - -
> >
> > It repeats, looks at index 5 in the second parent and finds 3. 3 is located in the first parent at index 8 and copied too the offspring.
> >
> > Offspring - - 6 - 1 - - 3 - -
> >
> > Now when the algorithm comes onto the next loop of the cycle, and checks index 8 of the second parent it finds 6 which is the value it started with. At this point the cycle halts and any empty spaces are filled with values from the second parent.
> >
> > Offspring: 2 5 6 7 1 8 4 3 10 9
> >
> > My entire function to achieve this:
> >
> > function offspring = cycle_crossover(firstParent, secondParent)
> >
> > % Initilisation of starting point and output matrix
> > seedNumber = randi(length(firstParent), 1, 1);
> > offspring = zeros(1, length(firstParent));
> > offspring(seedNumber) = firstParent(seedNumber);
> > currentNumber = seedNumber;
> > cycleIndex = 0;
> >
> > % Perform swapping until cycle is complete
> > [~, locations] = ismember(secondParent, firstParent);
> > while cycleIndex ~= seedNumber;
> > cycleIndex = locations(currentNumber);
> > offspring(cycleIndex) = firstParent(cycleIndex);
> > currentNumber = cycleIndex;
> > end
> >
> > % Keep remaining items in their original location
> > zeroLocations = ismember(offspring, 0);
> > offspring(zeroLocations) = secondParent(zeroLocations);
> > end
> >
> >
> > I changed the 'copy remaining items' section code also based on what Bruno provided, and that too seemed to speed up significantly.
> >
> > For the scenario laid out, is that how you believe the code should look?
> >
> > Thanks again!
> >
> > Jonny
> - - - - - - - - - -
> If you are sure that both firstParent and secondParent will be permutations on the same integers 1 to n, that would take care of the two difficulties I listed earlier. However, much more importantly than that, it should not be necessary to use ismember with its order n*log(n) algorithm in that case. By taking the inverse of one of these permutations applied to the other one you can easily provide the equivalent of the 'location' array that ismember would otherwise furnish, and inverses of permutations require only an order n algorithm. If n were large that could make an important difference in computation time required.
>
> I haven't worked out the details yet. If in general you will need to use firstParent and secondParent arrays that are not restricted in the above respect, please let me know before I go to all the trouble of figuring it out. My aging brain cells object to doing unnecessary work these days.
>
> Roger Stafford
- - - - - - -
  It didn't strain my brain cells as much as I feared. In place of

 [~, locations] = ismember(secondParent, firstParent);

just put

 locations = 1:n;
 locations(firstParent) = locations; % <-- Gets inverse of firstParent
 locations = locations(secondParent);

and voila, you have the same 'locations' array as with 'ismember'. (This is of course assuming that 'firstParent' and 'secondParent' are both permutations on 1:n.) Each of these three lines is an order n algorithm, so for large n they should be faster than 'ismember'. The 'ismember' function is intended for more difficult problems where the two arguments are not merely the same set in different orderings.

  Also it unnecessary to replace the zeros in 'offspring' with original values from 'secondParent'. Just copy 'secondParent' over to 'offspring' in the beginning and your while loop will later simply replace those that are encountered with 'cycleIndex' while the rest remain as the original 'secondParent' values. No need to do 'ismember' here either. That gives you an overall order n algorithm.

Roger Stafford

Subject: Vectorization Help...I think

From: Bruno Luong

Date: 18 Jul, 2010 07:47:03

Message: 8 of 9

Completely agree with Roger.

If firstParent and secondParent are independent uniform permutations as generated with RANDPERM, then firstParent and locations are also independent uniform permutations. So further simplification would be generating firstParent and locations with RANDPERM, then obtaining secondParent by indexing:

n = 100;
firstParent = randperm(n);
locations = randperm(n);
secondParent = firstParent(locations);

% then pass firstParent and locations (instead of secondParent) on your algorithm.

Bruno

Subject: Vectorization Help...I think

From: Jonny O'Connell

Date: 18 Jul, 2010 14:46:06

Message: 9 of 9

Thanks very much to you both!

I've made the changes and they all work great! I can't unfortunately use the first and second parent generation you referred to Bruno as only on the first loop of my entire algorithm are they truly random, after that they evolve.

I'll be honest it was a bit confusing what you were both saying, but stepping through it with the debugger for a while and I got there in the end!

I appreciate you both going to the time to not only help me, but to also explain what you were doing. I'm a firm believer in teach a man to fish etc.. :)

Thanks again,

Jonny

Tags for this Thread

No tags are associated with 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