@@ -3930,29 +3930,31 @@ def set_coordinate_section_and_such(mesh, coordinates):
39303930 plex = mesh.topology_dm
39313931 coordinate_element = coordinates.ufl_element()
39323932 gdim = plex.getCoordinateDim()
3933+
3934+ # the existing section has the correct numbering but is scalar, expand to gdim
3935+ scalar_dm = coordinates.function_space().dm
3936+ scalar_section = scalar_dm.getLocalSection()
3937+ p_start, p_end = scalar_section.getChart()
3938+
3939+ vector_section = PETSc.Section().create(comm = scalar_section.comm)
3940+ vector_section.setNumFields(1 )
3941+ vector_section.setFieldComponents(0 , gdim)
3942+ vector_section.setChart(p_start, p_end)
3943+ vector_section.setPermutation(scalar_section.getPermutation())
3944+
3945+ for p in range (p_start, p_end):
3946+ scalar_ndof = scalar_section.getDof(p)
3947+ CHKERR(PetscSectionSetDof(vector_section.sec, p, scalar_ndof* gdim))
3948+ CHKERR(PetscSectionSetFieldDof(vector_section.sec, p, 0 , scalar_ndof* gdim))
3949+ vector_section.setUp()
3950+
39333951 if coordinate_element.family() == " Lagrange" :
3934- # the existing section has the correct numbering but is scalar, expand to gdim
3935- scalar_dm = coordinates.function_space().dm
3936- scalar_section = scalar_dm.getLocalSection()
3937- p_start, p_end = scalar_section.getChart()
3938-
3939- vector_section = PETSc.Section().create(comm = scalar_section.comm)
3940- vector_section.setNumFields(1 )
3941- vector_section.setFieldComponents(0 , gdim)
3942- vector_section.setChart(p_start, p_end)
3943- vector_section.setPermutation(scalar_section.getPermutation())
3944-
3945- for p in range (p_start, p_end):
3946- scalar_ndof = scalar_section.getDof(p)
3947- CHKERR(PetscSectionSetDof(vector_section.sec, p, scalar_ndof* gdim))
3948- CHKERR(PetscSectionSetFieldDof(vector_section.sec, p, 0 , scalar_ndof* gdim))
3949- vector_section.setUp()
3950-
3951- # set coordinates
39523952 # apparently gdim ignored in this call and set explicitly below
39533953 CHKERR(DMSetCoordinateSection(plex.dm, gdim, vector_section.sec))
39543954 plex.setCoordinateDim(gdim)
3955-
39563955 plex.setCoordinatesLocal(coordinates.dat._vec)
39573956 else :
3958- raise NotImplementedError (" Need to do DG" )
3957+ plex.setCellCoordinateDM(plex.getCoordinateDM().clone())
3958+ plex.setCellCoordinateSection(gdim, vector_section) # gdim ignored
3959+ plex.setCoordinateDim(gdim)
3960+ plex.setCellCoordinatesLocal(coordinates.dat._vec)
0 commit comments