diff --git a/src/nos/data/helmholtz/solver/util.py b/src/nos/data/helmholtz/solver/util.py
index cd4606d3..14f47e5e 100644
--- a/src/nos/data/helmholtz/solver/util.py
+++ b/src/nos/data/helmholtz/solver/util.py
@@ -8,7 +8,7 @@
 
 def get_mesh(
     mesh_file: pathlib.Path, comm: MPI.Comm = MPI.COMM_SELF
-) -> Tuple[dolfinx.mesh.Mesh, dolfinx.mesh.MeshTags, dolfinx.mesh.MeshTags]:
+) -> (Tuple)[dolfinx.mesh.Mesh, dolfinx.mesh.MeshTags, dolfinx.mesh.MeshTags]:
     """Loads mesh from a msh file and returns it as a dolfinx mesh. Initializes mesh on rank 0 of the given
     communicator.
 
@@ -19,5 +19,8 @@ def get_mesh(
     Returns:
         Tuple containing the dolfinx mesh, mesh-tags and boundary tags.
     """
+    gmsh.initialize()
     gmsh.open(mesh_file)
-    return dolfinx.io.gmshio.model_to_mesh(gmsh.model, comm, 0)
+    mesh, cell_tags, boundary_tags = dolfinx.io.gmshio.model_to_mesh(gmsh.model, comm, 0)
+    gmsh.finalize()
+    return mesh, cell_tags, boundary_tags
diff --git a/test/data/helmholtz/solver/__init__.py b/test/data/helmholtz/solver/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/test/data/helmholtz/solver/test_solver_util.py b/test/data/helmholtz/solver/test_solver_util.py
new file mode 100644
index 00000000..54769965
--- /dev/null
+++ b/test/data/helmholtz/solver/test_solver_util.py
@@ -0,0 +1,52 @@
+import pathlib
+import tempfile
+
+import dolfinx
+import gmsh
+from mpi4py import MPI
+
+from nos.data.helmholtz.solver import get_mesh
+
+
+def test_get_mesh():
+    # create mesh
+    gmsh.initialize()
+    gmsh.model.add("test_get_mesh")
+
+    rect = gmsh.model.occ.add_rectangle(0.0, 0.0, 0.0, 1.0, 1.0)
+    circ = gmsh.model.occ.add_disk(0.0, 0.0, 0.0, 0.5, 0.5)
+
+    objs, _ = gmsh.model.occ.fragment([(2, rect)], [(2, circ)])
+
+    gmsh.model.occ.synchronize()
+
+    surfs = gmsh.model.occ.get_entities(2)
+    surfs = [surf[1] for surf in surfs]
+    lines = gmsh.model.occ.get_entities(1)
+    lines = [line[1] for line in lines]
+
+    for i, surf in enumerate(surfs):
+        gmsh.model.add_physical_group(2, [surf], i)
+    for i, line in enumerate(lines):
+        gmsh.model.add_physical_group(1, [line], i)
+
+    gmsh.model.occ.synchronize()
+    gmsh.model.mesh.generate(2)
+
+    with tempfile.TemporaryDirectory() as tmp:
+        # write file
+        file_path = pathlib.Path(tmp).joinpath("test_mesh.msh")
+        gmsh.write(file_path)
+        gmsh.finalize()
+
+        # load
+        mesh, cell_tags, boundary_tags = get_mesh(file_path, MPI.COMM_WORLD)
+
+    # test types
+    assert isinstance(mesh, dolfinx.mesh.Mesh)
+    assert isinstance(cell_tags, dolfinx.mesh.MeshTags)
+    assert isinstance(boundary_tags, dolfinx.mesh.MeshTags)
+
+    # test size
+    assert len(lines) == len(set(boundary_tags.values))
+    assert len(surfs) == len(set(cell_tags.values))