Issue Importing XDMF Mesh in FESTIM

Hello,

I generated a 2D mesh for my case study (1 cm Tungsten, 3 cm EUROFER, and 2 cm FLiBe) using DolfinX and exported it into two files: mesh_domains.xdmf and mesh_boundaries.xdmf.

Here is the code I used to load the mesh in FESTIM:
import festim as F
import numpy as np

my_model = F.Simulation()

my_model.mesh = F.MeshFromXDMF(volume_file=“mesh_domains.xdmf”, boundary_file=“mesh_boundaries.xdmf”)

However, I get the following error:

*** Error: Unable to open XDMF file.
*** Reason: XDMF file “mesh_domains.xdmf” does not exist.
*** Where: This error was encountered inside XDMFFile.cpp.
*** Process: 0


*** DOLFIN version: 2019.1.0
*** Git changeset: ec7503cef5f328d4927ba1050abe55c3519ac4e9
*** -------------------------------------------------------------------------

ERROR conda.cli.main_run:execute(127): conda run python /home/ndrea_i_aio/festim/TESI/MODEL_W_EURO_FLiBe failed. (See above for error)

I have confirmed that the files exist in my TESI directory.
Could you please advise how to correctly import XDMF meshes into FESTIM, or suggest what might be causing this error?

Thank you very much for your help.

Best regards,
Andrea

Hi @Andrea! this error message indicates that the file mesh_domains.xdmf isn’t found by your script. The path that is passed to volume_file needs to be the relative path from where you execute your script.

Thanks for the quick answer. I have now fixed the paths and they are correct. However, I am still facing an issue when importing the mesh. This is the error I get:

  File "/home/ndrea_i_aio/festim/TESI/MODEL_W_EURO_FLiBe", line 6, in <module>
    my_model.mesh = F.MeshFromXDMF(volume_file="/home/ndrea_i_aio/festim/TESI/mesh_domains.xdmf", boundary_file="/home/ndrea_i_aio/festim/TESI/mesh_boundaries.xdmf")
                    ~~~~~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/ndrea_i_aio/miniconda3/envs/festim_env_2/lib/python3.13/site-packages/festim/meshing/mesh_from_xdmf.py", line 28, in __init__
    self.define_markers()
    ~~~~~~~~~~~~~~~~~~~^^
  File "/home/ndrea_i_aio/miniconda3/envs/festim_env_2/lib/python3.13/site-packages/festim/meshing/mesh_from_xdmf.py", line 43, in define_markers
    f.XDMFFile(self.boundary_file).read(surface_markers, "f")
    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^^^
RuntimeError: 

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
***     fenics-support@googlegroups.com
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to recognise cell type.
*** Reason:  Unknown value "".
*** Where:   This error was encountered inside XDMFFile.cpp.
*** Process: 0
*** 
*** DOLFIN version: 2019.1.0
*** Git changeset:  ec7503cef5f328d4927ba1050abe55c3519ac4e9
*** -------------------------------------------------------------------------


ERROR conda.cli.main_run:execute(127): `conda run python /home/ndrea_i_aio/festim/TESI/MODEL_W_EURO_FLiBe` failed. (See above for error)

Could this error be due to an incompatibility between FESTIM (Fenics 2019) and DolfinX?

I am also including the code I use to generate the mesh:

import numpy as np
import gmsh
from mpi4py import MPI
import os
from dolfinx.io import gmshio
import matplotlib.pyplot as plt
import pyvista
from dolfinx import plot

# Geometry
wi_tun = 1e-2   # tungsten thickness
wi_euro = 3e-2  # EUROFER thickness
wi_FLiBe = 2e-2 # FLiBe thickness
Le = 1

def create_three_region_mesh():
    width_tung = wi_tun
    width_euro = wi_euro
    width_FLiBe = wi_FLiBe
    length = Le
    gmsh.initialize()
    gmsh.model.add("three_regions")

    # Points
    p1 = gmsh.model.geo.addPoint(0, 0, 0)
    p2 = gmsh.model.geo.addPoint(width_tung, 0, 0)
    p3 = gmsh.model.geo.addPoint(0, length, 0)
    p4 = gmsh.model.geo.addPoint(width_tung, length, 0)

    p5 = gmsh.model.geo.addPoint(width_tung + width_euro, 0, 0)
    p6 = gmsh.model.geo.addPoint(width_tung + width_euro, length, 0)

    p7 = gmsh.model.geo.addPoint(width_tung + width_euro + width_FLiBe, 0, 0)
    p8 = gmsh.model.geo.addPoint(width_tung + width_euro + width_FLiBe, length, 0)

    # Lines Tungsten
    l1 = gmsh.model.geo.addLine(p1, p2)
    l2 = gmsh.model.geo.addLine(p3, p1)
    l3 = gmsh.model.geo.addLine(p4, p3)
    l4 = gmsh.model.geo.addLine(p2, p4)

    # Lines EUROFER
    l5 = gmsh.model.geo.addLine(p2, p5)
    l6 = gmsh.model.geo.addLine(p5, p6)
    l7 = gmsh.model.geo.addLine(p6, p4)

    # Lines FLiBe
    l8 = gmsh.model.geo.addLine(p5, p7)
    l9 = gmsh.model.geo.addLine(p7, p8)
    l10 = gmsh.model.geo.addLine(p8, p6)
    l11 = gmsh.model.geo.addLine(p6, p5)

    # Curve Loops
    tung_loop = gmsh.model.geo.addCurveLoop([l1, l4, l3, l2])
    euro_loop = gmsh.model.geo.addCurveLoop([l5, l6, l7, -l4])
    flibe_loop = gmsh.model.geo.addCurveLoop([l8, l9, l10, l11])

    # Surfaces
    tung_surf = gmsh.model.geo.addPlaneSurface([tung_loop])
    euro_surf = gmsh.model.geo.addPlaneSurface([euro_loop])
    flibe_surf = gmsh.model.geo.addPlaneSurface([flibe_loop])
    gmsh.model.geo.synchronize()

    # Mesh size
    gmsh.option.setNumber("Mesh.MeshSizeMin", 0.02)
    gmsh.option.setNumber("Mesh.MeshSizeMax", 0.02)

    # Physical groups
    gmsh.model.addPhysicalGroup(2, [tung_surf], 1)
    gmsh.model.setPhysicalName(2, 1, "Tungsten")
    gmsh.model.addPhysicalGroup(2, [euro_surf], 2)
    gmsh.model.setPhysicalName(2, 2, "EuroFer")
    gmsh.model.addPhysicalGroup(2, [flibe_surf], 3)
    gmsh.model.setPhysicalName(2, 3, "FLiBe")

    gmsh.model.addPhysicalGroup(1, [l2], 4)   # Left
    gmsh.model.addPhysicalGroup(1, [l9], 5)   # Right
    gmsh.model.addPhysicalGroup(1, [l3, l7, l10], 6)  # Top
    gmsh.model.addPhysicalGroup(1, [l1, l5, l8], 7)   # Bottom
    gmsh.model.addPhysicalGroup(1, [l4], 8)   # Interface W-Euro
    gmsh.model.addPhysicalGroup(1, [l11], 9)  # Interface Euro-FLiBe

    gmsh.model.mesh.generate(2)

    # Conversion to dolfinx
    mesh_2d, cell_tags, facet_tags = gmshio.model_to_mesh(
        gmsh.model, MPI.COMM_WORLD, 0, gdim=2
    )
    gmsh.finalize()
    total_width = width_tung + width_euro + width_FLiBe
    return mesh_2d, total_width, width_tung, width_euro, width_FLiBe, cell_tags, facet_tags

mesh_2d, total_width, width_tung, width_euro, width_FLiBe, cell_tags, facet_tags = create_three_region_mesh()

# Export
from dolfinx.io import XDMFFile
output_dir = "TESI"
os.makedirs(output_dir, exist_ok=True)

domains_file = os.path.join(output_dir, "mesh_domains.xdmf")
with XDMFFile(MPI.COMM_WORLD, domains_file, "w") as xdmf:
    xdmf.write_mesh(mesh_2d)
    xdmf.write_meshtags(cell_tags, mesh_2d.geometry, geometry_xpath="/Mesh/Geometry")

boundaries_file = os.path.join(output_dir, "mesh_boundaries.xdmf")
with XDMFFile(MPI.COMM_WORLD, boundaries_file, "w") as xdmf:
    xdmf.write_mesh(mesh_2d)
    xdmf.write_meshtags(facet_tags, mesh_2d.geometry, geometry_xpath="/Mesh/Geometry")

Thanks in advance!

I don’t think mixing FESTIM1 and dolfinx is a good idea. You may want to use FESTIM2 - which relies on dolfinx.

Our FESTIM2 tutorial is building up here

Got it, thanks for the clarification! Sorry for the trouble.
I’ll use FESTIM1 with SALOME.