<html><head><meta http-equiv="content-type" content="text/html; charset=utf-8"></head><body style="overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;"><br id="lineBreakAtBeginningOfMessage"><div>
<meta charset="UTF-8"><div dir="auto" style="caret-color: rgb(0, 0, 0); color: rgb(0, 0, 0); letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;"><div>Re:</div><div><br></div><div>Thank you Pierre, I really appreciate. I’am testing it right now to access the improvements. </div><div><br></div><div>BR,</div><div>Medane</div><div><br></div><div><br></div></div></div><div><blockquote type="cite"><div>On 27 Jan 2025, at 20:19, Pierre Jolivet <pierre@joliv.et> wrote:</div><br class="Apple-interchange-newline"><div><meta http-equiv="content-type" content="text/html; charset=utf-8"><div style="overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;">Please always keep the list in copy.<div>The way you create A is not correct, I’ve attached a fixed code.</div><div>If you want to keep your own distribution for A (and not the one associated to R_part), you’ll need to first call <a href="https://urldefense.us/v3/__https://petsc.org/main/manualpages/Mat/MatCreateSubMatrix/__;!!G_uCfscf7eWS!Ye7A4fqD5xLobOPjgvYkh9cj1-JExIqX_EJIHFm-NHw5rEk2PU5kvs3GfKlJd2TZPorWhvb0Jh7eTcKii9t7Z7tgYSIoeSHTchrF1snH$">https://petsc.org/main/manualpages/Mat/MatCreateSubMatrix/</a> to redistribute A and then do a MatCopy() of the resulting Mat into R_part</div><div><br></div><div>Thanks,</div><div>Pierre</div><div><br></div><div>$ /Volumes/Data/repositories/petsc/arch-darwin-c-debug-real/bin/mpirun -n 4 ./ex1234</div><div><div>Mat Object: 4 MPI processes</div><div> type: mpidense</div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div>Mat Object: 4 MPI processes</div><div> type: mpidense</div><div> 2.6219599187040323e+00 1.9661197867318445e+00 1.5218640363910978e+00 </div><div> 3.5202261875977947e+00 3.6311893358251384e+00 2.2279492868785069e+00 </div><div> 2.7505403755038014e+00 3.1546072728892756e+00 1.8416294994524489e+00 </div><div> 2.4676055638467314e+00 2.3185625557889602e+00 2.0401666986599833e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div><div> 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 </div></div><div><br></div><div><br></div><div></div></div><span id="cid:762AAAA5-F332-45DD-8798-85AB922A9453"><ex1234.c></span><meta http-equiv="content-type" content="text/html; charset=utf-8"><div style="overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;"><div><br id="lineBreakAtBeginningOfMessage"><div><br><blockquote type="cite"><div>On 27 Jan 2025, at 6:53 PM, medane.tchakorom@univ-fcomte.fr wrote:</div><br class="Apple-interchange-newline"><div><meta http-equiv="content-type" content="text/html; charset=utf-8"><div style="overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;">Re:<br id="lineBreakAtBeginningOfMessage"><div>
<meta charset="UTF-8"><div dir="auto" style="text-align: start; text-indent: 0px; overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;"><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">This is a small reproductible example using MatDenseGetSubMatrix</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">Command: petscmpiexec -n 4 ./example</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">==========================================================</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt nlines = 8; // lines</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt ncolumns = 3; // columns</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt random_size = 12;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt rank;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt size;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Initialize PETSc</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInitialize(&argc, &args, NULL, NULL);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MPI_Comm_rank(MPI_COMM_WORLD, &rank);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MPI_Comm_size(MPI_COMM_WORLD, &size);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // R_full with all values to zero</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> Mat R_full;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, nlines, ncolumns, NULL, &R_full);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatZeroEntries(R_full);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatView(R_full, PETSC_VIEWER_STDOUT_WORLD);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Creating and setting A and S to rand values</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> Mat A, S;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, nlines / 2, random_size, NULL, &A);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, random_size, ncolumns, NULL, &S);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatSetRandom(A, NULL);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatSetRandom(S, NULL);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Computing R_part</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> Mat R_part;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseGetSubMatrix(R_full, PETSC_DECIDE, nlines / 2, PETSC_DECIDE, PETSC_DECIDE, &R_part);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatMatMult(A, S, MAT_REUSE_MATRIX, PETSC_DECIDE, &R_part);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Visualizing R_full</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseRestoreSubMatrix(R_full, &R_part);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatView(R_full, PETSC_VIEWER_STDOUT_WORLD);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Destroying matrices</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDestroy(&R_part);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDestroy(&R_full);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscFinalize();</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> return 0;</span></font></div></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">==========================================================</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">Part of the error output contains….:</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">"Cannot change/reset row sizes to 1 local 4 global after previously setting them to 2 local 4 global ….”</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">==========================================================</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div><div><font><span style="caret-color: rgb(0, 0, 0);">PetscInt nlines = 8; // lines</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt ncolumns = 3; // columns</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt random_size = 12;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt rank;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInt size;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Initialize PETSc</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscInitialize(&argc, &args, NULL, NULL);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MPI_Comm_rank(MPI_COMM_WORLD, &rank);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MPI_Comm_size(MPI_COMM_WORLD, &size);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // R_full with all values to zero</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> Mat R_full;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, nlines, ncolumns, NULL, &R_full);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatZeroEntries(R_full);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatView(R_full, PETSC_VIEWER_STDOUT_WORLD);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Creating and setting A and S to rand values</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> Mat A, S;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, nlines / 2, random_size, NULL, &A);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, random_size, ncolumns, NULL, &S);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatSetRandom(A, NULL);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatSetRandom(S, NULL);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Computing R_part</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> Mat R_part;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatMatMult(A, S, MAT_INITIAL_MATRIX, PETSC_DECIDE, &R_part);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatView(R_part, PETSC_VIEWER_STDOUT_WORLD);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> Mat R_sub;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseGetSubMatrix(R_full, PETSC_DECIDE, nlines / 2, PETSC_DECIDE, PETSC_DECIDE, &R_sub);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscScalar *storage = NULL;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseGetArray(R_part, &storage);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscScalar *storage_sub = NULL;</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseGetArray(R_sub, &storage_sub);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscArraycpy(storage_sub, storage, (nlines / 2) * ncolumns);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseRestoreArray(R_part, &storage);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseRestoreArray(R_sub, &storage_sub);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDenseRestoreSubMatrix(R_full, &R_sub);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatView(R_full, PETSC_VIEWER_STDOUT_WORLD);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> // Destroying matrices</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDestroy(&R_part);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> MatDestroy(&R_full);</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"><br></span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> PetscFinalize();</span></font></div><div><font><span style="caret-color: rgb(0, 0, 0);"> return 0;</span></font></div></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">==========================================================</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">Now Using MatDenseGetArray</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">Please let me know if I need to clarify something.</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">Thanks</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;">Medane</div><div style="caret-color: rgb(0, 0, 0); letter-spacing: normal; text-transform: none; white-space: normal; word-spacing: 0px; text-decoration: none; -webkit-text-stroke-width: 0px;"><br></div></div>
</div>
<div><br><blockquote type="cite"><div>On 27 Jan 2025, at 16:26, Pierre Jolivet <pierre@joliv.et> wrote:</div><br class="Apple-interchange-newline"><div><meta charset="UTF-8"><br class="Apple-interchange-newline"><br style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;"><blockquote type="cite" style="font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;"><div>On 27 Jan 2025, at 3:52 PM, Matthew Knepley <knepley@gmail.com> wrote:</div><br class="Apple-interchange-newline"><div><div dir="ltr"><div dir="ltr">On Mon, Jan 27, 2025 at 9:23 AM<span class="Apple-converted-space"> </span><a href="mailto:medane.tchakorom@univ-fcomte.fr">medane.tchakorom@univ-fcomte.fr</a><span class="Apple-converted-space"> </span><<a href="mailto:medane.tchakorom@univ-fcomte.fr">medane.tchakorom@univ-fcomte.fr</a>> wrote:</div><div class="gmail_quote gmail_quote_container"><blockquote class="gmail_quote" style="margin: 0px 0px 0px 0.8ex; border-left-width: 1px; border-left-style: solid; border-left-color: rgb(204, 204, 204); padding-left: 1ex;"><div>Re:<div><br></div><div><b>MatDenseGetSubMatrix</b> in fact could be the best alternative (cleaner code), but as mentioned earlier, I would like to use the smaller matrix R_part to get the result of a MatMatMult operation.</div><div><br></div><div> MatMatMult(A, S, MAT_INITIAL_MATRIX, PETSC_DECIDE, &R_part);</div><div><br></div><div>When trying to use MAT_REUSE_MATRIX, this gave an error (expected error I think). On the other side, I mentionned on <b>MatDenseGetSubMatrix<span class="Apple-converted-space"> </span></b>page, "The output matrix is not redistributed by PETSc”. So will R_part be a valid output matrix for MatMatMult?</div></div></blockquote><div><br></div><div>I believe it is, but I am not the expert.</div><div><br></div><div>Pierre, can the SubMatrix be used as the output of a MatMatMult()?</div></div></div></div></blockquote><div style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;"><br></div><div style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;">It should be, but there may be some limitations due to leading dimensions and what not.</div><div style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;">By looking at just the single line of code we got from you, I can see at least one issue: it should be MAT_REUSE_MATRIX, not MAT_INITIAL_MATRIX (assuming you got R_part from MatDenseGetSubMatrix).</div><div style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;">Feel free to share a (minimal) reproducer.</div><div style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;"><br></div><div style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;">Thanks,</div><div style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;">Pierre</div><br style="caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;"><blockquote type="cite" style="font-family: Helvetica; font-size: 12px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;"><div><div dir="ltr"><div class="gmail_quote gmail_quote_container"><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin: 0px 0px 0px 0.8ex; border-left-width: 1px; border-left-style: solid; border-left-color: rgb(204, 204, 204); padding-left: 1ex;"><div><div>Thanks</div><div><br></div><div><br></div><div><div><blockquote type="cite"><div>On 27 Jan 2025, at 14:52, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> wrote:</div><br><div><div dir="ltr"><div dir="ltr">On Mon, Jan 27, 2025 at 8:42 AM Pierre Jolivet <<a href="mailto:pierre@joliv.et" target="_blank">pierre@joliv.et</a>> wrote:</div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin: 0px 0px 0px 0.8ex; border-left-width: 1px; border-left-style: solid; border-left-color: rgb(204, 204, 204); padding-left: 1ex;">> On 27 Jan 2025, at 2:23 PM,<span class="Apple-converted-space"> </span><a href="mailto:medane.tchakorom@univ-fcomte.fr" target="_blank">medane.tchakorom@univ-fcomte.fr</a><span class="Apple-converted-space"> </span>wrote:<br>><span class="Apple-converted-space"> </span><br>> Dear PETSc users,<br>><span class="Apple-converted-space"> </span><br>> I hope this message finds you well. I don’t know If my question is relevant, but I’am currently working with DENSE type matrix, and would like to copy one matrix R_part [ n/2 x m] (resulted from a MatMatMult operation) into another dense matrix R_full [n x m]. <span class="Apple-converted-space"> </span><br>> Both matrices being on the same communicator, I would like to efficiently copy R_part in the first half of R_full.<span class="Apple-converted-space"> </span><br>> I have being using MatSetValues, but for large matrices, the subsequent assembling operation is costly.<br><br>Could you please share the output of -log_view as well as a single file that will be generated with -info dump (e.g., dump.0, the file associated to process #0)?<br>This shouldn’t be that costly, so there may be an option missing, like MAT_NO_OFF_PROC_ENTRIES.<br>Anyway, if you want to optimize this, the fastest way would be to call MatDenseGetArray[Read,Write]() and then do the necessary PetscArraycpy().<br></blockquote><div><br></div><div>The other alternative (which I think makes cleaner code) is to use</div><div><br></div><div> <a href="https://urldefense.us/v3/__https://petsc.org/main/manualpages/Mat/MatDenseGetSubMatrix/__;!!G_uCfscf7eWS!Ye7A4fqD5xLobOPjgvYkh9cj1-JExIqX_EJIHFm-NHw5rEk2PU5kvs3GfKlJd2TZPorWhvb0Jh7eTcKii9t7Z7tgYSIoeSHTchex8FLh$" target="_blank">https://petsc.org/main/manualpages/Mat/MatDenseGetSubMatrix/</a></div><div><br></div><div>to create your R_part matrix. Then you are directly acting on the memory you want when assemble the smaller matrix.</div><div><br></div><div> THanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin: 0px 0px 0px 0.8ex; border-left-width: 1px; border-left-style: solid; border-left-color: rgb(204, 204, 204); padding-left: 1ex;">Thanks,<br>Pierre<br><br>> Please could you suggest me some strategies or functions to do this efficiently.<span class="Apple-converted-space"> </span><br>><span class="Apple-converted-space"> </span><br>> Thank you for your time and assistance.<br>><span class="Apple-converted-space"> </span><br>> Best regards,<br>> Tchakorom Medane<br>><span class="Apple-converted-space"> </span><br><br></blockquote></div><div><br clear="all"></div><div><br></div><span class="gmail_signature_prefix">--<span class="Apple-converted-space"> </span></span><br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><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="https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!Ye7A4fqD5xLobOPjgvYkh9cj1-JExIqX_EJIHFm-NHw5rEk2PU5kvs3GfKlJd2TZPorWhvb0Jh7eTcKii9t7Z7tgYSIoeSHTcjygY4TU$" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div></div></blockquote></div><br></div></div></blockquote></div><div><br clear="all"></div><div><br></div><span class="gmail_signature_prefix">--<span class="Apple-converted-space"> </span></span><br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><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="https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!Ye7A4fqD5xLobOPjgvYkh9cj1-JExIqX_EJIHFm-NHw5rEk2PU5kvs3GfKlJd2TZPorWhvb0Jh7eTcKii9t7Z7tgYSIoeSHTcjygY4TU$" target="_blank">https://www.cse.buffalo.edu/~knepley/</a></div></div></div></div></div></div></div></div></div></blockquote></div></blockquote></div><br></div></div></blockquote></div><br></div></div></div></blockquote></div><br></body></html>