messy mesh obtained when using spline/bsline

[yujingll]

Hi, I have problems with generate desired mesh when using spline/bspline. When using spline/bspline one gives very messy mesh grid:

while the desired mesh from python using dolfinx and gmsh is like

Here attached the code for production of the messy mesh in Julia and desired mesh in python

using Ferrite
using FerriteGmsh
using FerriteViz

gmsh.initialize()
gdim = 2
p1 = gmsh.model.occ.addPoint(0.0, 0.0, 0, 0)
p2 = gmsh.model.occ.addPoint(1.0, 0.0, 0, 0)
p3 = gmsh.model.occ.addPoint(1.0, 0.5, 0, 0)
p4 = gmsh.model.occ.addPoint(1.0, 1.0, 0, 0)
s1 = gmsh.model.occ.addBSpline([p1, p2, p3, p4])

p2 = gmsh.model.occ.addPoint(0.0, 1.0, 0, 0)
p3 = gmsh.model.occ.addPoint(0.5, 1.0, 0, 0)
s2 = gmsh.model.occ.addSpline([p4, p3, p2, p1])

ll = gmsh.model.occ.addCurveLoop([s1, s2])
pl = gmsh.model.occ.addPlaneSurface([ll])

gmsh.model.occ.synchronize()
gmsh.model.addPhysicalGroup(2, [1], 1)

meshSize = 0.1
gmsh.option.setNumber("Mesh.CharacteristicLengthMin",meshSize)
gmsh.option.setNumber("Mesh.CharacteristicLengthMax",meshSize)
gmsh.model.mesh.generate(gdim)

dim = Int64(gmsh.model.getDimension())
facedim = dim - 1 

nodes = tonodes()
elements, gmsh_elementidx = toelements(dim)
cellsets = tocellsets(dim, gmsh_elementidx)

domaincellset = cellsets["1"]
elements = elements[collect(domaincellset)]

boundarydict = toboundary(facedim)
facesets = tofacesets(boundarydict, elements)
gmsh.finalize()
grid = Grid(elements, nodes, facesets=facesets, cellsets=cellsets)
import WGLMakie 

FerriteViz.wireframe(grid,markersize=10,strokewidth=2)

import gmsh
from dolfinx.io import gmshio
import pyvista
from mpi4py import MPI
from dolfinx import fem
from dolfinx import plot

gmsh.initialize()
p1 = gmsh.model.occ.addPoint(0.0, 0.0, 0, 0)
p2 = gmsh.model.occ.addPoint(1.0, 0.0, 0, 0)
p3 = gmsh.model.occ.addPoint(1.0, 0.5, 0, 0)
p4 = gmsh.model.occ.addPoint(1.0, 1.0, 0, 0)
s1 = gmsh.model.occ.addBSpline([p1, p2, p3, p4])

p2 = gmsh.model.occ.addPoint(0.0, 1.0, 0, 0)
p3 = gmsh.model.occ.addPoint(0.5, 1.0, 0, 0)
s2 = gmsh.model.occ.addSpline([p4, p3, p2, p1])

ll = gmsh.model.occ.addCurveLoop([s1, s2])
pl = gmsh.model.occ.addPlaneSurface([ll])

gmsh.model.occ.synchronize()
gmsh.model.addPhysicalGroup(2, [1], 1)

meshSize = 0.1
gdim = 2
gmsh.option.setNumber("Mesh.CharacteristicLengthMin",meshSize)
gmsh.option.setNumber("Mesh.CharacteristicLengthMax",meshSize)
gmsh.model.mesh.generate(gdim)

gmsh_model_rank = 0
mesh_comm = MPI.COMM_WORLD
domain, cell_markers, facet_markers = gmshio.model_to_mesh(gmsh.model, mesh_comm, gmsh_model_rank, gdim=gdim)

V = fem.functionspace(domain, ("Lagrange", 1))
cells, types, x = plot.vtk_mesh(V)
grid = pyvista.UnstructuredGrid(cells, types, x)
plotter = pyvista.Plotter()
plotter.add_mesh(grid, show_edges=True)
plotter.show()
gmsh.finalize()

[koehlerson]

Hello,

seems to be the same issue as https://github.com/Ferrite-FEM/FerriteGmsh.jl/issues/20 I thought that I fixed it with the saveall flag but it seems to be still the same. There is some difference in how the API behaves when reading from file and when reading from memory, see:

using Ferrite
using FerriteGmsh
using FerriteViz

gmsh.initialize()
gdim = 2
p1 = gmsh.model.occ.addPoint(0.0, 0.0, 0, 0)
p2 = gmsh.model.occ.addPoint(1.0, 0.0, 0, 0)
p3 = gmsh.model.occ.addPoint(1.0, 0.5, 0, 0)
p4 = gmsh.model.occ.addPoint(1.0, 1.0, 0, 0)
s1 = gmsh.model.occ.addBSpline([p1, p2, p3, p4])

p2 = gmsh.model.occ.addPoint(0.0, 1.0, 0, 0)
p3 = gmsh.model.occ.addPoint(0.5, 1.0, 0, 0)
s2 = gmsh.model.occ.addSpline([p4, p3, p2, p1])

ll = gmsh.model.occ.addCurveLoop([s1, s2])
pl = gmsh.model.occ.addPlaneSurface([ll])

gmsh.model.occ.synchronize()
gmsh.model.addPhysicalGroup(2, [1], 1)

meshSize = 0.1
gmsh.option.setNumber("Mesh.CharacteristicLengthMin",meshSize)
gmsh.option.setNumber("Mesh.CharacteristicLengthMax",meshSize)
gmsh.model.mesh.generate(gdim)
gmsh.write("test.msh")

dim = Int64(gmsh.model.getDimension())
facedim = dim - 1 

nodes = tonodes()
elements, gmsh_elementidx = toelements(dim)
cellsets = tocellsets(dim, gmsh_elementidx)

domaincellset = cellsets["1"]
elements = elements[collect(domaincellset)]

boundarydict = toboundary(facedim)
facesets = tofacesets(boundarydict, elements)
gmsh.finalize()
grid_bad = Grid(elements, nodes, facesets=facesets, cellsets=cellsets)
grid_good = togrid("test.msh")
import GLMakie 

FerriteViz.wireframe(grid_good,markersize=10,strokewidth=2)

gives me:

[koehlerson]

mh yeah it's related to the SaveAll flag. In this snippet I flipped the flag and obtained the correct mesh for grid_bad:

using Ferrite
using FerriteGmsh
using FerriteViz

gmsh.initialize()
gdim = 2
saveall_flag = Bool(gmsh.option.getNumber("Mesh.SaveAll"))
if !saveall_flag
    gmsh.option.setNumber("Mesh.SaveAll",1)
end
p1 = gmsh.model.occ.addPoint(0.0, 0.0, 0, 0)
p2 = gmsh.model.occ.addPoint(1.0, 0.0, 0, 0)
p3 = gmsh.model.occ.addPoint(1.0, 0.5, 0, 0)
p4 = gmsh.model.occ.addPoint(1.0, 1.0, 0, 0)
s1 = gmsh.model.occ.addBSpline([p1, p2, p3, p4])

p2 = gmsh.model.occ.addPoint(0.0, 1.0, 0, 0)
p3 = gmsh.model.occ.addPoint(0.5, 1.0, 0, 0)
s2 = gmsh.model.occ.addSpline([p4, p3, p2, p1])

ll = gmsh.model.occ.addCurveLoop([s1, s2])
pl = gmsh.model.occ.addPlaneSurface([ll])

gmsh.model.occ.synchronize()
gmsh.model.addPhysicalGroup(2, [1], 1)

meshSize = 0.1
gmsh.option.setNumber("Mesh.CharacteristicLengthMin",meshSize)
gmsh.option.setNumber("Mesh.CharacteristicLengthMax",meshSize)
gmsh.model.mesh.generate(gdim)
gmsh.write("test.msh")

dim = Int64(gmsh.model.getDimension())
facedim = dim - 1 

nodes = tonodes()
elements, gmsh_elementidx = toelements(dim)
cellsets = tocellsets(dim, gmsh_elementidx)

domaincellset = cellsets["1"]
elements = elements[collect(domaincellset)]

boundarydict = toboundary(facedim)
facesets = tofacesets(boundarydict, elements)
gmsh.finalize()
grid_bad = Grid(elements, nodes, facesets=facesets, cellsets=cellsets)
grid_good = togrid("test.msh")
import GLMakie 

FerriteViz.wireframe(grid_bad,markersize=10,strokewidth=2)

[koehlerson]

I don't really understand why this changes the API's behavior

[koehlerson]

If you translate from memory instead of from disk and enable the SaveAll flag, then it seems that the control nodes of the BSpline are exported as well. I'm quite unsure how to handle it. My advice is to follow as closely as possible the togrid operations (the auxillary ones here: https://github.com/Ferrite-FEM/FerriteGmsh.jl/blob/master/src/FerriteGmsh.jl#L195C1-L201C39) and I need to update the README and docs