<div dir="ltr"><div dir="ltr">On Wed, May 17, 2023 at 11:20 AM Berend van Wachem <<a href="mailto:berend.vanwachem@ovgu.de">berend.vanwachem@ovgu.de</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">Dear Matt,<br>
<br>
Is there a way to 'redo' the DMLocalizeCoordinates() ? Or to undo it?<br>
Alternatively, can we make the calling of DMLocalizeCoordinates() in the  DMPlexCreate...() routines optional?<br>
<br>
Otherwise, we would have to copy all arrays of coordinates from DMGetCoordinatesLocal() and DMGetCellCoordinatesLocal() before <br>
scaling them.<br></blockquote><div><br></div><div>I am likely not being clear. I think all you have to do is the following:</div><div><br></div><div>  DMGetCoordinatesLocal(dm, &xl);</div><div>  VecScale(xl, scale);</div><div>  DMSetCoordinatesLocal(dm, xl);</div>  DMGetCellCoordinatesLocal(dm, &xl);<br>  VecScale(xl, scale);<br>  DMSetCellCoordinatesLocal(dm, xl);<br><div><br></div><div>Does this not work?</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</div><div><br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
Best regards, Berend.<br>
<br>
On 5/17/23 16:35, Matthew Knepley wrote:<br>
> On Wed, May 17, 2023 at 10:21 AM Berend van Wachem <<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a> <mailto:<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a>>> wrote:<br>
> <br>
>     Dear Matt,<br>
> <br>
>     Thanks for getting back to me so quickly.<br>
> <br>
>     If I scale each of the coordinates of the mesh (say, I want to cube each<br>
>     co-ordinate), and I do this for both:<br>
> <br>
>     DMGetCoordinatesLocal();<br>
>     DMGetCellCoordinatesLocal();<br>
> <br>
>     How do I know I am not cubing one coordinate multiple times?<br>
> <br>
> <br>
> Good question. Right now, the only connection between the two sets of coordinates is DMLocalizeCoordinates(). Since sometimes <br>
> people want to do non-trivial things to<br>
> coordinates, I prefer not to push in an API for "just" scaling, but I could be convinced<br>
> the other way.<br>
> <br>
>    Thanks,<br>
> <br>
>       Matt<br>
> <br>
>     Thanks, Berend.<br>
> <br>
>     On 5/17/23 16:10, Matthew Knepley wrote:<br>
>      > On Wed, May 17, 2023 at 10:02 AM Berend van Wachem<br>
>      > <<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a> <mailto:<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a>> <mailto:<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a><br>
>     <mailto:<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a>>>> wrote:<br>
>      ><br>
>      >     Dear PETSc Team,<br>
>      ><br>
>      >     We are using DMPlex, and we create a mesh using<br>
>      ><br>
>      >     DMPlexCreateBoxMesh (.... );<br>
>      ><br>
>      >     and get a uniform mesh. The mesh is periodic.<br>
>      ><br>
>      >     We typically want to "scale" the coordinates (vertices) of the mesh,<br>
>      >     and<br>
>      >     to achieve this, we call<br>
>      ><br>
>      >     DMGetCoordinatesLocal(dm, &coordinates);<br>
>      ><br>
>      >     and scale the entries in the Vector coordinates appropriately.<br>
>      ><br>
>      >     and then<br>
>      ><br>
>      >     DMSetCoordinatesLocal(dm, coordinates);<br>
>      ><br>
>      ><br>
>      >     After this, we localise the coordinates by calling<br>
>      ><br>
>      >     DMLocalizeCoordinates(dm);<br>
>      ><br>
>      >     This worked fine up to PETSc 3.18, but with versions after this, the<br>
>      >     coordinates we get from the call<br>
>      ><br>
>      >     DMPlexGetCellCoordinates(dm, CellID, &isDG, &CoordSize,<br>
>      >     &ArrayCoordinates, &Coordinates);<br>
>      ><br>
>      >     are no longer correct if the mesh is periodic. A number of the<br>
>      >     coordinates returned from calling DMPlexGetCellCoordinates are wrong.<br>
>      ><br>
>      >     I think, this is because DMLocalizeCoordinates is now automatically<br>
>      >     called within the routine DMPlexCreateBoxMesh.<br>
>      ><br>
>      >     So, my question is: How should we scale the coordinates from a periodic<br>
>      >     DMPlex mesh so that they are reflected correctly when calling both<br>
>      >     DMGetCoordinatesLocal and DMPlexGetCellCoordinates, with PETSc versions<br>
>      >       >= 3.18?<br>
>      ><br>
>      ><br>
>      > I think we might have to add an API function. For now, when you scale<br>
>      > the coordinates,<br>
>      > can you scale both copies?<br>
>      ><br>
>      >    DMGetCoordinatesLocal()<br>
>      >    DMGetCellCoordinatesLocal();<br>
>      ><br>
>      > and then set them back.<br>
>      ><br>
>      >    Thanks,<br>
>      ><br>
>      >       Matt<br>
>      ><br>
>      >     Many thanks, Berend.<br>
>      ><br>
>      > --<br>
>      > What most experimenters take for granted before they begin their<br>
>      > experiments is infinitely more interesting than any results to which<br>
>      > their experiments lead.<br>
>      > -- Norbert Wiener<br>
>      ><br>
>      > <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a> <<a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a>> <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a><br>
>     <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>>><br>
> <br>
> <br>
> <br>
> -- <br>
> What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to <br>
> which their experiments lead.<br>
> -- Norbert Wiener<br>
> <br>
> <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a> <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
</blockquote></div><br clear="all"><div><br></div><span class="gmail_signature_prefix">-- </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="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>