On Thu, Jun 16, 2022 at 12:06 PM Zongze Yang <[email protected]> wrote:
> > > 在 2022年6月16日，23:22，Matthew Knepley <[email protected]> 写道： > > > On Thu, Jun 16, 2022 at 11:11 AM Zongze Yang <[email protected]> wrote: > >> Hi, if I load a `gmsh` file with second-order elements, the coordinates >> will be stored in a DG-P2 space. After obtaining the coordinates of a cell, >> how can I map the coordinates to vertex and edge? >> > > By default, they are stored as P2, not DG. > > > I checked the coordinates vector, and found the dogs only defined on cell > other than vertex and edge, so I said they are stored as DG. > Then the function DMPlexVecGetClosure > <https://petsc.org/main/docs/manualpages/DMPLEX/DMPlexVecGetClosure/> seems > return > the coordinates in lex order. > > Some code in reading gmsh file reads that > > > 1756: if (isSimplex) continuity = PETSC_FALSE > <https://petsc.org/main/docs/manualpages/Sys/PETSC_FALSE/>; /* XXX FIXME > Requires DMPlexSetClosurePermutationLexicographic() */ > > > 1758: GmshCreateFE(comm, NULL, isSimplex, continuity, nodeType, dim, > coordDim, order, &fe) > > > The continuity is set to false for simplex. > Oh, yes. That needs to be fixed. For now, you can just project it to P2 if you want using https://petsc.org/main/docs/manualpages/DM/DMProjectCoordinates/ Thanks, Matt > Thanks, > Zongze > > You can ask for the coordinates of a vertex or an edge directly using > > https://petsc.org/main/docs/manualpages/DMPLEX/DMPlexPointLocalRead/ > > by giving the vertex or edge point. You can get all the coordinates on a > cell, in the closure order, using > > https://petsc.org/main/docs/manualpages/DMPLEX/DMPlexVecGetClosure/ > > Thanks, > > Matt > > >> Below is some code load the gmsh file, I want to know the relation >> between `cl` and `cell_coords`. >> >> ``` >> import firedrake as fd >> import numpy as np >> >> # Load gmsh file (2rd) >> plex = fd.mesh._from_gmsh('test-fd-load-p2-rect.msh') >> >> cs, ce = plex.getHeightStratum(0) >> >> cdm = plex.getCoordinateDM() >> csec = dm.getCoordinateSection() >> coords_gvec = dm.getCoordinates() >> >> for i in range(cs, ce): >> cell_coords = cdm.getVecClosure(csec, coords_gvec, i) >> print(f'coordinates for cell {i} :\n{cell_coords.reshape([-1, 3])}') >> cl = dm.getTransitiveClosure(i) >> print('closure:', cl) >> break >> ``` >> >> Best wishes, >> Zongze >> > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ > <http://www.cse.buffalo.edu/~knepley/> > > -- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
