Input 2nd order quadrilateral mesh (.msh) to doflinx

Hello,

Can I generate custom quadrilateral mesh of second order in dolfinx (using nodes and connectivity) ?

import numpy as np
import meshio

# Define nodal coordinates (8-node second-order quad elements)
points = np.array([
    [0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [2.0, 0.0, 0.0],  # Bottom row (Corner nodes)
    [0.5, 0.0, 0.0], [1.5, 0.0, 0.0],  # Mid-side nodes on bottom row

    [0.0, 1.0, 0.0], [1.0, 1.0, 0.0], [2.0, 1.0, 0.0],  # Middle row (Corner nodes)
    [0.5, 1.0, 0.0], [1.5, 1.0, 0.0]  # Mid-side nodes in the middle row
])

# Define second-order quadrilateral connectivity (8 nodes per quad)
cells = [np.array([
    [0, 1, 6, 5, 3, 8, 7, 4],  # First quadrilateral element
    [1, 2, 7, 6, 4, 9, 8, 5]   # Second quadrilateral element
])]
el= basix.ufl.element("CG", "quadrilateral", 2, shape=(3, ))
domain = ufl.Mesh(el)
# Create second-order quadrilateral mesh
quad_mesh = create_mesh(MPI.COMM_WORLD, points, cells, domain)

The above code restarts the kernel. How can I generate the 2nd order quad mesh. I need the 2nd order to use curved elements.
mesh, subdomains, boundaries = gmshio.read_from_msh("2ndorder_quad.msh", MPI.COMM_WORLD,0, gdim=3) does not take quadratic elements as input.

Thanks

Comparing to: Mesh generation — FEniCS Workshop and DefElement I am noticing quite a few inconsistencies:

  • you’re passing the connectivity and the nodes to create_mesh in the wrong order
  • Your quad has 8 nodes, but should have 9
  • your list of nodes misses the entire intermediate row (y=0.5)
  • in your text you indicate 2D, but your nodes and gdim are 3d
  • you’re missing the dtype indicators
  • You’re using MPI.COMM_WORLD, which should be MPI.COMM_SELF

Not sure which of those factors ended up being decisive, but here is a running code:

import basix
import ufl
import numpy as np
from mpi4py import MPI
import dolfinx

# Define nodal coordinates (8-node second-order quad elements)
points = np.array([
    [0.0, 0.0], #0 
    [1.0, 0.0], #1
    [2.0, 0.0], #2
    [0.5, 0.0], #3
    [1.5, 0.0], #4
    [0.0, 1.0], #5
    [1.0, 1.0], #6
    [2.0, 1.0], #7
    [0.5, 1.0], #8
    [1.5, 1.0], #9
    [0.0, 0.5], #10
    [1.0, 0.5], #11
    [2.0, 0.5], #12
    [0.5, 0.5], #13
    [1.5, 0.5], #14
], dtype=np.float64)

# Visualization of the ordering:
# 5  8  6  9  7
# 10 13 11 14 12
# 0  3  1  4  2

connectivity = np.array([
    [0, 1, 5, 6, 3, 10, 11, 8, 13],  # First quadrilateral element
    [1, 2, 6, 7, 4, 11, 12, 9, 14]   # Second quadrilateral element
], dtype=np.int64)

el= basix.ufl.element("CG", "quadrilateral", 2, shape=(2, ))
domain_element = ufl.Mesh(el)

# Create second-order quadrilateral mesh
quad_mesh = dolfinx.mesh.create_mesh(MPI.COMM_SELF, connectivity, points, domain_element)
1 Like

See for instance
Https://jsdokken.com/FEniCS23-tutorial/src/mesh_generation.html#higher-order-meshes

If the mesh is made quadratic in gmsh, it will be read in as a quadratic mesh.

1 Like