<div dir="ltr">The problem isn't only with 1 process, but it seems to be with all non-square matrices?<div><br></div><div>For example, here I have a two process program which tries to set a 6x7 matrix to (each process owns 3 rows):</div><div><br></div><div>0 0 0 1 1 1 1</div><div>0 0 0 1 1 1 1</div><div>0 0 0 1 1 1 1</div><div><div>0 0 0 1 1 1 1</div><div>0 0 0 1 1 1 1</div><div>0 0 0 1 1 1 1</div></div><div><div><br></div><div><div>#include <petsc.h></div><div>#include <mpi.h></div><div><br></div><div>int main(int argc, char** argv)</div><div>{</div><div> PetscErrorCode err;</div><div> err = PetscInitialize(&argc, &argv, NULL, "help");</div><div> CHKERRQ(err);</div><div><br></div><div> int rank, size;</div><div> MPI_Comm_rank(MPI_COMM_WORLD, &rank);</div><div> MPI_Comm_size(MPI_COMM_WORLD, &size);</div><div> if(size != 2)</div><div> {</div><div> printf("must run with 2 processes");</div><div> MPI_Abort(MPI_COMM_WORLD, -1);</div><div> }</div><div><br></div><div> // create a sparse AIJ matrix distributed across MPI</div><div> Mat A;</div><div> err = MatCreate(MPI_COMM_WORLD, &A);</div><div> CHKERRQ(err);</div><div> err = MatSetType(A, MATMPIAIJ);</div><div> CHKERRQ(err);</div><div> // setup pre-allocation for matrix space</div><div> {</div><div> err =</div><div> MatSetSizes(A, 3, PETSC_DECIDE, 6, 7);</div><div> CHKERRQ(err);</div><div> if(rank == 0)</div><div> {</div><div> PetscInt d_nnz[] = {0, 0, 0};</div><div> PetscInt o_nnz[] = {4, 4, 4};</div><div> err = MatMPIAIJSetPreallocation(A, 0, d_nnz, 0, o_nnz);</div><div> CHKERRQ(err);</div><div> }</div><div> else</div><div> {</div><div> PetscInt d_nnz[] = {3, 3, 3};</div><div> PetscInt o_nnz[] = {1, 1, 1};</div><div> err = MatMPIAIJSetPreallocation(A, 0, d_nnz, 0, o_nnz);</div><div> CHKERRQ(err);</div><div> }</div><div> }</div><div> err = MatSetUp(A);</div><div> CHKERRQ(err);</div><div><br></div><div> // set values inside the matrix</div><div> for (PetscInt row = 0; row < 3; ++row)</div><div> {</div><div> for (PetscInt col = 3; col < 7; ++col)</div><div> {</div><div> err = MatSetValue(A, 3 * rank + row, col, 1, INSERT_VALUES);</div><div> CHKERRQ(err);</div><div> }</div><div> }</div><div><br></div><div> err = MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);</div><div> CHKERRQ(err);</div><div> err = MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);</div><div> CHKERRQ(err);</div><div><br></div><div> err = MatView(A, PETSC_VIEWER_STDOUT_WORLD);</div><div> CHKERRQ(err);</div><div><br></div><div> // free memory</div><div> err = MatDestroy(&A);</div><div> CHKERRQ(err);</div><div><br></div><div> // cleanup any internal PETSc data at end of program</div><div> err = PetscFinalize();</div><div> CHKERRQ(err);</div><div>}</div></div><div><br></div></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Tue, Feb 14, 2017 at 5:26 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:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><span class="">On Tue, Feb 14, 2017 at 6:41 AM, Andrew Ho <span dir="ltr"><<a href="mailto:andrewh0@uw.edu" target="_blank">andrewh0@uw.edu</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div>I have a 3x6 matrix, and I want to set it to (just as an example):</div><div><br></div><div>0 0 0 1 1 1</div><div>0 0 0 1 1 1</div><div>0 0 0 1 1 1<br><br>From what I can tell, according to the documentation the MPIAIJ sparsity of this matrix is:</div></div></blockquote><div><br></div></span><div>I believe this is an inconsistency in PETSc for rectangular matrices. We divide matrices into diagonal and</div><div>off-diagonal parts mainly to separate communication from computation. Thus on 1 proc, we never divide</div><div>them, even if the matrix is rectangular. Therefore we misinterpret your preallocation.</div><div><br></div><div>Barry, do you think we want to try and change this?</div><div><br></div><div> Matt</div><div><div class="h5"><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div>d_nnz = [0, 0, 0]</div><div>o_nnz = [3, 3, 3]</div><div><br></div><div>However, when I do this, I get the following error:</div><div><br></div><div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">[0]PETSC ERROR: --------------------- Error Message ------------------------------<wbr>------------------------------<wbr>-- <br>[0]PETSC ERROR: Argument out of range<br>[0]PETSC ERROR: New nonzero at (0,3) caused a malloc<br>Use MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR<wbr>, PETSC_FALSE) to turn off this check<br>[0]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/d<wbr>ocumentation/faq.html</a> for trouble shooting.<br>[0]PETSC ERROR: Petsc Release Version 3.7.5, unknown <br>[0]PETSC ERROR: ./ex3 on a arch-linux2-c-opt Tue Feb 14 04:23:56 2017<br>[0]PETSC ERROR: Configure options --with-debugging=0 --COPTFLAGS="-O3 -march=native" --CXXOPTFLAGS="-O3 -march=native" --FOPTFLAGS="-O3 -march=native" <br>[0]PETSC ERROR: #1 MatSetValues_MPIAIJ() line 582 in petsc/src/mat/impls/aij/mpi/mp<wbr>iaij.c<br>[0]PETSC ERROR: #2 MatSetValues() line 1190 in petsc/src/mat/interface/matrix<wbr>.c <br>[0]PETSC ERROR: #3 main() line 36 in ex3.c<br>[0]PETSC ERROR: No PETSc Option Table entries<br>[0]PETSC ERROR: ----------------End of Error Message -------send entire error message to petsc-maint@mcs.anl.gov-------<wbr>---</blockquote><div><br></div><div>Here's a working test code:</div><div><br></div><div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">#include <petsc.h><br>#include <mpi.h><br>int main(int argc, char** argv)<br>{<br> PetscErrorCode err;<br> err = PetscInitialize(&argc, &argv, NULL, "help");<br> CHKERRQ(err);<br> // create a sparse AIJ matrix distributed across MPI<br> PetscInt global_width = 6;<br> PetscInt global_height = 3;<br> Mat A;<br> err = MatCreate(MPI_COMM_WORLD, &A);<br> CHKERRQ(err);<br> err = MatSetType(A, MATMPIAIJ);<br> CHKERRQ(err);<br> // setup pre-allocation for matrix space<br> {<br> err =<br> MatSetSizes(A, global_height, PETSC_DECIDE, global_height, global_width);<br> CHKERRQ(err);<br> PetscInt d_nnz[] = {0, 0, 0};<br> PetscInt o_nnz[] = {3, 3, 3};<br> err = MatMPIAIJSetPreallocation(A, 0, d_nnz, 0, o_nnz);<br> CHKERRQ(err);<br> }<br> err = MatSetUp(A);<br> CHKERRQ(err);<br> // set values inside the matrix<br> for (PetscInt row = 0; row < global_height; ++row)<br> {<br> for (PetscInt col = global_height; col < global_width; ++col)<br> {<br> err = MatSetValue(A, row, col, 1, INSERT_VALUES);<br> CHKERRQ(err);<br> }<br> }<br> err = MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);<br> CHKERRQ(err);<br> err = MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);<br> CHKERRQ(err);<br> err = MatView(A, PETSC_VIEWER_STDOUT_WORLD);<br> CHKERRQ(err);<br> // free memory<br> err = MatDestroy(&A);<br> CHKERRQ(err);<br> // cleanup any internal PETSc data at end of program<br> err = PetscFinalize();<br> CHKERRQ(err);<br>}</blockquote></div><div><br></div><div>Am I mis-understanding what the d_nnz and o_nnz parameter are supposed to mean?</div><span class="m_-5553325249802122313HOEnZb"><font color="#888888">-- <br><div class="m_-5553325249802122313m_3224953740535681211gmail_signature"><div dir="ltr">Andrew Ho</div></div></font></span></div></div>
</blockquote></div></div></div><span class="HOEnZb"><font color="#888888"><br><br clear="all"><div><br></div>-- <br><div class="m_-5553325249802122313gmail_signature" data-smartmail="gmail_signature">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>
</font></span></div></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr">Andrew Ho</div></div>
</div>