I would strongly suggest to load the boundary markers used for boundary conditions from file (as you can define them in your mesh generation software). There are plenty of posts describing such a pipeline on the forum.
Thanks! can I try the following option:
import dolfinx
import numpy as np
from mpi4py import MPI
from dolfinx.cpp.mesh import CellType
import meshio
def create_mesh(mesh, cell_type, prune_z=False):
cells = mesh.get_cells_type(cell_type)
# This line needs to get cell data from what ever the markers have been called in the vtk file.
# Expore this by printing mesh.cell_data
cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
if prune_z:
points = mesh.points[:,:2]
else:
points = mesh.points
out_mesh = meshio.Mesh(points=points, cells={cell_type: cells})
cell_data={"name_to_read":[cell_data]})
return out_mesh
in_mesh = meshio.read("sphere.vtk")
out_mesh = create_mesh(in_mesh, "tetra")
meshio.write("sphere.xdmf", out_mesh)
would this preserve the name tag of the mesh? Because I do not use the gmsh for the mesh creation, I use another tool to create a CAD and mesh. And I can save the mesh either in .cgns or .vtk format.
As written here, you need to explore mesh.cell_data to find the appropriate markers.
Please note that you need to store boundary markers in a separate mesh, as shown in:
Thank you! I will try that.
I have tried two options (one with .cgns and .vtk):
import meshio
import numpy as np
def create_mesh(mesh: meshio.Mesh, cell_type: str, data_name: str = "name_to_read",
prune_z: bool = False):
cells = mesh.get_cells_type(cell_type)
cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
points = mesh.points[:, :2] if prune_z else mesh.points
out_mesh = meshio.Mesh(points=points, cells={cell_type: cells}, cell_data={
data_name: [cell_data]})
return out_mesh
# import and convert msh file to xdmf
msh = meshio.read('mesh.cgns')
line_mesh = create_mesh(msh, "line", data_name="Boundary", prune_z=True)
meshio.write("mf.xdmf", line_mesh)
triangle_mesh = create_mesh(msh, "tetra", prune_z=True)
meshio.write("mesh.xdmf", triangle_mesh)
Then I am getting the following error:
Traceback (most recent call last):
File "mesh.py", line 17, in <module>
line_mesh = create_mesh(msh, "line", data_name="Boundary", prune_z=True)
File "mesh.py", line 8, in create_mesh
cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
File "/python-venv/lib/python3.8/site-packages/meshio/_mesh.py", line 228, in get_cell_data
[d for c, d in zip(self.cells, self.cell_data[name]) if c.type == cell_type]
KeyError: 'gmsh:physical'
and I also tried another option to get .msh file using meshio-convert mesh.cgns mesh.gmsh
Traceback (most recent call last):
File "/python-venv/bin/meshio-convert", line 8, in <module>
sys.exit(convert())
File "/python-venv/lib/python3.8/site-packages/meshio/_cli/_convert.py", line 15, in convert
mesh = read(args.infile, file_format=args.input_format)
File "/python-venv/lib/python3.8/site-packages/meshio/_helpers.py", line 69, in read
return reader_map[file_format](filename)
File "/python-venv/lib/python3.8/site-packages/meshio/cgns/_cgns.py", line 16, in read
f = h5py.File(filename, "r")
File "/python-venv/lib/python3.8/site-packages/h5py/_hl/files.py", line 507, in __init__
fid = make_fid(name, mode, userblock_size, fapl, fcpl, swmr=swmr)
File "/python-venv/lib/python3.8/site-packages/h5py/_hl/files.py", line 220, in make_fid
fid = h5f.open(name, flags, fapl=fapl)
File "h5py/_objects.pyx", line 54, in h5py._objects.with_phil.wrapper
File "h5py/_objects.pyx", line 55, in h5py._objects.with_phil.wrapper
File "h5py/h5f.pyx", line 106, in h5py.h5f.open
OSError: Unable to open file (file signature not found)
Could you please tell me how to solve this issue?
For the mesh creation, I am using the ANSYS, and I believe that name tags will be saved in .cgns file. Or is there any other format that I can export from ANSYS and that can be compatible with FEniCS-X mesh conversion.
Again, you should look at what mesh.cell_data contains. There should be a dictionary there with markers, you need to use the correct key. “Gmsh:physical” is the key used by Gmsh to indicate markers.
I would submit an issue at:
For your second issue.
Meshio supports ANSYS msh, according to the readme at: GitHub - nschloe/meshio: 🕸 input/output for many mesh formats
Thanks! that is the thing, I do not know how to get those tags (where to see or how to get those names). In ANSYS I create the name tags during the CAD modeling and later will be exported to during mesh creation. Do you mean, I need to replace the name in here Gmsh:physical with, for example, I create a name surfaces with inlet, outlet, walls and sphere during the meshing process.
Simply do print(mesh.cell_data). This will give you what data meshio has been able to read in.
it is giving some blank {}, I will check one more time, in ANSYS if I am exporting the mesh with name tags. Because when I try to open the file in gmsh, I do not see any names either.
I have just used
meshio-info mesh.vtk
<meshio mesh object>
Number of points: 369976
Number of cells:
tetra: 2127126
triangle: 52414
Would this info be useful?
Or from the beginning of the thread of this problem, it works, the only problem is defining the BC for the obstacle. Is there any example of how to define the BC for obstacle when there is no naming tag is in the mesh?
For example, in FEniCS 2019, the following was working:
class Sphere(SubDomain):
def inside(self, x, on_boundary):
result = on_boundary\
and x[0] > 0.8 - DOLFIN_EPS and x[0] < 2.0 + DOLFIN_EPS\
and x[1] > 0.2 - DOLFIN_EPS and x[1] < 1.8 + DOLFIN_EPS\
and x[2] > 0.05 - DOLFIN_EPS and x[2] < 0.3 + DOLFIN_EPS
return result
But when we try the same approach in FEniCS-X it does not work:
def Sphere(x):
result = np.logical_and(x[0] > 2.0 - DOLFIN_EPS, x[0] < 4.0 - DOLFIN_EPS)
result2 = np.logical_and(x[1] > 2.0 - DOLFIN_EPS, x[1] < 4.0 - DOLFIN_EPS)
result3 = np.logical_and(x[2] > 2.0 - DOLFIN_EPS, x[2] < 4.0 - DOLFIN_EPS)
return np.logical_and(np.logical_and(result, result2), result3)
In FEniCSx we do not have a DOLFIN_EPS defined. You can define your own tolerance (for instance eps=1e-14).
Thanks! changed it and running the simulation.
Before, it was like that DOLFIN_EPS=0.000001 now I have changed to eps=1e-14
DOLFIN_EPS is not the issue, I think the issue comes from how we define the logical condition or something it has to do with FEniCSx.
As you have not supplied the mesh, i cant really help you any further.
Please supply the mesh and a minimal code to reproduce the issue
Many thanks! Here is the code and mesh
You need to upload the mesh .h5 file as well.
Sorry my bad, it is now uploaded.
I would suggest using locate_entities_boundary and locate_dofs_topological to create your bc, as follows:
facets = dolfinx.mesh.locate_entities_boundary(
mesh, mesh.topology.dim-1, Sphere)
mt = dolfinx.mesh.MeshTags(mesh, mesh.topology.dim-1,
facets, np.full(len(facets), 1, dtype=np.int32))
with XDMFFile(mesh.comm, "mt.xdmf", "w") as xdmf:
xdmf.write_mesh(mesh)
xdmf.write_meshtags(mt)
sphere_boundary_dofs = fem.locate_dofs_topological(
V, mesh.topology.dim-1, facets)
Here you can inspect mt.xdmf to see that the correct boundary is tagged.
However, as your problem requires more than 16 GB RAM I cannot use my laptop to verify the solution, so that is up to you.
Thanks for trying out. If possible, could you please try with a few time (longer time would require lots of memory) steps? I have uploaded a small mesh (1m W(y-axis), 1m H(z-axis), 1.5 L(x-axis), and sphere radius is .1 m and located at .5m (from x,y,z-axis)) (62799 nodes and 344194 elements). Since it is an incompressible flow, just checking an initial few time steps solution is enough to find an error.
Because, with the tag, I am not sure even sure if I am properly exporting those properties along with mesh. But anyway I will try in parallel as you suggested.
It is not the time steps that makes it memory bound, it is the fact that the original mesh you posted has 2.1 million dofs in the velocity space. As I already showed you with the code above, you will get the correct tagged sphere using locate_entities_boundary.
You can simply test this yourself with a few time steps.