Question about multithreading a function

Hello, I have a question about a function that by itself, creates a lot of overhead.

The function in question generates a point cloud in a plane defined by three points. The point cloud is a 51 x 51 point cloud, meaning I have 51 points along one direction vector for an edge of the plane, then I have 51 points connecting every of the 51 points on the direction vector. I attached a picture of what this looks like with 11 x 11 points.

Now, what I do that takes up a lot of time is connect each of the 51 x 51 points to all of the three points, making a binary tree that bifurcates. I have the three points that define the edges of the plane, then the 51 x 51 points that make up the point of bifurcation.

The part that adds the most overhead is this right here:

    for i=1:length(newPointCloud)
        subtreePtMx(:,:,i) = [ccoPoints; randomxyz; newPointCloud(i,:)];
        subPropMx(:,1,i) = sqrt(sum((bsxfun(@minus, subtreePtMx(subFaceMx(:,1),:,i),subtreePtMx(subFaceMx(:,2),:,i))).^2, 2));
        subPropMx(:,2,i) = [1; (1+((subPropMx(2,1,i)/subPropMx(3,1,i))^(1/4))^(-3))^(-1/3); (1+((subPropMx(2,1,i)/subPropMx(3,1,i))^(1/4))^(3))^(-1/3)];
        subVolumes(i) = pi*(subPropMx(1,1,i)+subPropMx(2,1,i)*subPropMx(2,2,i)^2+subPropMx(3,1,i)*subPropMx(3,2,i)^2);
    end

In the first part, I create a (4,3,N) matrix where I have the first three points, and add one point in the point cloud as the forth point. I then in the next line calculate the distance between the points for each for each (4,3,i) sub matrix, and this takes roughly 70% of the time for the function. The next couple of parts do some math calculating the eventual volume of the subtree, I am trying to find the min volume.

Now, here are my questions: is it possible to change this line:

        subPropMx(:,1,i) = sqrt(sum((bsxfun(@minus, subtreePtMx(subFaceMx(:,1),:,i),subtreePtMx(subFaceMx(:,2),:,i))).^2, 2));

to make it less computationally expensive? Its basically just a distance function of two points, done three times iteratively. Barring that, is it possible to have this function run in parallel using matlabpool and parfor?

Thanks, any advice would be greatly appreciated!

Answer 1

[Matt J]

I don't quite understand your current implementation, but the pattern in the image you've posted should be very simple to compute using a few builtin function calls (which, of course, will be very well multi-threaded internally). Let A be a 2xN matrix whose columns are the (x,y) points along edge #1, and B the opposing points along edge #2 so that A(:,i) and B(:,i) correspond. Then the set of all points can be generated as follows

   n=11; %or n=51

   c=reshape(linspace(0,1,n),1,1,[]);

   AllPoints=bsxfun(@plus ,  bsxfun(@times,c,B-A),   A);

The result should be that AllPoints(:,i,j) will be the j-th point connecting edge points A(:,i) and B(:,i)