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)