<div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div><br></div><div><br><div class="gmail_quote"><div dir="ltr">On Tue, 6 Nov 2018 at 16:01, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@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"><div class="gmail_quote"><div dir="ltr">On Tue, Nov 6, 2018 at 9:45 AM Francesco Magaletti <<a href="mailto:francesco_magaletti@fastwebnet.it" target="_blank">francesco_magaletti@fastwebnet.it</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 style="overflow-wrap: break-word;">Hi Matt, <div>thanks for the reply, but I didn’t catch exactly your point. What do you mean with “first bin the points into cells, and then use the</div><div>interpolant from the discretization” ? To me it sounds like my “naive” approach, but maybe I missed something in your suggestion. </div></div></blockquote><div><br></div><div>I did not read the whole mail before. Your naive approach is exactly right. However, in parallel, you need to keep track of the bounds of each process, send the points to the right process, and then do local location and interpolation. It is some bookkeeping.</div><div><br></div><div>Explaining this gave me another idea. I think that DMSwarm might do this automatically. You create a DMSwarm, stick in</div><div>all the points you need to interpolate at, set the cellDM to the DMDA used for interpolation. Then you call DMSwarmMigrate()</div><div>to put the points on the correct process, and then interpolate. Tell me if this does not work. I am Cc'ing the author Dave May</div><div>in case I have made an error.</div></div></div></blockquote><div dir="auto"><br></div><div dir="auto">It will work using dmswarm but not exactly as Matt describes (although I might have missed some details in your naive approach). </div><div dir="auto"><br></div><div dir="auto">Here is why I think Matts exact definition is not quite what you want. When you set the dmda on the swarm object, the swarm uses it for point location and selects a particular communication pattern for the data exchange. The point location indicates if a point was contained in the sub-domain. If the answer is no, the point is scattered to all ranks which are neighbours of the current sub-domain. The point location does not identify which sub-domain contains each point.</div><div dir="auto">Hence I think this comm pattern may not be sufficient for what you want to do (in general)</div><div dir="auto"><br></div><div dir="auto">The alternative approach which will work is to not attach the dmda to the swarm. In this case, you are responsible to determine where the point should go and assign the target rank value to an internal swarm variable (textual name is "DMSwarm_rank"). In this mode you define the communication pattern, cf. inferring it from the dmda decomposition. </div><div dir="auto"><br></div><div>Here is what I would do:</div><div>* I'll assume you want to interpolate from dmda2 to dmda1 and you have fields Vec F2 and Vec F1 defined on dmda2 and dmda1.</div><div>* I'll assume that you ultimately want the interpolated field values to live in F1. If that is not the case the procedure described below will be slightly different.</div><div>* Below is some pseudo code. I'll write things swarm[i]->"xxx" to indicate registered entires in the swam named "xxx" which live in the i-th slot of the dmswarm's internal storage.</div><div dir="auto"><br></div><div>[1] </div><div>Create a swarm and register the following fields</div><div>"coor" (PetscReal)</div><div>"dmda_index" (PetscInt)</div><div><br></div><div>[2] </div><div>Get the sub-domain bounding box from dmda2 and broadcast these to all ranks, so that every rank knows the bounds of all sub-domains associated with dmda2</div><div><br></div><div>[3] </div><div>Set the local size of the swarm to match the number of locally owned points in dmda1. Do not include the ghost points.</div><div>Copy all the coordinates from dmda1 into the registered swarm field "coor".</div><div>Copy the global index of each point in dmda1 into the registered field "dmda_index"</div><div><br></div><div>[4] </div><div>Get the internally registered dmswarm field named "DMSwarm_rank".</div><div>for all points in the swarm, i</div><div> swarm[i]->"DMSwarm_rank" = comm.rank</div><div> for all bounding boxes of dmda2, b[k]</div><div> if swarm[i]->"coor" is contained in b[k]</div><div> swarm[i]->"DMSwarm_rank" = k</div><div> break</div><div><br></div><div>[5]</div><div dir="auto">Call DMSwarmMigrate(swarm,PETSC_TRUE) and you data will get shipped to where you indicated it should go. </div><div>Now all the points you want to interpolate to should be located within the correct sub-domains of dmda2.</div><div>Note we will remove the sent points from the local dmswarm storage</div><div><br></div><div>[6]</div><div>Get the registered swarm field "coor" and "dmda_index"</div><div>VecZeroEntries(F1)</div><div><div>for all points in the swarm, i</div></div><div> locate which cell in dmda2 contains swarm[i]->"coor"</div><div> /* interpolate F2 defined on dmda2 to swarm[i]->"coor" - call this F_interp */</div><div> F_interp = .....</div><div> VecSetValue(F1,swarm[i]->"dmda_index",F_interp,INSERT_VALUES);</div><div dir="auto"><br></div><div>VecAssemblyBegin(F1)</div><div dir="auto"><div>VecAssemblyBegin(F2)</div><div><br></div><div><br></div><div>Hope this helps.</div><div><br>Thanks,</div><div> Dave</div></div><div dir="auto"><br></div><div dir="auto"><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 class="gmail_quote"><div></div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div></div></div><div dir="ltr"><div class="gmail_quote"><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 style="overflow-wrap: break-word;"><div>My problem comes from the solution of a 1D time evolving pde coupled with a hand-made mesh adapting algorithm. Therefore I’m already using a DMDA to manage all the communication issues among processors. Every n timesteps I move the grid and then I need to evaluate the old solution on the new grid points, this is why I need an efficient way to do it. Is it feasible to let DMDA’s objects speak with DMPlexes, in order to keep the previous code structure unaltered and use DMPlex only for the interpolation routine?</div><div><br></div><div>Thanks</div><div><br></div><div><br></div><div><br><div><div><blockquote type="cite"><div>Il giorno 06/nov/2018, alle ore 15:14, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> ha scritto:</div><br class="gmail-m_-741461629713182391m_5266060687362831173m_-6324208327439674341Apple-interchange-newline"><div><div dir="ltr" style="font-family:Helvetica;font-size:12px;font-style:normal;font-variant:normal;font-weight:normal;letter-spacing:normal;line-height:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px"><div class="gmail_quote"><div dir="ltr">On Tue, Nov 6, 2018 at 5:11 AM Francesco Magaletti via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</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">Dear all,<span class="gmail-m_-741461629713182391m_5266060687362831173m_-6324208327439674341Apple-converted-space"> </span><br><br>I would like to ask you if is there an easy and efficient way in parallel (by using PetsC functions) to interpolate a DMDA vector associated with a nonuniform 1D grid to another DMDA vector with the same length but associated with a different nonuniform grid.<span class="gmail-m_-741461629713182391m_5266060687362831173m_-6324208327439674341Apple-converted-space"> </span><br><br>Let me rephrase it to be as clearer as I can:<br><br>I have two structured nonuniform 1D grids with coordinate vectors x[i] and y[i]. Both the domains have been discretized with the same number of points, but the coordinate vectors x and y are different. I have a discretized field u[i] = u(x[i]) and I would like to use these point values to evaluate the values u(y[i]) in the points of the second grid.<br><br>I read on the manual pages that functions like DMCreateInterpolation or similar work only with different but uniform DMDAs. Did I understand correctly?<span class="gmail-m_-741461629713182391m_5266060687362831173m_-6324208327439674341Apple-converted-space"> </span><br><br>A naive approach, with a serial code, could be to find the points x[i] and x[i+1] that surround the point y[j] for every j and then simply linear interpolating the values u[i] and u[i+1]. I suspect that this is not the most efficient way to do it. Moreover it won’t work in parallel since, in principle, I do not know beforehand how many ghost nodes could be necessary to perform all the interpolations.<span class="gmail-m_-741461629713182391m_5266060687362831173m_-6324208327439674341Apple-converted-space"> </span><br><br>Thank you in advance for your help!<br></blockquote><div><br></div><div>This has not been written, but is not that hard. You would first bin the points into cells, and then use the</div><div>interpolant from the discretization. This is how we do it in the unstructured case. Actually, if you wrote</div><div>your nonuniform 1D meshes as DMPlexes, I think it would work right now (have not tested).</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">Francesco</blockquote></div><br clear="all"><div><br></div>--<span class="gmail-m_-741461629713182391m_5266060687362831173m_-6324208327439674341Apple-converted-space"> </span><br><div dir="ltr" class="gmail-m_-741461629713182391m_5266060687362831173m_-6324208327439674341gmail_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></div></div></div></div></div></div></div></div></div></blockquote></div><br></div></div></div></blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail-m_-741461629713182391m_5266060687362831173gmail_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>
</blockquote></div></div>
</div></div></div></div></div></div></div></div>