FEniCS / dolfinx

Next generation FEniCS problem solving environment
https://fenicsproject.org
GNU Lesser General Public License v3.0
776 stars 181 forks source link

2D XDMFFile vector output broken #3094

Closed jorgensd closed 8 months ago

jorgensd commented 8 months ago

Summarize the issue

No longer padding xdmf-file vectors causes them to not be renderable in Paraview (say with the glyphs filter).

How to reproduce the bug

Run script below. Get:

cat mesh_0.7.3.xdmf 
<?xml version="1.0"?>
<!DOCTYPE Xdmf SYSTEM "Xdmf.dtd" []>
<Xdmf Version="3.0" xmlns:xi="https://www.w3.org/2001/XInclude">
  <Domain>
    <Grid Name="mesh" GridType="Uniform">
      <Topology TopologyType="Triangle" NumberOfElements="40" NodesPerElement="3">
        <DataItem Dimensions="40 3" NumberType="Int" Format="HDF">mesh_0.7.3.h5:/Mesh/mesh/topology</DataItem>
      </Topology>
      <Geometry GeometryType="XY">
        <DataItem Dimensions="33 2" Format="HDF">mesh_0.7.3.h5:/Mesh/mesh/geometry</DataItem>
      </Geometry>
    </Grid>
    <Grid Name="f" GridType="Collection" CollectionType="Temporal">
      <Grid Name="f" GridType="Uniform">
        <xi:include xpointer="xpointer(/Xdmf/Domain/Grid[@GridType='Uniform'][1]/*[self::Topology or self::Geometry])" />
        <Time Value="0" />
        <Attribute Name="f" AttributeType="Vector" Center="Node">
          <DataItem Dimensions="33 3" Format="HDF">mesh_0.7.3.h5:/Function/f/0</DataItem>
        </Attribute>
      </Grid>
    </Grid>
  </Domain>
</Xdmf>

and

root@5f5d5ecfb0e9:~/shared# cat mesh_0.8.0.0.xdmf 
<?xml version="1.0"?>
<!DOCTYPE Xdmf SYSTEM "Xdmf.dtd" []>
<Xdmf Version="3.0" xmlns:xi="https://www.w3.org/2001/XInclude">
  <Domain>
    <Grid Name="mesh" GridType="Uniform">
      <Topology TopologyType="Triangle" NumberOfElements="40" NodesPerElement="3">
        <DataItem Dimensions="40 3" NumberType="Int" Format="HDF">mesh_0.8.0.0.h5:/Mesh/mesh/topology</DataItem>
      </Topology>
      <Geometry GeometryType="XY">
        <DataItem Dimensions="33 2" Format="HDF">mesh_0.8.0.0.h5:/Mesh/mesh/geometry</DataItem>
      </Geometry>
    </Grid>
    <Grid Name="f" GridType="Collection" CollectionType="Temporal">
      <Grid Name="f" GridType="Uniform">
        <xi:include xpointer="xpointer(/Xdmf/Domain/Grid[@GridType='Uniform'][1]/*[self::Topology or self::Geometry])" />
        <Time Value="0" />
        <Attribute Name="f" AttributeType="Vector" Center="Node">
          <DataItem Dimensions="33 2" Format="HDF">mesh_0.8.0.0.h5:/Function/f/0</DataItem>
        </Attribute>
      </Grid>
    </Grid>
  </Domain>
</Xdmf>

Minimal Example (Python)

from mpi4py import MPI

import dolfinx
import numpy as np

mesh = dolfinx.mesh.create_unit_square(
    MPI.COMM_WORLD, 10, 2, ghost_mode=dolfinx.mesh.GhostMode.shared_facet)

with dolfinx.io.XDMFFile(mesh.comm, f"mesh_{dolfinx.__version__}.xdmf", "w") as xdmf:
    xdmf.write_mesh(mesh)
    V = dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (2,)))
    u = dolfinx.fem.Function(V)
    u.interpolate(lambda x: (x[1], np.zeros(x.shape[1], dtype=x.dtype)))
    xdmf.write_function(u)

Output (Python)

No response

Version

main branch

DOLFINx git commit

1ccffb0aa545634ea6e160236a7b043dc2a7fdd2

Installation

Using docker nightly image

Additional information

No response

jorgensd commented 8 months ago

Seems to have been introduced in: https://github.com/FEniCS/dolfinx/pull/2848 I think we need to have special handling for rank 1 tensors. @IgorBaratta what do you think?