<div dir="ltr"><div dir="ltr">On Fri, May 20, 2022 at 4:45 PM Mike Michell <<a href="mailto:mi.mike1021@gmail.com">mi.mike1021@gmail.com</a>> wrote:<br></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 dir="ltr"><div dir="ltr">Thanks for the reply. <br><br>> "What I want to do is to exchange data (probably just MPI_Reduce)" which confuses me, because halo exchange is a point-to-point exchange and not a reduction.  Can you clarify?<br>PetscSFReduceBegin/End seems to be the function that do reduction for PetscSF object. What I intended to mention was either reduction or exchange, not specifically intended "reduction".<br><br>As a follow-up question: <br>Assuming that the code has its own local solution arrays (not Petsc type), and if the plex's DAG indices belong to the halo region are the only information that I want to know (not the detailed section description, such as degree of freedom on vertex, cells, etc.). I have another PetscSection for printing out my solution. <br>Also if I can convert that DAG indices into my local cell/vertex index, can I just use the PetscSF object created from DMGetPointSF(), instead of "creating PetscSection + DMGetSectionSF()"? In other words, can I use the PetscSF object declared from DMGetPointSF() for the halo communication?<br></div></div></blockquote><div><br></div><div>No, because that point SF will index information by point number. You would need to build a new SF that indexes your dofs. The steps you would</div><div>go through are exactly the same as you would if you just told us what the Section is that indexes your data.</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:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div dir="ltr">Thanks,<br>Mike<br></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr"><br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>The PetscSF that is created automatically is the "point sf" (<a href="https://petsc.org/main/docs/manualpages/DM/DMGetPointSF/" target="_blank">https://petsc.org/main/docs/manualpages/DM/DMGetPointSF/</a>): it says which mesh points (cells, faces, edges and vertices) are duplicates of others.</div><div><br></div><div>In a finite volume application we typically want to assign degrees of freedom just to cells: some applications may only have one degree of freedom, others may have multiple.</div><div><br></div><div>You encode where you want degrees of freedom in a PetscSection and set that as the section for the DM in DMSetLocalSection() (<a href="https://petsc.org/release/docs/manualpages/DM/DMSetLocalSection.html" target="_blank">https://petsc.org/release/docs/manualpages/DM/DMSetLocalSection.html</a>)</div><div><br></div><div>(A c example of these steps that sets degrees of freedom for *vertices* instead of cells is `src/dm/impls/plex/tutorials/ex7.c`)<br></div><div><br></div><div>After that you can call DMGetSectionSF() (<a href="https://petsc.org/main/docs/manualpages/DM/DMGetSectionSF/" target="_blank">https://petsc.org/main/docs/manualpages/DM/DMGetSectionSF/</a>) to the the PetscSF that you want for halo exchange: the one for your solution variables.</div><div><br></div><div>After that, the only calls you typically need in a finite volume code is PetscSFBcastBegin() to start a halo exchange and PetscSFBcastEnd() to complete it.</div><div><br></div><div>You say</div><div><br></div><div>> What I want to do is to exchange data (probably just MPI_Reduce)</div><div><br></div><div>which confuses me, because halo exchange is a point-to-point exchange and not a reduction.  Can you clarify?<br> </div><div><br></div><div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Fri, May 20, 2022 at 8:35 PM Mike Michell <<a href="mailto:mi.mike1021@gmail.com" target="_blank">mi.mike1021@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Dear PETSc developer team,<br><br>Hi, I am using DMPlex for a finite-volume code and trying to understand the usage of PetscSF. What is a typical procedure for doing halo data exchange at parallel boundary using PetscSF object on DMPlex? Is there any example that I can refer to usage of PetscSF with distributed DMPlex? <br><br>Assuming to use the attached mock-up code and mesh, if I give "-dm_distribute_overlap 1 -over_dm_view" to run the code, I can see a PetscSF object is already created, although I have not called "PetscSFCreate" in the code. How can I import & use that PetscSF already created by the code to do the halo data exchange?<br><br>What I want to do is to exchange data (probably just MPI_Reduce) in a parallel boundary region using PetscSF and its functions. I might need to have an overlapping layer or not. <br><br>Thanks,<br>Mike<br></div>
</blockquote></div>
</blockquote></div></div>
</blockquote></div><br clear="all"><div><br></div>-- <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="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>