MeshValueCollection in dolfinx

dolfinx
1

Dear community,

I’m currently transitioning my scripts from dolfin to dolfinx and wondered where the MeshValueCollection has gone. Unfortunately, all the supplied demos only work with simple meshes and there are no examples where a mesh is read from an XDMF file including mesh values.

My code snippet

import numpy as np
from mpi4py import MPI
from petsc4py import PETSc
from dolfinx import *
from dolfinx.mesh import *
from dolfinx.cpp.mesh import *
from dolfinx.io import *
from dolfinx.la import *
from ufl import *

with XDMFFile(MPI.COMM_WORLD,'path/file.xdmf','r',encoding=dolfinx.cpp.io.XDMFFile.Encoding.ASCII) as infile:
    mesh = infile.read_mesh("Grid")

mvc_s = MeshValueCollection("size_t", mesh, 2)
with XDMFFile(MPI.COMM_WORLD,'path/file_boundary.xdmf','r',encoding=dolfinx.cpp.io.XDMFFile.Encoding.ASCII) as infile:
    infile.read(mvc_s, "boundaries_surf")

This produces
NameError: name ‘MeshValueCollection’ is not defined

Does anyone have an idea what happend to this function (and the MeshFunction as well)?

I’m using the latest dolfinx from today. The code snippet works fine for dolfin (with the modifications on the XDMF read…).

Maybe someone has an idea or a different solution how to read the boundary values?

Thanks!

Best,
Marc

2

The new syntax is using MeshTags, which is a sparse representation of a MesFunction, combining the benifits is MeshFunction and MeshValueCollection.
See these tests or this demo for a he new syntax.

3

Thank you! I made some progress, however I’m struggling to read line and point data… I have the following boundary XDMF file:

<Xdmf Version="3.0">
  <Domain>
    <Grid Name="Grid">
      <Topology NodesPerElement="3" NumberOfElements="24" TopologyType="triangle">
        <DataItem DataType="Int" Dimensions="24 3" Format="XML" Precision="4">
          1 0 8
          0 2 8
          3 1 8
          2 3 8
          0 1 10
          4 0 10
          1 5 10
          5 4 10
          1 3 12
          5 1 12
          3 7 12
          7 5 12
          5 9 4
          4 9 6
          7 9 5
          6 9 7
          2 11 3
          6 11 2
          3 11 7
          7 11 6
          0 13 2
          4 13 0
          2 13 6
          6 13 4
        </DataItem>
      </Topology>
      <Attribute AttributeType="Scalar" Center="Face" Name="boundaries_surf">
        <DataItem DataType="Int" Dimensions="24" Format="XML" Precision="4">
          1
          1
          1
          1
          2
          2
          2
          2
          3
          3
          3
          3
          4
          4
          4
          4
          5
          5
          5
          5
          6
          6
          6
          6
        </DataItem>
      </Attribute>
    </Grid>
    <Grid Name="Grid2">
      <Topology NodesPerElement="2" NumberOfElements="2" TopologyType="Polyline">
        <DataItem DataType="Int" Dimensions="2 3" Format="XML" Precision="4">
          1 3
          1 0
        </DataItem>
      </Topology>
      <Attribute AttributeType="Scalar" Center="Edge" Name="boundaries_edge">
        <DataItem DataType="Int" Dimensions="2" Format="XML" Precision="4">
          1
          2
        </DataItem>
      </Attribute>
    </Grid>
    <Grid Name="Grid3">
      <Topology NodesPerElement="1" NumberOfElements="1" TopologyType="Polyvertex">
        <DataItem DataType="Int" Dimensions="1 3" Format="XML" Precision="4">
          1
        </DataItem>
      </Topology>
      <Attribute AttributeType="Scalar" Center="Node" Name="boundaries_point">
        <DataItem DataType="Int" Dimensions="1" Format="XML" Precision="4">
          1
        </DataItem>
      </Attribute>
    </Grid>
  </Domain>
</Xdmf>

I can manage to read the face element data with

mesh.topology.create_connectivity(mesh.topology.dim-1, mesh.topology.dim)
with XDMFFile(comm,'.../input_boundary.xdmf','r',encoding=encoding) as infile:
    mt = infile.read_meshtags(mesh, "Grid")

When I want to read the edge/point data

mesh.topology.create_connectivity(mesh.topology.dim-2, mesh.topology.dim-1)
with XDMFFile(comm,'.../input_boundary.xdmf','r',encoding=encoding) as infile:
    mt = infile.read_meshtags(mesh, "Grid2")

I always get
RuntimeError: Number of entities and values must match

I’m not sure what the line
mesh.topology.create_connectivity(mesh.topology.dim-2, mesh.topology.dim-1)
exactly does and if setting it like this for the edge data is right (?).

I could read this data easily with dolfin via the Attribute names (boundaries_surf, boundaries_edge, boundaries_point), but now these do not seem to be accessible anymore?

Thanks!

Bw,
Marc

P.S.: Maybe, if it helps for reproducibility, this would be the corresponding domain XDMF file:

<Xdmf Version="3.0">
  <Domain>
    <Grid Name="Grid">
      <Geometry GeometryType="XYZ">
        <DataItem DataType="Float" Dimensions="14 3" Format="XML" Precision="8">
          0.0 0.0 1.0
          0.0 0.0 0.0
          0.0 1.0 1.0
          0.0 1.0 0.0
          1.0 0.0 1.0
          1.0 0.0 0.0
          1.0 1.0 1.0
          1.0 1.0 0.0
          0.0 0.5 0.5
          1.0 0.5 0.5
          0.5 0.0 0.5
          0.5 1.0 0.5
          0.5 0.5 0.0
          0.5 0.5 1.0
        </DataItem>
      </Geometry>
      <Topology NodesPerElement="4" NumberOfElements="24" TopologyType="tetrahedron">
        <DataItem DataType="Int" Dimensions="24 4" Format="XML" Precision="4">
          9 11 10 12
          8 13 11 10
          11 10 13 9
          8 10 11 12
          1 0 8 10
          0 2 8 13
          10 0 13 4
          3 11 8 2
          1 8 3 12
          11 13 2 6
          4 13 9 6
          6 11 9 7
          10 4 9 5
          11 7 3 12
          12 9 7 5
          12 1 10 5
          13 8 0 10
          8 2 11 13
          10 8 1 12
          11 3 8 12
          7 11 9 12
          13 9 10 4
          11 9 13 6
          12 10 9 5
        </DataItem>
      </Topology>
    </Grid>
  </Domain>
</Xdmf>
4 
<DataItem DataType="Int" Dimensions="2 3" Format="XML" Precision="4"> 1 3 1 0 </DataItem>

That is because this is wrong, the dimensions should be 2 2.
Dolfinx does a more throughout check of your data than old dolfin does.

You need to create the connectivity between the entity you are loading and the cell. For the grid2 case, it is

mesh.topology.create_connectivity(mesh.topology.dim-2, mesh.topology.dim)

You have the same dimension issue for the point data, where it should be 1 1 and not 1 3.

5

Great! Thanks, this worked! :slight_smile:

Best wishes,
Marc