Problem with fem.form()

dolfinx
1

Hello everyone,

I’m working with a 2D mesh generated by Gmsh and encountered an error while writing the following code. Any help would be greatly appreciated.

element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.triangle, 1)
e = basix.ufl_wrapper.BasixElement(element, 3)
V = fem.FunctionSpace(domain, e)
Et = ufl.TrialFunction(V)
W = ufl.TestFunction(V)
a = fem.form(ufl.inner(ufl.curl(Et), ufl.curl(W)) * ufl.dx)

Error:

Expecting pulled back expression with shape '(2,)', got '(3,)'
ERROR:UFL:Expecting pulled back expression with shape '(2,)', got '(3,)'
---------------------------------------------------------------------------
UFLException                              Traceback (most recent call last)
foo.ipynb Cell 3 line 9
      6 Et = ufl.TrialFunction(V)
      7 W = ufl.TestFunction(V)
----> 9 a=fem.form(ufl.inner(ufl.curl(Et), ufl.curl(W)) * ufl.dx)
2

Why are you sending in 3 as the second argument? This is the geometric dimension of the mesh, which is 2 if you are solving a non-manifold 2D problem

3

Thank you for your response. When I send 2 as the second argument, I encountered the following error:

Non-matching cell of finite element and domain.
ERROR:UFL:Non-matching cell of finite element and domain.

Possibly, this is due to domain.ufl_cell() = Cell('triangle', 3). I generate the domain using the following code:

from dolfinx.io.gmshio import model_to_mesh
import gmsh

R = 500

gmsh.initialize()
gmsh.model.add("model_1")
   
meshsize = 18     
x0, y0 = 0, 0

p1 = gmsh.model.geo.addPoint(x0  , y0, 0, meshsize)    
p2 = gmsh.model.geo.addPoint(x0 + R, y0, 0, meshsize)
p3 = gmsh.model.geo.addPoint(x0 - R, y0, 0, meshsize)   

arc1 = gmsh.model.geo.addCircleArc(p2, p1, p3)   # arch top
arc2 = gmsh.model.geo.addCircleArc(p3, p1, p2)   # arch bottom
c_out = gmsh.model.geo.addCurveLoop([arc1, arc2])

s_bg  =  gmsh.model.geo.addPlaneSurface([c_out])        
gmsh.model.addPhysicalGroup(2, [s_bg], tag=0)

gmsh.model.geo.synchronize()
gdim = 2

gmsh.model.mesh.generate(gdim)

model_rank = 0
domain, cell_tags, _ = model_to_mesh(gmsh.model, MPI.COMM_WORLD, model_rank)

I’m curious about why domain.geometry.x.shape yields (2975, 3) . Does this imply that the domain is in 3D? If so, how can I modify it to be in 2D?

4

Do you want to solve on a manifold? If you simply want to solve on a plane 2D geometry, you need to change

to
domain, cell_tags, _ = model_to_mesh(gmsh.model, MPI.COMM_WORLD, model_rank,gdim=gdim)

1 Like
5

Thank you so much! This really helped me.