Thread Subject:

parfor help

Hi, I have this problem

pa and pm are integer matrices, that I use as indexes to other
matrices
c is an 8X8 that is different in every loop


for i=1:100
    for j=1:80
        ia=pa(j); im=pm(j);
        sum2=0;
        for n=1:num_of_elements(i)
            el=jtmp(n,i);
            c=kernel(:,:,el);
            tmp=cm.icon(1:8,el);
            a1=pot(ia,tmp);
            a2=pot(im,tmp);
            S2=((conj(a1))'*a2).*c;
            jam=sum(S2(:));
            sum2=sum2+jam
        end
        jac(i,j)=sum2;
    end
end


The first for loop is independent from the other's so I want to
transform it to
parfor i=1:100

but I can't store the value of sum2 to jac, because is is not valid
indexes are restricted.


How can i parallelize this?

Marios Karaoulis <marios.karaoulis@gmail.com> wrote in message <e194a14d-4ff3-46d5-892e-5d58581d55dc@s24g2000pri.googlegroups.com>...
> Hi, I have this problem
>
> pa and pm are integer matrices, that I use as indexes to other
> matrices
> c is an 8X8 that is different in every loop
>
>
> for i=1:100
> for j=1:80
> ia=pa(j); im=pm(j);
> sum2=0;
> for n=1:num_of_elements(i)
> el=jtmp(n,i);
> c=kernel(:,:,el);
> tmp=cm.icon(1:8,el);
> a1=pot(ia,tmp);
> a2=pot(im,tmp);
> S2=((conj(a1))'*a2).*c;
> jam=sum(S2(:));
> sum2=sum2+jam
> end
> jac(i,j)=sum2;
> end
> end
>
>
> The first for loop is independent from the other's so I want to
> transform it to
> parfor i=1:100
>
> but I can't store the value of sum2 to jac, because is is not valid
> indexes are restricted.
>
>
> How can i parallelize this?

Give this a try:
%%%%
[ii, jj] = meshgrid(1:100,1:80); %come up with all indices ahead of time
jac = zeros(100,80); %preallocate: listen to M-LINT!!!
for idx = 1:numel(ii)
    ia=pa(jj(idx)); im=pm(jj(idx));
    sum2=0;
    jtmp_this_iter = jtmp(:,ii(idx));
    parfor n=1:num_of_elements(ii(idx))
        el=jtmp_this_iter(n);
        c=kernel(:,:,el); %I don't think there's a way to slice Kernel. At least not with the minimal amount of data you've given us. There may be a better way to run this whole thing if we had sample matrices.
        tmp=cm.icon(1:8,el);
        a1=pot(ia,tmp);
        a2=pot(im,tmp);
        S2=((conj(a1))'*a2).*c;
        jam=sum(S2(:));
        sum2=sum2+jam;
    end
    jac(ii(idx),jj(idx))=sum2;
end

With the knowledge you've provided the parallelization is best done inside the regular for-loop. There might be another way around this. Okay there is another way around this; but we'd have to know more about your matrices: What we would need is two more matrices with direct values i.e. not num_of_elements(ii(idx)) but rather a list of those elements that could then be associated with the single loop (i.e. the idx refers to the number directly. jtmp would also have to be indexable this way.

The result would be this:
%single parfor loop
parfor idx = 1:numel(ii)
    ia=pa(jj(idx)); im=pm(jj(idx));
    sum2=0;
    el = jtmp(idx);
    c=kernel(:,:,el);
    tmp=cm.icon(1:8,el);
    a1=pot(ia,tmp);
    a2=pot(im,tmp);
    S2=((conj(a1))'*a2).*c;
    jam=sum(S2(:));
    sum2=sum2+jam;
    jac(ii(idx),jj(idx))=sum2;
end

%This shouldn't be too hard to accomplish, Good Luck!
-Sean

And why not do parfor on j=1:80 ?

"Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i6re2h$bq0$1@fred.mathworks.com>...
> And why not do parfor on j=1:80 ?

That would work just as well as the first example most likely. However, depending on how big 'kernel' is there could be quite a bit of wasted overhead since it won't be sliced.

On Sep 15, 5:33 pm, "Matt J " <mattjacREM...@THISieee.spam> wrote:
> And why not do parfor on j=1:80  ?

I tried this, it works, but the speed is way slower.

The same loor with parfor took 300 seconds,
with just for took 150 sec.

So I need another technique...

On Sep 15, 3:30 pm, Marios Karaoulis <marios.karaou...@gmail.com>
wrote:
> Hi, I have this problem
>
> pa and pm are integer matrices, that I use as indexes to other
> matrices
> c is an 8X8 that is different in every loop
>
> for i=1:100
>     for j=1:80
>         ia=pa(j); im=pm(j);
>         sum2=0;
>         for n=1:num_of_elements(i)
>             el=jtmp(n,i);
>             c=kernel(:,:,el);
>             tmp=cm.icon(1:8,el);
>             a1=pot(ia,tmp);
>             a2=pot(im,tmp);
>             S2=((conj(a1))'*a2).*c;
>             jam=sum(S2(:));
>             sum2=sum2+jam
>         end
>         jac(i,j)=sum2;
>     end
> end
>
> The first for loop is independent from the other's so I want to
>



I will try to describe what are each matrix.


We are talking for the finite element method.


for every element (el), i have already calculate 8x8 matrix c. So I
sored is to 3d matrix (kernel), where in each 3rd dimension, there is
this 8x8 matrix.

for every loop, I need to take some elements from the pot matrix, in
the indexes defined in matrices pa,pm and tmp.

On Sep 15, 4:57 pm, "Sean "
<sean.dewol...@nospamplease.umit.maine.edu> wrote:
> Marios Karaoulis <marios.karaou...@gmail.com> wrote in message <e194a14d-4ff3-46d5-892e-5d58581d5...@s24g2000pri.googlegroups.com>...
> > Hi, I have this problem
>
> > pa and pm are integer matrices, that I use as indexes to other
> > matrices
> > c is an 8X8 that is different in every loop
>
> > for i=1:100
> >     for j=1:80
> >         ia=pa(j); im=pm(j);
> >         sum2=0;
> >         for n=1:num_of_elements(i)
> >             el=jtmp(n,i);
> >             c=kernel(:,:,el);
> >             tmp=cm.icon(1:8,el);
> >             a1=pot(ia,tmp);
> >             a2=pot(im,tmp);
> >             S2=((conj(a1))'*a2).*c;
> >             jam=sum(S2(:));
> >             sum2=sum2+jam
> >         end
> >         jac(i,j)=sum2;
> >     end
> > end
>
> > The first for loop is independent from the other's so I want to
> > transform it to
> > parfor i=1:100
>
> > but I can't store the value of sum2 to  jac, because is is not valid
> > indexes are restricted.
>
> > How can i parallelize this?
>
> Give this a try:
> %%%%
> [ii, jj] = meshgrid(1:100,1:80); %come up with all indices ahead of time
> jac = zeros(100,80); %preallocate: listen to M-LINT!!!
> for idx = 1:numel(ii)
>     ia=pa(jj(idx)); im=pm(jj(idx));
>     sum2=0;
>     jtmp_this_iter = jtmp(:,ii(idx));
>     parfor n=1:num_of_elements(ii(idx))
>         el=jtmp_this_iter(n);
>         c=kernel(:,:,el); %I don't think there's a way to slice Kernel.  At least not with the minimal amount of data you've given us.  There may be a better way to run this whole thing if we had sample matrices.
>         tmp=cm.icon(1:8,el);
>         a1=pot(ia,tmp);
>         a2=pot(im,tmp);
>         S2=((conj(a1))'*a2).*c;
>         jam=sum(S2(:));
>         sum2=sum2+jam;
>     end
>     jac(ii(idx),jj(idx))=sum2;
> end
>


By the way, by using this without the parfor loop, reduced the
computation speed for almost 40%.
Cool.


But with the parfor loop, as you mentioned it, it takes for ever...

Marios Karaoulis <marios.karaoulis@gmail.com> wrote in message <b464d524-e1de-48ed-ac04-86dc06e1e153@x24g2000pro.googlegroups.com>...
> On Sep 15, 4:57 pm, "Sean "
> <sean.dewol...@nospamplease.umit.maine.edu> wrote:
> > Marios Karaoulis <marios.karaou...@gmail.com> wrote in message <e194a14d-4ff3-46d5-892e-5d58581d5...@s24g2000pri.googlegroups.com>...
> > > Hi, I have this problem
> >
> > > pa and pm are integer matrices, that I use as indexes to other
> > > matrices
> > > c is an 8X8 that is different in every loop
> >
> > > for i=1:100
> > >     for j=1:80
> > >         ia=pa(j); im=pm(j);
> > >         sum2=0;
> > >         for n=1:num_of_elements(i)
> > >             el=jtmp(n,i);
> > >             c=kernel(:,:,el);
> > >             tmp=cm.icon(1:8,el);
> > >             a1=pot(ia,tmp);
> > >             a2=pot(im,tmp);
> > >             S2=((conj(a1))'*a2).*c;
> > >             jam=sum(S2(:));
> > >             sum2=sum2+jam
> > >         end
> > >         jac(i,j)=sum2;
> > >     end
> > > end
> >
> > > The first for loop is independent from the other's so I want to
> > > transform it to
> > > parfor i=1:100
> >
> > > but I can't store the value of sum2 to  jac, because is is not valid
> > > indexes are restricted.
> >
> > > How can i parallelize this?
> >
> > Give this a try:
> > %%%%
> > [ii, jj] = meshgrid(1:100,1:80); %come up with all indices ahead of time
> > jac = zeros(100,80); %preallocate: listen to M-LINT!!!
> > for idx = 1:numel(ii)
> >     ia=pa(jj(idx)); im=pm(jj(idx));
> >     sum2=0;
> >     jtmp_this_iter = jtmp(:,ii(idx));
> >     parfor n=1:num_of_elements(ii(idx))
> >         el=jtmp_this_iter(n);
> >         c=kernel(:,:,el); %I don't think there's a way to slice Kernel.  At least not with the minimal amount of data you've given us.  There may be a better way to run this whole thing if we had sample matrices.
> >         tmp=cm.icon(1:8,el);
> >         a1=pot(ia,tmp);
> >         a2=pot(im,tmp);
> >         S2=((conj(a1))'*a2).*c;
> >         jam=sum(S2(:));
> >         sum2=sum2+jam;
> >     end
> >     jac(ii(idx),jj(idx))=sum2;
> > end
> >
>
>
> By the way, by using this without the parfor loop, reduced the
> computation speed for almost 40%.
> Cool.
>
>
> But with the parfor loop, as you mentioned it, it takes for ever...
>

First: Are the two results equal?
Second: How many workers do you have open for parfor? If you only have one worker it'll be much slower. I would recommend closing your old workers before starting this job and then reopening them (just to be sure!). (matlabpool close force;matlabpool open N)

I think you can avoid the inner loop entirely if you restructure (combination of reshape/replicate) your pot, jtmp, and (maybe) cm.icon to be 3-dimensional. Thus each slice corresponds to one element.
You could then S2 and sum2 the same way.

On Sep 17, 9:42 am, "Sean "
<sean.dewol...@nospamplease.umit.maine.edu> wrote:
> Marios Karaoulis <marios.karaou...@gmail.com> wrote in message <b464d524-e1de-48ed-ac04-86dc06e1e...@x24g2000pro.googlegroups.com>...
> > On Sep 15, 4:57 pm, "Sean "
> > <sean.dewol...@nospamplease.umit.maine.edu> wrote:
> > > Marios Karaoulis <marios.karaou...@gmail.com> wrote in message <e194a14d-4ff3-46d5-892e-5d58581d5...@s24g2000pri.googlegroups.com>...
> > > > Hi, I have this problem
>
> > > > pa and pm are integer matrices, that I use as indexes to other
> > > > matrices
> > > > c is an 8X8 that is different in every loop
>
> > > > for i=1:100
> > > >     for j=1:80
> > > >         ia=pa(j); im=pm(j);
> > > >         sum2=0;
> > > >         for n=1:num_of_elements(i)
> > > >             el=jtmp(n,i);
> > > >             c=kernel(:,:,el);
> > > >             tmp=cm.icon(1:8,el);
> > > >             a1=pot(ia,tmp);
> > > >             a2=pot(im,tmp);
> > > >             S2=((conj(a1))'*a2).*c;
> > > >             jam=sum(S2(:));
> > > >             sum2=sum2+jam
> > > >         end
> > > >         jac(i,j)=sum2;
> > > >     end
> > > > end
>
> > > > The first for loop is independent from the other's so I want to
> > > > transform it to
> > > > parfor i=1:100
>
> > > > but I can't store the value of sum2 to  jac, because is is not valid
> > > > indexes are restricted.
>
> > > > How can i parallelize this?
>
> > > Give this a try:
> > > %%%%
> > > [ii, jj] = meshgrid(1:100,1:80); %come up with all indices ahead of time
> > > jac = zeros(100,80); %preallocate: listen to M-LINT!!!
> > > for idx = 1:numel(ii)
> > >     ia=pa(jj(idx)); im=pm(jj(idx));
> > >     sum2=0;
> > >     jtmp_this_iter = jtmp(:,ii(idx));
> > >     parfor n=1:num_of_elements(ii(idx))
> > >         el=jtmp_this_iter(n);
> > >         c=kernel(:,:,el); %I don't think there's a way to slice Kernel.  At least not with the minimal amount of data you've given us.  There may be a better way to run this whole thing if we had sample matrices.
> > >         tmp=cm.icon(1:8,el);
> > >         a1=pot(ia,tmp);
> > >         a2=pot(im,tmp);
> > >         S2=((conj(a1))'*a2).*c;
> > >         jam=sum(S2(:));
> > >         sum2=sum2+jam;
> > >     end
> > >     jac(ii(idx),jj(idx))=sum2;
> > > end
>
> > By the way, by using this without the parfor loop, reduced the
> > computation speed for almost 40%.
> > Cool.
>
> > But with the parfor loop, as you mentioned it, it takes for ever...
>
> First:  Are the two results equal?
> Second: How many workers do you have open for parfor?  If you only have one worker it'll be much slower.  I would recommend closing your old workers before starting this job and then reopening them (just to be sure!).  (matlabpool close force;matlabpool open N)
>
> I think you can avoid the inner loop entirely if you restructure (combination of reshape/replicate) your pot, jtmp, and (maybe) cm.icon to be 3-dimensional.  Thus each slice corresponds to one element.
> You could then S2 and sum2 the same way.



Yes, results using the first algorithm, are the same. I haven't
succeeded with the second algorithm.
I have a core i7, and matlab sees 8 workers.

I will try again with your idea, but any help would be appreciate it.

Thanks.

P.S. By the way, can I send you the matrices?

Marios Karaoulis <marios.karaoulis@gmail.com> wrote in message <3e05ef26-b03b-
> P.S. By the way, can I send you the matrices?

Go for it.
I'm coincidentally taking an FEA course right now.

On Sep 17, 10:45 am, "Sean "
<sean.dewol...@nospamplease.umit.maine.edu> wrote:
> Marios Karaoulis <marios.karaou...@gmail.com> wrote in message <3e05ef26-b03b-
> > P.S. By the way, can I send you the matrices?
>
> Go for it.
> I'm coincidentally taking an FEA course right now.

How can I send it to you? I do not have your email