<html><head><meta http-equiv="Content-Type" content="text/html; charset=utf-8"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><br class=""><div><blockquote type="cite" class=""><div class="">On 22 Sep 2019, at 6:03 PM, Smith, Barry F. <<a href="mailto:bsmith@mcs.anl.gov" class="">bsmith@mcs.anl.gov</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div class=""><br class=""><br class=""><blockquote type="cite" class="">On Sep 22, 2019, at 10:14 AM, Pierre Jolivet via petsc-dev <<a href="mailto:petsc-dev@mcs.anl.gov" class="">petsc-dev@mcs.anl.gov</a>> wrote:<br class=""><br class="">FWIW, I’ve fixed MatMatMult and MatTransposeMatMult here <a href="https://gitlab.com/petsc/petsc/commit/93d7d1d6d29b0d66b5629a261178b832a925de80" class="">https://gitlab.com/petsc/petsc/commit/93d7d1d6d29b0d66b5629a261178b832a925de80</a> (with MAT_INITIAL_MATRIX).<br class="">I believe there is something not right in your MR (2032) with MAT_REUSE_MATRIX (without having called MAT_INITIAL_MATRIX first), cf. <a href="https://gitlab.com/petsc/petsc/merge_requests/2069#note_220269898" class="">https://gitlab.com/petsc/petsc/merge_requests/2069#note_220269898</a>.<br class="">Of course, I’d love to be proved wrong!<br class=""></blockquote><br class=""> I don't understand the context.<br class=""><br class=""> MAT_REUSE_MATRIX requires that the C matrix has come from a previous call with MAT_INITIAL_MATRIX, you cannot just put any matrix in the C location.</div></div></blockquote><div><br class=""></div><div>1) It was not the case before the MR, I’ve used that “feature” (which may be specific for MatMatMult_MPIAIJ_MPIDense) for as long as I can remember</div><div>2) If it is not the case anymore, I think it should be mentioned somewhere (and not only in the git log, because I don’t think all users will go through that)</div><div>3) This comment should be removed from the code as well: <a href="https://www.mcs.anl.gov/petsc/petsc-dev/src/mat/impls/aij/mpi/mpimatmatmult.c.html#line398" class="">https://www.mcs.anl.gov/petsc/petsc-dev/src/mat/impls/aij/mpi/mpimatmatmult.c.html#line398</a></div><br class=""><blockquote type="cite" class=""><div class=""><div class="">This is documented in the manual page. We should have better error checking that this is the case so the code doesn't crash at memory access but instead produces a very useful error message if the matrix was not obtained with MAT_INITIAL_MATRIX. <br class=""><br class=""> Is this the issue or do I not understand?<br class=""></div></div></blockquote><div><br class=""></div><div>This is exactly the issue.</div><br class=""><blockquote type="cite" class=""><div class=""><div class=""> Barry<br class=""><br class="">BTW: yes MAT_REUSE_MATRIX has different meanings for different matrix operations in terms of where the matrix came from, this is suppose to be all documented in each methods manual page but some may be missing or incomplete, and error checking is probably not complete for all cases. Perhaps the code should be changed to have multiple different names for each reuse case for clarity to user?<br class=""></div></div></blockquote><div><br class=""></div><div>Definitely, cf. above.</div><div><br class=""></div><div>Thanks,</div><div>Pierre</div><br class=""><blockquote type="cite" class=""><div class=""><div class=""><blockquote type="cite" class=""><br class="">Thanks,<br class="">Pierre<br class=""><br class=""><blockquote type="cite" class="">On 22 Sep 2019, at 5:04 PM, Zhang, Hong <<a href="mailto:hzhang@mcs.anl.gov" class="">hzhang@mcs.anl.gov</a>> wrote:<br class=""><br class="">I'll check it tomorrow.<br class="">Hong<br class=""><br class="">On Sun, Sep 22, 2019 at 1:04 AM Pierre Jolivet via petsc-dev <<a href="mailto:petsc-dev@mcs.anl.gov" class="">petsc-dev@mcs.anl.gov</a>> wrote:<br class="">Jed,<br class="">I’m not sure how easy it is to put more than a few lines of code on GitLab, so I’ll just send the (tiny) source here, as a follow-up of our discussion <a href="https://gitlab.com/petsc/petsc/merge_requests/2069#note_220229648" class="">https://gitlab.com/petsc/petsc/merge_requests/2069#note_220229648</a>.<br class="">Please find attached a .cpp showing the brokenness of C=A*B with A of type MPIAIJ and B of type MPIDense when the LDA of B is not equal to its number of local rows.<br class="">It does [[1,1];[1,1]] * [[0,1,2,3];[0,1,2,3]]<br class="">C should be equal to 2*B, but it’s not, unless lda = m (= 1).<br class="">Mat Object: 2 MPI processes<br class=""> type: mpidense<br class="">0.0000000000000000e+00 1.0000000000000000e+00 2.0000000000000000e+00 3.0000000000000000e+00<br class="">0.0000000000000000e+00 1.0000000000000000e+00 2.0000000000000000e+00 3.0000000000000000e+00<br class=""><br class="">If you change Bm here <a href="https://www.mcs.anl.gov/petsc/petsc-dev/src/mat/impls/aij/mpi/mpimatmatmult.c.html#line549" class="">https://www.mcs.anl.gov/petsc/petsc-dev/src/mat/impls/aij/mpi/mpimatmatmult.c.html#line549</a> to the LDA of B, you’ll get the correct result.<br class="">Mat Object: 2 MPI processes<br class=""> type: mpidense<br class="">0.0000000000000000e+00 2.0000000000000000e+00 4.0000000000000000e+00 6.0000000000000000e+00<br class="">0.0000000000000000e+00 2.0000000000000000e+00 4.0000000000000000e+00 6.0000000000000000e+00<br class=""><br class="">Unfortunately, w.r.t. MR 2069, I still don’t get the same results with a plain view LDA > m (KO) and a view + duplicate LDA = m (OK).<br class="">So there might be something else to fix (or this might not even be a correct fix), but the only reproducer I have right now is the full solver.<br class=""><br class="">Thanks,<br class="">Pierre<br class=""><br class=""></blockquote><br class=""></blockquote><br class=""></div></div></blockquote></div><br class=""></body></html>