Internal mesh vs gmsh

Hi,
I’m new to fenicx. I would like to use it to simulate heat transport under laser irradiation at the surface.
I did a test with internal mesh and it seems to work. However, when I try to use gmsh to add more materials I start to get artifacts even with simple cube.
I believe it is related to the mesh topology but could not find the reason.
I will really appreciate some help.

Here are the results and the 2 codes:


Internal mesh:

import numpy as np
import ufl
from dolfinx.fem.petsc import LinearProblem
from petsc4py import PETSc
from mpi4py import MPI
from dolfinx import fem, mesh, plot
from dolfinx.io import XDMFFile
from petsc4py.PETSc import ScalarType

L = 100.0
H = 100.0
dL = 30
dH = 30
domain = mesh.create_box(MPI.COMM_WORLD,[[0.0,0.0,0.0], [L, L, H]], [dL, dL, dH], mesh.CellType.tetrahedron)
V = fem.FunctionSpace(domain, ("CG", 1))

x = ufl.SpatialCoordinate(domain)
rpump = 25.0
q =  ((H-x[2])/ H) * x[0] * (x[0]-L) * x[1] * (x[1]-L) *  (2 / (np.pi * rpump * rpump) )  * ufl.exp( - 2 * ((x[0]-L/2.0)**2 + (x[1]-L/2.0)**2) / (rpump * rpump ))

omega = 1e2
c=0.4
k=20
eta = 2*np.pi*1j*c/k

u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)


a = ( ufl.inner(ufl.grad(u), ufl.grad(v))  + eta * omega * ufl.inner(u,v) ) * ufl.dx 
L =  (omega / k) * ufl.inner(q, v) * ufl.ds

problem = LinearProblem(a, L, petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
uh = problem.solve()

import dolfinx.io
with dolfinx.io.XDMFFile(domain.comm, "outputreg.xdmf", "w") as xdmf:
    xdmf.write_mesh(domain)
    xdmf.write_function(uh)

Using gmsh

import numpy as np
import ufl
from dolfinx.fem.petsc import LinearProblem
from petsc4py import PETSc
from mpi4py import MPI
from dolfinx import fem, mesh, plot
from dolfinx.io import XDMFFile
from petsc4py.PETSc import ScalarType
import gmsh
import sys

L = 100.0
H = 100.0
dL = 40
dH = 20
interfaces = [H]
lc = 8
gmsh.initialize()
proc = MPI.COMM_WORLD.rank 

if proc == 0:

    gmsh.model.add("Model")

    gmsh.model.occ.addBox(0, 0, 0, L, L, H, 1)
    gmsh.model.occ.synchronize()
    gmsh.model.addPhysicalGroup(3, [1], 1)
    gmsh.option.setNumber("Mesh.CharacteristicLengthMin",2)
    gmsh.option.setNumber("Mesh.CharacteristicLengthMax",5)

    #gmsh.option.setNumber("Mesh.RecombineAll", 1)
    gmsh.option.setNumber("Mesh.RecombinationAlgorithm", 1)
    gmsh.option.setNumber("Mesh.Algorithm", 8)
    gmsh.option.setNumber("Mesh.Algorithm3D", 1)
    #gmsh.option.setNumber("Mesh.SubdivisionAlgorithm", 2)

    gmsh.model.mesh.generate(3)

    gmsh.option.setNumber("Mesh.SaveAll", 0)
    gmsh.model.geo.synchronize()
    gmsh.write("t2.msh")

    if '-nopopup' not in sys.argv:
        gmsh.fltk.run()
gmsh.finalize()
from dolfinx.io.gmshio import read_from_msh
model_rank = 0
domain, cell_tags, facet_tags = read_from_msh("t2.msh", MPI.COMM_WORLD, model_rank, gdim=3)
V = fem.FunctionSpace(domain, ("CG", 1))
x = ufl.SpatialCoordinate(domain)
rpump = 25.0
q =  ((H-x[2])/H) * x[0] * (x[0]-L) * x[1] * (x[1]-L) * (2 / (np.pi * rpump * rpump) )  * ufl.exp( - 2 * ((x[0]-L/2.0)**2 + (x[1]-L/2.0)**2) / (rpump * rpump ))

omega = 1e2
c=0.4
k=20
eta = 2*np.pi*1j*c/k

u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)

V = fem.FunctionSpace(domain, ("CG", 1))
x = ufl.SpatialCoordinate(domain)
rpump = 25.0
q =  ((H-x[2])/H) * x[0] * (x[0]-L) * x[1] * (x[1]-L) * (2 / (np.pi * rpump * rpump) )  * ufl.exp( - 2 * ((x[0]-L/2.0)**2 + (x[1]-L/2.0)**2) / (rpump * rpump ))

omega = 1e3
c=0.4
k=20
eta = 2*np.pi*1j*c/k

u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)
a = ( ufl.inner(ufl.grad(u), ufl.grad(v))  + eta * omega * ufl.inner(u,v) ) * ufl.dx 
L =  (omega / k) * ufl.inner(q, v) * ufl.ds

problem = LinearProblem(a, L,  petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
uh = problem.solve()

import dolfinx.io
with dolfinx.io.XDMFFile(domain.comm, "outputocc.xdmf", "w") as xdmf:
           xdmf.write_mesh(domain)
           xdmf.write_function(uh)

First of all, please format each code with 3x`, such as

```python
import dolfinx
#add more code here
```

Done, I’m sorry about it.