<div dir="ltr"><br><div class="gmail_extra"><br><div class="gmail_quote">On Wed, Oct 18, 2017 at 8:18 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:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);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 storage format.<br>
<br>
1)So when I was creating sequentially , I just used  MatCreateSeqAIJWithArrays , but the same for MPI version is quite confusing to use. I don't understand how to decide on the local rows(it would be really helpful if there is an example) .<br></blockquote><div><br></div><div>You don't use local row indices (you but you don't want to). The code does not change. As Jed says don't use <span style="font-size:12.8px">MatCreateSeqAIJWithArrays. Just use MatCreate. This is from ksp ex56.c:</span></div><div><span style="font-size:12.8px"><br></span></div><div><div><span style="font-size:12.8px">    /* create stiffness matrix */</span></div><div><span style="font-size:12.8px">    ierr = MatCreate(comm,&Amat);CHKERRQ(ierr);</span></div><div><span style="font-size:12.8px">    ierr = MatSetSizes(Amat,m,m,M,M);CHKERRQ(ierr);</span></div><div><span style="font-size:12.8px">    if (!test_late_bs) {</span></div><div><span style="font-size:12.8px">      ierr = MatSetBlockSize(Amat,3);CHKERRQ(ierr);</span></div><div><span style="font-size:12.8px">    }</span></div><div><span style="font-size:12.8px">    ierr = MatSetType(Amat,MATAIJ);CHKERRQ(ierr);</span></div><div><span style="font-size:12.8px">    ierr = MatSeqAIJSetPreallocation(Amat,0,d_nnz);CHKERRQ(ierr);</span></div><div><span style="font-size:12.8px">    ierr = MatMPIAIJSetPreallocation(Amat,0,d_nnz,0,o_nnz);CHKERRQ(ierr);</span></div></div><div><span style="font-size:12.8px"><br></span></div><div><span style="font-size:12.8px">You can preallocate with an estimate (upper bound). See the documentation for  </span><span style="font-size:12.8px">MatMPIAIJSetPreallocation (Google it).</span></div><div><span style="font-size:12.8px"><br></span></div><div><span style="font-size:12.8px">You then run your code on one processor and PETSc will distribute it.  Now you just add values with (i,j,value) with MatSetValues (Google it). You will find that it is very simple.</span></div><div><span style="font-size:12.8px"><br></span></div><div><span style="font-size:12.8px">Now, this simple way will just chop your domain up in a dumb way. If you have a regular grid then you will get a 1D partitioning, which will work for a while. Otherwise you can partition this matrix. Bat that is another story. You want to start with this simple way anyway.</span></div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<br>
2)When I also tried using MatSetValues it doesn't seem to use the same indexing as compressed row storage format.What type of indexing should be used when MatSetValues are used and called from rank 0 for CRS Matrices?<div class="gmail-HOEnZb"><div class="gmail-h5"><br>
<br>
On 2017-10-18 13:33, Jed Brown wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);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>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
Hello,<br>
<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>
</blockquote></blockquote>
</div></div></blockquote></div><br></div></div>