<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Sun, Apr 27, 2014 at 7:47 AM, Justin Dong <span dir="ltr"><<a href="mailto:jsd1@rice.edu" target="_blank">jsd1@rice.edu</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div dir="ltr">I think I've got the right idea on how to divide up the for loop:<div><br></div><div>    MPI_Comm_size(MPI_COMM_WORLD,&numprocs); </div><div>    MPI_Comm_rank(MPI_COMM_WORLD,&myid); </div><div><br>

</div><div><div>    mystart = (nelements / numprocs) * myid;</div><div>    if (nelements % numprocs > myid)</div><div>    {</div><div>    <span style="white-space:pre-wrap">  </span>mystart += myid;</div><div>    <span style="white-space:pre-wrap">   </span>myend = mystart + (nelements / numprocs) + 1;</div>

<div>    }</div><div>  <span style="white-space:pre-wrap">  </span>else</div><div>  <span style="white-space:pre-wrap">  </span>{</div><div>    <span style="white-space:pre-wrap">      </span>mystart += nelements % numprocs;</div>
<div>
    <span style="white-space:pre-wrap">       </span>myend = mystart + (nelements / numprocs);</div><div>  <span style="white-space:pre-wrap">  </span>}</div></div></div></blockquote><div><br></div><div>We have a function that does exactly this: <a href="http://www.mcs.anl.gov/petsc/petsc-dev/docs/manualpages/Sys/PetscSplitOwnership.html">http://www.mcs.anl.gov/petsc/petsc-dev/docs/manualpages/Sys/PetscSplitOwnership.html</a></div>
<div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div>But then, do I do this?</div>
</div></blockquote><div><br></div><div>Yes.</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div dir="ltr">
<div><div style="font-family:arial,sans-serif;font-size:13px">for (k = mystart+1; k < myend; ++k)</div><div style="font-family:arial,sans-serif;font-size:13px">{</div><div style="font-family:arial,sans-serif;font-size:13px">

    get_index(k,ie,je); /* ie and je are arrays that indicate where to place A_local */</div><div style="font-family:arial,sans-serif;font-size:13px">    compute_local_matrix(A_local,...);</div><div style="font-family:arial,sans-serif;font-size:13px">

<br></div><div style="font-family:arial,sans-serif;font-size:13px">    MatSetValues(A_global, nlocal, ie, nlocal, je, A_local, ADD_VALUES);</div><div style="font-family:arial,sans-serif;font-size:13px">}</div></div><div>
<br>
</div><div>The indices are global indices and I'm doing </div><div><br></div><div><div>MatCreateSeqAIJ(PETSC_COMM_WORLD, nglobal, nglobal, 3*nlocal, </div><div><span style="white-space:pre-wrap">                           </span>    PETSC_NULL, &A_global);</div>
</div></div></blockquote><div><br></div><div>This should be MatCreateMPIAIJ() for parallel execution.</div><div><br></div><div>   Matt</div><div><br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div dir="ltr"><div>to create the matrix. Running the program seems to give me multiple errors, but mainly</div><div>







<p>[3]PETSC ERROR: Object is in wrong state!</p>
<p>[3]PETSC ERROR: Mat object's type is not set: Argument # 1!</p><p><br></p></div><div><br></div></div><div class="gmail_extra"><br><br><div class="gmail_quote">On Sun, Apr 27, 2014 at 6:54 AM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br>

<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">
<div>On Sun, Apr 27, 2014 at 6:25 AM, Justin Dong <span dir="ltr"><<a href="mailto:jsd1@rice.edu" target="_blank">jsd1@rice.edu</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div dir="ltr">Hi Matt,<div><br></div><div>For option 1), that would be using MPI and not any special functions in PETSc? I ask since I've never done parallel programming before. I tried using OpenMP to parallelize that for loop but it resulted in some errors and I didn't investigate further, but I'm assuming it's because each process wasn't communicating properly with the others during MatSetValues?</div>


</div></blockquote><div><br></div></div><div>Yes, using MPI, so each process owns a set of elements that it loops over. The Mat object manages the global</div><div>values as long as you use global indices for the (row, column). There are of course many refinements for this,</div>


<div>but I think the thing to do is get something working fast, and then optimize it.</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</div><div><div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">


<div dir="ltr"><div>Sincerely,</div><div>Justin</div><div class="gmail_extra"><br><br><div class="gmail_quote">On Sun, Apr 27, 2014 at 5:57 AM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br>



<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">
<div><div>On Sun, Apr 27, 2014 at 4:14 AM, Justin Dong <span dir="ltr"><<a href="mailto:jsd1@rice.edu" target="_blank">jsd1@rice.edu</a>></span> wrote:<br>


<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div dir="ltr">Hi all,<div><br></div><div>I'm currently coding a finite element problem in PETSc. I'm computing all of the matrices by myself and would prefer to do it that way. I want to parallelize the assembly of the global matrices. This is a simplified version of the assembly routine (pseudocode):</div>





<div><br></div><div>for (k = 0; k < nelements; ++k)</div><div>{</div><div>    get_index(k,ie,je); /* ie and je are arrays that indicate where to place A_local */</div><div>    compute_local_matrix(A_local,...);</div><div>





<br></div><div>    MatSetValues(A_global, nlocal, ie, nlocal, je, A_local, ADD_VALUES);</div><div>}</div><div><br></div><div>This is for DG finite elements and I know the matrix easily be assembled in parallel. Even on my laptop, this would allow me to significantly speed up my code. The only problem is, I'm not sure how to do this in PETSc. I'd assume I need to use MatCreateMPIAIJ instead of MatCreateSeqAIJ, but wasn't able to find any using tutorials on how I might do this. </div>




</div></blockquote><div><br></div></div></div><div>1) If you just split this loop across processes, it would work immediately. However, that can be non-optimal in terms</div><div>     of communication.</div><div><br></div>



<div>2) PETSc matrices are distributed by contiguous blocks of rows (see manual section on matrices), so you would like</div>
<div>    those row blocks to correspond roughly to your element blocks.</div><div><br></div><div>3) You will also have to preallocate, but that should be the same as you do now for the sequential case, except you</div><div>




     check whether a column is inside the diagonal block.</div><div><br></div><div>  Thanks,</div><div><br></div><div>      Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">




<div dir="ltr"><div>Sincerely,</div><div>Justin</div></div><span><font color="#888888">
</font></span></blockquote></div><span><font color="#888888"><br><br clear="all"><span class=""><font color="#888888"><span><font color="#888888"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>




-- Norbert Wiener
</font></span></font></span></font></span></div></div><span class=""><font color="#888888">
</font></span></blockquote></div><span class=""><font color="#888888"><br></font></span></div></div><span class=""><font color="#888888">
</font></span></blockquote></div></div></div><span class=""><font color="#888888"><div><div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>


-- Norbert Wiener
</div></div></font></span></div></div>
</blockquote></div><br></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener
</div></div>