<div dir="ltr">I don't understand how the matrices are partitioned then. The documentation online all shows each mpi rank owning entire rows, with no partitioning along columns.<div><br></div><div>I've attached an example where I try to partition a 6x7 matrix into 4 ranks:</div><div><br></div><div>rank 0 should own a 3x4 section</div><div>rank 1 should own a 3x3 section</div><div>rank 2 should own a 3x4 section</div><div>rank 3 should own a 3x3 section</div><div><br></div><div>However, when I run the example and print out the ownership range and ownership column range, I get:</div><div><br></div><div><div>rank 0 owns rows/cols: [0,3), [0, 4)</div><div>rank 1 owns rows/cols: [3,6), [4, 7)</div><div>rank 2 owns rows/cols: [6,9), [7, 11)</div><div>rank 3 owns rows/cols: [9,12), [11, 14)</div></div><div><br></div><div>Which makes no sense since these ranges extend beyond the actual size of the matrix.</div></div><div class="gmail_extra"><br><div class="gmail_quote">On Tue, Feb 14, 2017 at 6:05 PM, Barry Smith <span dir="ltr"><<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><br>
1) on a single process all the columns are in the diagonal block and none in the off diagonal block (note that the column partitioning corresponding to a partitioning of the vector in the product A *x<br>
<br>
2) on two processes you guessed wrong how many columns PETSc would put on the first process. You guessed it would put the first three on the first process and the last three on the second process.<br>
<br>
a) it cannot do that, every column has to been owned by some process (in the vector x above) so it cannot be 3 and 3, it has to be 4 and 3 or 3 and 4.<br>
<br>
b) PETSc puts the "extra" columns on the earlier processes not the later processes.<br>
<br>
For rectangular matrices you really cannot get away with using "PETSC_DECIDE" for local columns when using preallocation<br>, since it may decide something different than what you assume.<br>
<br>
I've attached the parallel code that behaves as expected.<br>
<br>
Barry<br>
<br>
> On Feb 14, 2017, at 1:06 PM, Andrew Ho <<a href="mailto:andrewh0@uw.edu">andrewh0@uw.edu</a>> wrote:<br>
><br>
> The problem isn't only with 1 process, but it seems to be with all non-square matrices?<br>
><br>
> For example, here I have a two process program which tries to set a 6x7 matrix to (each process owns 3 rows):<br>
><br>
> 0 0 0 1 1 1 1<br>
> 0 0 0 1 1 1 1<br>
> 0 0 0 1 1 1 1<br>
> 0 0 0 1 1 1 1<br>
> 0 0 0 1 1 1 1<br>
> 0 0 0 1 1 1 1<br>
><br>
> #include <petsc.h><br>
> #include <mpi.h><br>
><br>
> int main(int argc, char** argv)<br>
> {<br>
> PetscErrorCode err;<br>
> err = PetscInitialize(&argc, &argv, NULL, "help");<br>
> CHKERRQ(err);<br>
><br>
> int rank, size;<br>
> MPI_Comm_rank(MPI_COMM_WORLD, &rank);<br>
> MPI_Comm_size(MPI_COMM_WORLD, &size);<br>
> if(size != 2)<br>
> {<br>
> printf("must run with 2 processes");<br>
> MPI_Abort(MPI_COMM_WORLD, -1);<br>
> }<br>
><br>
> // create a sparse AIJ matrix distributed across MPI<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, 3, PETSC_DECIDE, 6, 7);<br>
> CHKERRQ(err);<br>
> if(rank == 0)<br>
> {<br>
> PetscInt d_nnz[] = {0, 0, 0};<br>
> PetscInt o_nnz[] = {4, 4, 4};<br>
> err = MatMPIAIJSetPreallocation(A, 0, d_nnz, 0, o_nnz);<br>
> CHKERRQ(err);<br>
> }<br>
> else<br>
> {<br>
> PetscInt d_nnz[] = {3, 3, 3};<br>
> PetscInt o_nnz[] = {1, 1, 1};<br>
> err = MatMPIAIJSetPreallocation(A, 0, d_nnz, 0, o_nnz);<br>
> CHKERRQ(err);<br>
> }<br>
> }<br>
> err = MatSetUp(A);<br>
> CHKERRQ(err);<br>
><br>
> // set values inside the matrix<br>
> for (PetscInt row = 0; row < 3; ++row)<br>
> {<br>
> for (PetscInt col = 3; col < 7; ++col)<br>
> {<br>
> err = MatSetValue(A, 3 * rank + row, col, 1, INSERT_VALUES);<br>
> CHKERRQ(err);<br>
> }<br>
> }<br>
><br>
> err = MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);<br>
> CHKERRQ(err);<br>
> err = MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);<br>
> CHKERRQ(err);<br>
><br>
> err = MatView(A, PETSC_VIEWER_STDOUT_WORLD);<br>
> CHKERRQ(err);<br>
><br>
> // free memory<br>
> err = MatDestroy(&A);<br>
> CHKERRQ(err);<br>
><br>
> // cleanup any internal PETSc data at end of program<br>
> err = PetscFinalize();<br>
> CHKERRQ(err);<br>
> }<br>
><br>
><br>
> On Tue, Feb 14, 2017 at 5:26 AM, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> wrote:<br>
> On Tue, Feb 14, 2017 at 6:41 AM, Andrew Ho <<a href="mailto:andrewh0@uw.edu">andrewh0@uw.edu</a>> wrote:<br>
> I have a 3x6 matrix, and I want to set it to (just as an example):<br>
><br>
> 0 0 0 1 1 1<br>
> 0 0 0 1 1 1<br>
> 0 0 0 1 1 1<br>
><br>
> From what I can tell, according to the documentation the MPIAIJ sparsity of this matrix is:<br>
><br>
> I believe this is an inconsistency in PETSc for rectangular matrices. We divide matrices into diagonal and<br>
> off-diagonal parts mainly to separate communication from computation. Thus on 1 proc, we never divide<br>
> them, even if the matrix is rectangular. Therefore we misinterpret your preallocation.<br>
><br>
> Barry, do you think we want to try and change this?<br>
><br>
> Matt<br>
><br>
> d_nnz = [0, 0, 0]<br>
> o_nnz = [3, 3, 3]<br>
><br>
> However, when I do this, I get the following error:<br>
><br>
> [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_<wbr>ERR, PETSC_FALSE) to turn off this check<br>
> [0]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/<wbr>documentation/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/<wbr>mpiaij.c<br>
> [0]PETSC ERROR: #2 MatSetValues() line 1190 in petsc/src/mat/interface/<wbr>matrix.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>---<br>
><br>
> Here's a working test code:<br>
><br>
> #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>
> }<br>
><br>
> Am I mis-understanding what the d_nnz and o_nnz parameter are supposed to mean?<br>
> --<br>
> Andrew Ho<br>
><br>
><br>
><br>
> --<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<br>
><br>
><br>
><br>
> --<br>
> Andrew Ho<br>
<br>
<br></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>