<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Thu, Oct 19, 2017 at 10:52 AM, Jaganathan, Srikrishna <span dir="ltr"><<a href="mailto:srikrishna.jaganathan@fau.de" target="_blank">srikrishna.jaganathan@fau.de</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Thanks a lot for your detailed reply, it was very slow to insert values using MatSetValues , I guess its due to lack of preallocation. I still couldn't figure out how to calculate the parameters required to preallocate MatMPIAIJ. So I have the nnz array for serial preallocation and this is accurate. So is there any rule of thumb to arrive at decent numbers for d_nnz and o_nnz from nnz array ?<br></blockquote><div><br></div><div>Yes. The d_nnz array wants nonzero that occur in the diagonal block. Thus, if your process owns rows [rStart, rEnd) then</div><div>d_nnz is the number of nonzeros in a row that lie in columns [rStart, rEnd). The rest of the nonzeros in that row are counted</div><div>in o_nnz.</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
On 2017-10-19 14:41, Mark Adams wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
On Wed, Oct 18, 2017 at 8:18 AM, Jaganathan, Srikrishna <<br>
<a href="mailto:srikrishna.jaganathan@fau.de" target="_blank">srikrishna.jaganathan@fau.de</a>> wrote:<br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Thanks for your response, its helpful.<br>
<br>
I do have few more questions, most of my matrices are of compressed row<br>
storage format.<br>
<br>
1)So when I was creating sequentially , I just used<br>
MatCreateSeqAIJWithArrays , but the same for MPI version is quite confusing<br>
to use. I don't understand how to decide on the local rows(it would be<br>
really helpful if there is an example) .<br>
<br>
</blockquote>
<br>
You don't use local row indices (you but you don't want to). The code does<br>
not change. As Jed says don't use MatCreateSeqAIJWithArrays. Just use<br>
MatCreate. This is from ksp ex56.c:<br>
<br>
    /* create stiffness matrix */<br>
    ierr = MatCreate(comm,&Amat);CHKERRQ(<wbr>ierr);<br>
    ierr = MatSetSizes(Amat,m,m,M,M);CHKE<wbr>RRQ(ierr);<br>
    if (!test_late_bs) {<br>
      ierr = MatSetBlockSize(Amat,3);CHKERR<wbr>Q(ierr);<br>
    }<br>
    ierr = MatSetType(Amat,MATAIJ);CHKERR<wbr>Q(ierr);<br>
    ierr = MatSeqAIJSetPreallocation(Amat<wbr>,0,d_nnz);CHKERRQ(ierr);<br>
    ierr = MatMPIAIJSetPreallocation(Amat<wbr>,0,d_nnz,0,o_nnz);CHKERRQ(ierr<wbr>);<br>
<br>
You can preallocate with an estimate (upper bound). See the documentation<br>
for  MatMPIAIJSetPreallocation (Google it).<br>
<br>
You then run your code on one processor and PETSc will distribute it.  Now<br>
you just add values with (i,j,value) with MatSetValues (Google it). You<br>
will find that it is very simple.<br>
<br>
Now, this simple way will just chop your domain up in a dumb way. If you<br>
have a regular grid then you will get a 1D partitioning, which will work<br>
for a while. Otherwise you can partition this matrix. Bat that is another<br>
story. You want to start with this simple way anyway.<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
2)When I also tried using MatSetValues it doesn't seem to use the same<br>
indexing as compressed row storage format.What type of indexing should be<br>
used when MatSetValues are used and called from rank 0 for CRS Matrices?<br>
<br>
<br>
On 2017-10-18 13:33, Jed Brown wrote:<br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Easiest is to assemble into a distributed matrix from rank 0.  So<br>
instead of calling MatCreate using PETSC_COMM_SELF, use a parallel<br>
communicator (like PETSC_COMM_WORLD).  It is fine if only rank 0 calls<br>
MatSetValues, but all processes must call MatAssemblyBegin/End.<br>
<br>
"Jaganathan, Srikrishna" <<a href="mailto:srikrishna.jaganathan@fau.de" target="_blank">srikrishna.jaganathan@fau.de</a>> writes:<br>
<br>
Hello,<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
<br>
I have been trying to distribute a already existing stiffness matrix in<br>
my FEM code to petsc parallel matrix object , but I am unable to find<br>
any documentation regarding it. It was quite straightforward to create a<br>
sequential petsc matrix object and everything was working as intended.I<br>
have read some of the user comments in the mailing lists regarding<br>
similar situation and most of the times the solution suggested is to<br>
create stiffness matrix from the the mesh in distributed format. Since<br>
its a little difficult in my case to pass the mesh data in the code , is<br>
there anyway to distribute already existing stiffness matrix ?<br>
<br>
Thanks and Regards<br>
<br>
Srikrishna Jaganathan<br>
<br>
</blockquote>
<br>
</blockquote></blockquote></blockquote>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div>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><br></div><div><a href="http://www.caam.rice.edu/~mk51/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div>
</div></div>