<div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Fri, Oct 5, 2018 at 2:58 PM Tristan Konolige <<a href="mailto:Tristan.Konolige@colorado.edu">Tristan.Konolige@colorado.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hello PETSc Developers,<br>
<br>
I’d like to get an ISLocalToGlobalMapping from a DMPlex. But, currently, DMGetLocalToGlobalMapping returns an empty ISLocalToGlobalMapping. I can see that ISLocalToGlobalMappings are created in DMPlexDistribute, but they are never assigned to the DM. Is DMGetLocalToGlobalMapping supported with DMPlex or is this a bug?<br></blockquote><div><br></div><div>LocalToGlobal mappings are only meaningful for function spaces, not topology. DMPlex itself only deals with topology. You have to give it a PetscSection if you want to say something about the function space, using DMSetSection().</div><div><br></div><div>DMDA only has one function space (values at vertices), so it just uses this by default, and GetLocalToGlobal() works automatically.</div><div><br></div><div>  Thanks,</div><div><br></div><div>    Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
I’ve included a simple example below.<br>
<br>
```<br>
#include <petsc.h><br>
<br>
int main(int argc, char** argv) {<br>
  PetscErrorCode ierr;<br>
  PetscInt size;<br>
  ISLocalToGlobalMapping l2g;<br>
  DM dm, dm_dist;<br>
  PetscInt num_points[2] = {2, 2};<br>
  PetscInt cone_size[4] = {1, 1, 0, 0};<br>
  PetscInt cones[2] = {2,3};<br>
  PetscInt orientations[2] = {0, 0};<br>
  PetscInt local_points[2] = {0, 1};<br>
  PetscInt global_points[2];<br>
<br>
  PetscInitialize(&argc, &argv, NULL, NULL);<br>
<br>
  DMPlexCreate(PETSC_COMM_WORLD, &dm);<br>
  DMPlexCreateFromDAG(dm, 1, num_points, cone_size, cones, orientations, NULL);<br>
  DMPlexDistribute(dm, 0, NULL, &dm_dist);<br>
<br>
  DMGetLocalToGlobalMapping(dm_dist, &l2g);<br>
  ISLocalToGlobalMappingApply(l2g, 2, local_points, global_points);<br>
<br>
  PetscFinalize();<br>
}<br>
```<br>
<br>
Thank you,<br>
Tristan Konolige<br>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature" data-smartmail="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>