<html><head><meta http-equiv="Content-Type" content="text/html; charset=us-ascii"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On 4 Dec 2020, at 1:13 PM, Matthew Knepley <<a href="mailto:knepley@gmail.com" class="">knepley@gmail.com</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="ltr" class=""><div dir="ltr" class="">On Fri, Dec 4, 2020 at 6:51 AM Roland Richter <<a href="mailto:roland.richter@ntnu.no" class="">roland.richter@ntnu.no</a>> wrote:<br class=""></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div class=""><p class="">Hei,</p><p class="">is it possible to fill a dense distributed matrix from existing
data in rank 0, distribute it to all involved nodes and retrieve
it afterwards, such that it can be stored in a single matrix in
rank 0 again? The background behind this question is the following
thought process:</p>
<ul class="">
<li class="">I generate the matrix locally on node 0</li></ul></div></blockquote><div class="">I do not understand why you would ever want to do this. <br class=""></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><ul class="">
<li class="">I distribute it to all involved nodes</li></ul></div></blockquote><div class="">Just create your distributed matrix, but set the values from rank 0. <br class=""></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><ul class="">
<li class="">I use the matrix for data processing/process the matrix data
on all nodes</li></ul></div></blockquote><div class="">I assume you change the matrix here. </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><ul class="">
<li class="">For storing the data I need it back on node 0 in a single
matrix</li></ul></div></blockquote><div class="">You can use <a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateRedundantMatrix.html" class="">https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateRedundantMatrix.html</a></div><div class="">to easily create serial matrices from it. You can throw away any you do not want.</div><div class=""><br class=""></div><div class=""> Thanks,</div><div class=""><br class=""></div><div class=""> Matt </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><ul class="">
</ul><p class="">For testing I wrote the following code (and using armadillo for
generating the initial matrices within arma-namespace):</p><p class=""><i class=""> Mat C, F;</i><i class=""><br class="">
</i><i class=""> Vec x, y, z;</i><i class=""><br class="">
</i><i class=""> PetscViewer viewer;</i><i class=""><br class="">
</i><i class=""> PetscMPIInt rank, size;</i><i class=""><br class="">
</i><i class=""> PetscInitialize(&argc, &args, (char*) 0, help);</i><i class=""><br class="">
</i><i class=""><br class="">
</i><i class=""> MPI_Comm_size(PETSC_COMM_WORLD, &size);</i><i class=""><br class="">
</i><i class=""> MPI_Comm_rank(PETSC_COMM_WORLD, &rank);</i><i class=""><br class="">
</i><i class=""><br class="">
</i><i class=""> PetscPrintf(PETSC_COMM_WORLD,"Number of processors =
%d, rank = %d\n", size, rank);</i><i class=""><br class="">
</i><i class=""></i><i class=""><br class="">
</i><i class=""> PetscViewerCreate(PETSC_COMM_WORLD, &viewer);</i><i class=""><br class="">
</i><i class=""> PetscViewerSetType(viewer, PETSCVIEWERASCII);</i><i class=""><br class="">
</i><i class=""> arma::cx_mat local_mat, local_zero_mat;</i><i class=""><br class="">
</i><i class=""> const size_t matrix_size = 5;</i><i class=""><br class="">
</i><i class=""> if(rank == 0) {</i><i class=""><br class="">
</i><i class=""> local_mat =
arma::randu<arma::cx_mat>(matrix_size, matrix_size);</i><i class=""><br class="">
</i><i class=""> local_zero_mat =
arma::zeros<arma::cx_mat>(matrix_size, matrix_size);</i><i class=""><br class="">
</i><i class=""> MatCreateDense(PETSC_COMM_WORLD, matrix_size,
matrix_size, PETSC_DECIDE, PETSC_DECIDE, local_mat.memptr(),
&C);</i><i class=""><br class="">
</i><i class=""> MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY);</i><i class=""><br class="">
</i><i class=""> MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY);</i><i class=""><br class="">
</i><i class=""> }</i><i class=""><br class="">
</i><i class=""> MatCreateDense(PETSC_COMM_WORLD, matrix_size,
matrix_size, PETSC_DECIDE, PETSC_DECIDE, NULL, &F);</i><i class=""><br class="">
</i><i class=""> if(rank == 0)</i><i class=""><br class="">
</i><i class=""> MatCopy(C, F, SAME_NONZERO_PATTERN);</i><i class=""><br class="">
</i><i class=""> MatAssemblyBegin(F, MAT_FINAL_ASSEMBLY);</i><i class=""><br class="">
</i><i class=""> MatAssemblyEnd(F, MAT_FINAL_ASSEMBLY);</i><i class=""><br class="">
</i><i class=""> if(rank == 0) {</i><i class=""><br class="">
</i><i class=""> MatView(C, viewer);</i><i class=""><br class="">
</i><i class=""> std::cout << local_mat << '\n';</i><i class=""><br class="">
</i><i class=""> std::cout << local_zero_mat << '\n';</i><i class=""><br class="">
</i><i class=""> }</i><i class=""><br class="">
</i><i class=""> //How can I retrieve the distributed data from F to C
again when running on multiple nodes?</i><i class=""><br class="">
</i><i class=""><br class="">
</i><i class=""> if(rank == 0) {</i><i class=""><br class="">
</i><i class=""> MatDestroy(&C);</i><i class=""><br class="">
</i><i class=""> }</i><i class=""><br class="">
</i><i class=""> MatDestroy(&F);</i><i class=""><br class="">
</i><i class=""> PetscViewerDestroy(&viewer);</i></p><p class="">Running the code with mpirun -n 1 works fine, but when running it
with mpirun -n 2 it just stops after <i class="">PetscPrintf() </i>and
hangs.</p><p class="">Moreover, how can I retrieve the data from matrix F afterwards?
For vectors I have the function/class VecScatter, but can I apply
this for matrices (mpidense) too?</p></div></blockquote></div></div></div></blockquote><div>With respect to this very last question, this is currently an open issue <a href="https://gitlab.com/petsc/petsc/-/issues/693" class="">https://gitlab.com/petsc/petsc/-/issues/693</a>.</div><div><br class=""></div><div>Thanks,</div><div>Pierre</div><blockquote type="cite" class=""><div class=""><div dir="ltr" class=""><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><p class="">Thanks!</p><p class="">Cheers,</p><p class="">Roland<br class="">
</p>
</div>
</blockquote></div><br clear="all" class=""><div class=""><br class=""></div>-- <br class=""><div dir="ltr" class="gmail_signature"><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div class="">What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br class="">-- Norbert Wiener</div><div class=""><br class=""></div><div class=""><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank" class="">https://www.cse.buffalo.edu/~knepley/</a><br class=""></div></div></div></div></div></div></div></div>
</div></blockquote></div><br class=""></body></html>