Incomplete meshtags when using write_mesh_input_order in adios4dolfinx

General
1

Hello all,

I have been using adios4dolfinx for checkpointing function and meshtag data as described in this series of tutorials adios4dolfinx - A framework for checkpointing in dolfinx. It’s been working like a charm !

However, I would like to keep the mesh ordering across storing/reading operations. If I understood correctly, this requires using adios4dolfinx.write_mesh_input_order() along with adios4dolfinx.write_function_on_input_mesh(), at the cost of some computation overhead. I managed to get this working for function data, but it’s failing for meshtag data.

I put together a simpler version of this tutorial to test meshtag checkpointing with original mesh order in a single-thread context:

import dolfinx
import numpy as np
import adios4dolfinx

from mpi4py import MPI
from pathlib import Path

filename = Path("meshtags_mwe.bp")
mesh = dolfinx.mesh.create_unit_square(MPI.COMM_SELF, 10, 10)

meshtags = {}
adios4dolfinx.write_mesh_input_order(filename, mesh)
for i in range(mesh.topology.dim + 1):
    mesh.topology.create_entities(i)
    mesh.topology.create_connectivity(i, mesh.topology.dim)

    entity_map = mesh.topology.index_map(i)
    entities = np.arange(entity_map.size_local, dtype=np.int32)
    values = np.arange(entity_map.size_local, dtype=np.int32)
    meshtags[i] = dolfinx.mesh.meshtags(mesh, i, entities, values)
    adios4dolfinx.write_meshtags(filename, mesh, meshtags[i], meshtag_name=f"meshtags_{i}")

read_mesh = adios4dolfinx.read_mesh(filename, MPI.COMM_SELF)
for i in range(read_mesh.topology.dim + 1):
    read_meshtags = adios4dolfinx.read_meshtags(filename, read_mesh, meshtag_name=f"meshtags_{i}")

    # sort and verify that all tags are present
    np.testing.assert_allclose(np.sort(read_meshtags.values), np.sort(meshtags[i].values), err_msg=f"Meshtags values for dimension {i} are incomplete")
    np.testing.assert_allclose(np.sort(read_meshtags.indices), np.sort(meshtags[i].indices), err_msg=f"Meshtags indices for dimension {i} are incomplete")

which fails with this output:

AssertionError: 
Not equal to tolerance rtol=1e-07, atol=0
Meshtags values for dimension 1 are incomplete
(shapes (107,), (320,) mismatch)
 ACTUAL: array([  0,   4,   5,   7,  15,  30,  40,  44,  47,  61,  64,  67,  77,
        81,  84,  87, 100, 104, 107, 110, 113, 116, 126, 129, 130, 132,
       133, 135, 136, 138, 139, 141, 142, 144, 145, 147, 150, 152, 154,...
 DESIRED: array([  0,   1,   2,   3,   4,   5,   6,   7,   8,   9,  10,  11,  12,
        13,  14,  15,  16,  17,  18,  19,  20,  21,  22,  23,  24,  25,
        26,  27,  28,  29,  30,  31,  32,  33,  34,  35,  36,  37,  38,...

Some of the meshtag information is clearly missing. Would you have any hints for this one ?

I am running version v0.10.0.post2 of dolfinx on the dolfinx/dolfinx:stable docker container.

Best,

2

Right, I haven’t thought about the case where you would store a meshtag while writing the mesh in input order.

I guess the easiest way to go about this is to convert the meshtag into a DG-0 function, and write that on the input order mesh (would only work for cell tags).

If this is a feature you would like, I would raise an issue in the new repo:

that includes more IO than checkpointing (it can also write data to formats visualizable in Paraview/vtk).