IO (dolfinx::io)

namespace io

Support for file IO.

IO to files for checkpointing and visualisation.

Functions

std::tuple<std::vector<double>, std::array<std::size_t, 2>, std::vector<std::int64_t>, std::vector<std::uint8_t>, std::vector<std::int64_t>, std::array<std::size_t, 2>> vtk_mesh_from_space(const fem::FunctionSpace &V)

Given a FunctionSpace, create a topology and geometry based on the dof coordinates.

Parameters

V[in] The function space

Pre

V must be a (discontinuous) Lagrange space

Returns

Mesh data

  1. node coordinates (shape={num_nodes, 3}), row-major storage

  2. node coordinates shape

  3. unique global ID for each node (a node that appears on more than one rank will have the same global ID)

  4. ghost index for each node (0=non-ghost, 1=ghost)

  5. cells (shape={num_cells, nodes_per_cell)}), row-major storage

  6. cell shape (shape={num_cells, nodes_per_cell)})

std::pair<std::vector<std::int64_t>, std::array<std::size_t, 2>> extract_vtk_connectivity(const mesh::Mesh &mesh)

Extract the cell topology (connectivity) in VTK ordering for all cells the mesh. The ‘topology’ includes higher-order ‘nodes’. The index of a ‘node’ corresponds to the index of DOLFINx geometry ‘nodes’.

Note

The indices in the return array correspond to the point indices in the mesh geometry array

Note

Even if the indices are local (int32), both Fides and VTX require int64 as local input

Parameters

mesh[in] The mesh

Returns

The cell topology in VTK ordering and in term of the DOLFINx geometry ‘nodes’

class ADIOS2Writer
#include <ADIOS2Writers.h>

Base class for ADIOS2-based writers.

Subclassed by FidesWriter, VTXWriter

Public Types

using Fdr = fem::Function<double>

Typedefs.

using Fdc = fem::Function<std::complex<double>>

Typedefs.

using U = std::vector<std::variant<std::shared_ptr<const Fdr>, std::shared_ptr<const Fdc>>>

Typedefs.

Public Functions

void close()

Close the file.

class FidesWriter : public ADIOS2Writer
#include <ADIOS2Writers.h>

Output of meshes and functions compatible with the Fides Paraview reader, see https://fides.readthedocs.io/en/latest/paraview/paraview.html.

Public Functions

FidesWriter(MPI_Comm comm, const std::filesystem::path &filename, std::shared_ptr<const mesh::Mesh> mesh)

Create Fides writer for a mesh.

Note

The mesh geometry can be updated between write steps but the topology should not be changed between write steps

Parameters
  • comm[in] The MPI communicator to open the file on

  • filename[in] Name of output file

  • mesh[in] The mesh. The mesh must a degree 1 mesh.

FidesWriter(MPI_Comm comm, const std::filesystem::path &filename, const ADIOS2Writer::U &u)

Create Fides writer for list of functions.

Note

All functions in u must share the same Mesh

Parameters
  • comm[in] The MPI communicator

  • filename[in] Name of output file

  • u[in] List of functions. The functions must (1) share the same mesh (degree 1) and (2) be degree 1 Lagrange.

FidesWriter(FidesWriter &&file) = default

Move constructor.

~FidesWriter() = default

Destructor.

FidesWriter &operator=(FidesWriter&&) = default

Move assignment.

void write(double t)

Write data with a given time.

Parameters

t[in] The time step

class VTXWriter : public ADIOS2Writer
#include <ADIOS2Writers.h>

Writer for meshes and functions using the ADIOS2 VTX format, see https://adios2.readthedocs.io/en/latest/ecosystem/visualization.html#using-vtk-and-paraview.

The output files can be visualized using ParaView.

Public Functions

VTXWriter(MPI_Comm comm, const std::filesystem::path &filename, std::shared_ptr<const mesh::Mesh> mesh)

Create a VTX writer for a mesh. This format supports arbitrary degree meshes.

Note

This format support arbitrary degree meshes

Note

The mesh geometry can be updated between write steps but the topology should not be changed between write steps

Parameters
  • comm[in] The MPI communicator to open the file on

  • filename[in] Name of output file

  • mesh[in] The mesh to write

VTXWriter(MPI_Comm comm, const std::filesystem::path &filename, const U &u)

Create a VTX writer for list of functions.

Note

This format supports arbitrary degree meshes

Parameters
  • comm[in] The MPI communicator to open the file on

  • filename[in] Name of output file

  • u[in] List of functions. The functions must (1) share the same mesh and (2) be (discontinuous) Lagrange functions. The element family and degree must be the same for all functions.

VTXWriter(VTXWriter &&file) = default

Move constructor.

~VTXWriter() = default

Destructor.

VTXWriter &operator=(VTXWriter&&) = default

Move assignment.

void write(double t)

Write data with a given time.

Parameters

t[in] The time step

class HDF5Interface
#include <HDF5Interface.h>

This class provides an interface to some HDF5 functionality.

Public Static Functions

static hid_t open_file(MPI_Comm comm, const std::filesystem::path &filename, const std::string &mode, const bool use_mpi_io)

Open HDF5 and return file descriptor.

Parameters
  • comm[in] MPI communicator

  • filename[in] Name of the HDF5 file to open

  • mode[in] Mode in which to open the file (w, r, a)

  • use_mpi_io[in] True if MPI-IO should be used

static void close_file(const hid_t handle)

Close HDF5 file.

Parameters

handle[in] HDF5 file handle

static void flush_file(const hid_t handle)

Flush data to file to improve data integrity after interruption.

Parameters

handle[in] HDF5 file handle

static std::filesystem::path get_filename(hid_t handle)

Get filename.

Parameters

handle[in] HDF5 file handle return The filename

template<typename T>
static void write_dataset(const hid_t handle, const std::string &dataset_path, const T *data, const std::array<std::int64_t, 2> &range, const std::vector<std::int64_t> &global_size, bool use_mpi_io, bool use_chunking)

Write data to existing HDF file as defined by range blocks on each process.

Parameters
  • handle[in] HDF5 file handle

  • dataset_path[in] Path for the dataset in the HDF5 file

  • data[in] Data to be written, flattened into 1D vector (row-major storage)

  • range[in] The local range on this processor

  • global_size[in] The global shape shape of the array

  • use_mpi_io[in] True if MPI-IO should be used

  • use_chunking[in] True if chunking should be used

template<typename T>
static std::vector<T> read_dataset(const hid_t handle, const std::string &dataset_path, const std::array<std::int64_t, 2> &range)

Read data from a HDF5 dataset “dataset_path” as defined by range blocks on each process.

Parameters
  • handle[in] HDF5 file handle

  • dataset_path[in] Path for the dataset in the HDF5 file

  • range[in] The local range on this processor

Returns

Flattened 1D array of values. If range = {-1, -1}, then all data is read on this process.

static bool has_dataset(const hid_t handle, const std::string &dataset_path)

Check for existence of dataset in HDF5 file.

Parameters
  • handle[in] HDF5 file handle

  • dataset_path[in] Data set path

Returns

True if dataset_path is in the file

static std::vector<std::int64_t> get_dataset_shape(const hid_t handle, const std::string &dataset_path)

Get dataset shape (size of each dimension)

Parameters
  • handle[in] HDF5 file handle

  • dataset_path[in] Dataset path

Returns

The shape of the dataset (row-major)

static void set_mpi_atomicity(const hid_t handle, const bool atomic)

Set MPI atomicity. See https://support.hdfgroup.org/HDF5/doc/RM/RM_H5F.html#File-SetMpiAtomicity and https://www.open-mpi.org/doc/v2.0/man3/MPI_File_set_atomicity.3.php Writes must be followed by an MPI_Barrier on the communicator before any subsequent reads are guaranteed to return the same data.

static bool get_mpi_atomicity(const hid_t handle)

Get MPI atomicity. See https://support.hdfgroup.org/HDF5/doc/RM/RM_H5F.html#File-GetMpiAtomicity and https://www.open-mpi.org/doc/v2.0/man3/MPI_File_get_atomicity.3.php.

class VTKFile
#include <VTKFile.h>

Output of meshes and functions in VTK/ParaView format. Isoparametric meshes of arbitrary degree are supported. For finite element functions, cell-based (DG0) and Lagrange (point-based) functions can be saved. For vertex-based functions the output must be isoparametic, i.e. the geometry and the finite element functions must be defined using the same basis.

Warning

This format is not suitable to checkpointing

Public Functions

VTKFile(MPI_Comm comm, const std::filesystem::path &filename, const std::string &file_mode)

Create VTK file.

~VTKFile()

Destructor.

void close()

Close file.

void flush()

Flushes XML files to disk.

void write(const mesh::Mesh &mesh, double time = 0.0)

Write a mesh to file. Supports arbitrary order Lagrange isoparametric cells.

Parameters
  • mesh[in] The Mesh to write to file

  • time[in] Time parameter to associate with mesh

template<typename T>
inline void write(const std::vector<std::reference_wrapper<const fem::Function<T>>> &u, double t)

Write finite elements function with an associated timestep.

Parameters
  • u[in] List of functions to write to file

  • t[in] Time parameter to associate with u

Pre

All Functions in u must share the same mesh

Pre

All Functions in u with point-wise data must use the same element type (up to the block size) and the element must be (discontinuous) Lagrange

Pre

Functions in u cannot be sub-Functions. Interpolate sub-Functions before output

class XDMFFile
#include <XDMFFile.h>

Read and write mesh::Mesh, fem::Function and other objects in XDMF.

This class supports the output of meshes and functions in XDMF (http://www.xdmf.org) format. It creates an XML file that describes the data and points to a HDF5 file that stores the actual problem data. Output of data in parallel is supported.

XDMF is not suitable for higher order geometries, as their currently only supports 1st and 2nd order geometries.

Public Types

enum class Encoding

File encoding type.

Values:

enumerator HDF5
enumerator ASCII

Public Functions

XDMFFile(MPI_Comm comm, const std::filesystem::path &filename, const std::string file_mode, const Encoding encoding = default_encoding)

Constructor.

~XDMFFile()

Destructor.

void close()

Close the file.

This closes open underlying HDF5 file. In ASCII mode the XML file is closed each time it is written to or read from, so close() has no effect.

void write_mesh(const mesh::Mesh &mesh, const std::string xpath = "/Xdmf/Domain")

Save Mesh.

Parameters
  • mesh[in]

  • xpath[in] XPath where Mesh Grid will be written

void write_geometry(const mesh::Geometry &geometry, const std::string name, const std::string xpath = "/Xdmf/Domain")

Save Geometry.

Parameters
  • geometry[in]

  • name[in]

  • xpath[in] XPath of a node where Geometry will be inserted

mesh::Mesh read_mesh(const fem::CoordinateElement &element, const mesh::GhostMode &mode, const std::string name, const std::string xpath = "/Xdmf/Domain") const

Read in Mesh.

Parameters
  • element[in] Element that describes the geometry of a cell

  • mode[in] The type of ghosting/halo to use for the mesh when distributed in parallel

  • name[in]

  • xpath[in] XPath where Mesh Grid is located

Returns

A Mesh distributed on the same communicator as the XDMFFile

std::pair<std::vector<std::int64_t>, std::array<std::size_t, 2>> read_topology_data(const std::string name, const std::string xpath = "/Xdmf/Domain") const

Read Topology data for Mesh.

Parameters
  • name[in] Name of the mesh (Grid)

  • xpath[in] XPath where Mesh Grid data is located

Returns

(Cell type, degree), and cells topology (global node indexing)

std::pair<std::vector<double>, std::array<std::size_t, 2>> read_geometry_data(const std::string name, const std::string xpath = "/Xdmf/Domain") const

Read Geometry data for Mesh.

Parameters
  • name[in] Name of the mesh (Grid)

  • xpath[in] XPath where Mesh Grid data is located

Returns

points on each process

std::pair<mesh::CellType, int> read_cell_type(const std::string grid_name, const std::string xpath = "/Xdmf/Domain")

Read information about cell type.

Parameters
  • grid_name[in] Name of Grid for which cell type is needed

  • xpath[in] XPath where Grid is stored

void write_function(const fem::Function<double> &u, double t, const std::string &mesh_xpath = "/Xdmf/Domain/Grid[@GridType='Uniform'][1]")

Write Function.

Parameters
  • u[in] The Function to write to file

  • t[in] The time stamp to associate with the Function

  • mesh_xpath[in] XPath for a Grid under which Function will be inserted

void write_function(const fem::Function<std::complex<double>> &u, double t, const std::string &mesh_xpath = "/Xdmf/Domain/Grid[@GridType='Uniform'][1]")

Write Function.

Parameters
  • u[in] The Function to write to file

  • t[in] The time stamp to associate with the Function

  • mesh_xpath[in] XPath for a Grid under which Function will be inserted

void write_meshtags(const mesh::MeshTags<std::int32_t> &meshtags, const std::string &geometry_xpath, const std::string &xpath = "/Xdmf/Domain")

Write MeshTags.

Parameters
  • meshtags[in]

  • geometry_xpath[in] XPath where Geometry is already stored in file

  • xpath[in] XPath where MeshTags Grid will be inserted

mesh::MeshTags<std::int32_t> read_meshtags(const std::shared_ptr<const mesh::Mesh> &mesh, const std::string name, const std::string xpath = "/Xdmf/Domain")

Read MeshTags.

Parameters
  • mesh[in] The Mesh that the data is defined on

  • name[in]

  • xpath[in] XPath where MeshTags Grid is stored in file

void write_information(const std::string name, const std::string value, const std::string xpath = "/Xdmf/Domain/")

Write Information.

Parameters
  • name[in]

  • value[in] String to store into Information tag

  • xpath[in] XPath where Information will be inserted

std::string read_information(const std::string name, const std::string xpath = "/Xdmf/Domain/")

Read Information.

Parameters
  • name[in]

  • xpath[in] XPath where Information is stored in file

MPI_Comm comm() const

Get the MPI communicator.

Returns

The MPI communicator for the file object

Public Static Attributes

static const Encoding default_encoding = Encoding::HDF5

Default encoding type.

namespace cells

Functions for the re-ordering of input mesh topology to the DOLFINx ordering, and transpose orderings for file output.

Functions

std::vector<std::uint8_t> perm_vtk(mesh::CellType type, int num_nodes)

Permutation array to map from VTK to DOLFINx node ordering.

Parameters
  • type[in] The cell shape

  • num_nodes[in] The number of cell ‘nodes’

Returns

Permutation array for permuting from VTK ordering to DOLFINx ordering, i.e. a_dolfin[i] = a_vtk[p[i]] @details Ifp = [0, 2, 1, 3]anda = [10, 3, 4, 7], thena_p =[a[p[0]], a[p[1]], a[p[2]], a[p[3]]] = [10, 4, 3, 7]`

std::vector<std::uint8_t> perm_gmsh(mesh::CellType type, int num_nodes)

Permutation array to map from Gmsh to DOLFINx node ordering.

Parameters
  • type[in] The cell shape

  • num_nodes[in]

Returns

Permutation array for permuting from Gmsh ordering to DOLFINx ordering, i.e. a_dolfin[i] = a_gmsh[p[i]] @details Ifp = [0, 2, 1, 3]anda = [10, 3, 4, 7], thena_p =[a[p[0]], a[p[1]], a[p[2]], a[p[3]]] = [10, 4, 3, 7]`

std::vector<std::uint8_t> transpose(const std::vector<std::uint8_t> &map)

Compute the transpose of a re-ordering map.

Parameters

map[in] A re-ordering map

Returns

Transpose of the map. E.g., is map = {1, 2, 3, 0}, the transpose will be {3 , 0, 1, 2 }.

std::vector<std::int64_t> apply_permutation(const std::span<const std::int64_t> &cells, std::array<std::size_t, 2> shape, const std::span<const std::uint8_t> &p)

Permute cell topology by applying a permutation array for each cell.

Parameters
  • cells[in] Array of cell topologies, with each row representing a cell (row-major storage)

  • shape[in] The shape of the cells array

  • p[in] The permutation array that maps a_p[i] = a[p[i]], where a_p is the permuted array

Returns

Permuted cell topology, where for a cell v_new[i] = v_old[map[i]]. The storage is row-major and the shape is the same as cells.

std::int8_t get_vtk_cell_type(mesh::CellType cell, int dim)

Get VTK cell identifier.

Parameters
  • cell[in] The cell type

  • dim[in] The topological dimension of the cell

Returns

The VTK cell identifier

namespace xdmf_function

Low-level methods for reading/writing XDMF files.

Functions

void add_function(MPI_Comm comm, const fem::Function<double> &u, const double t, pugi::xml_node &xml_node, const hid_t h5_id)

TODO.

void add_function(MPI_Comm comm, const fem::Function<std::complex<double>> &u, const double t, pugi::xml_node &xml_node, const hid_t h5_id)

TODO.

namespace xdmf_mesh

Low-level methods for reading XDMF files.

Functions

void add_mesh(MPI_Comm comm, pugi::xml_node &xml_node, const hid_t h5_id, const mesh::Mesh &mesh, const std::string path_prefix)

Add Mesh to xml node.

Creates new Grid with Topology and Geometry xml nodes for mesh. In HDF file data is stored under path prefix.

void add_topology_data(MPI_Comm comm, pugi::xml_node &xml_node, const hid_t h5_id, const std::string path_prefix, const mesh::Topology &topology, const mesh::Geometry &geometry, int cell_dim, const std::span<const std::int32_t> &entities)

Add Topology xml node.

Parameters
  • comm[in]

  • xml_node[in]

  • h5_id[in]

  • path_prefix[in]

  • topology[in]

  • geometry[in]

  • cell_dim[in] Dimension of mesh entities to save

  • entities[in] Local-to-process indices of mesh entities whose topology will be saved. This is used to save subsets of Mesh.

void add_geometry_data(MPI_Comm comm, pugi::xml_node &xml_node, const hid_t h5_id, const std::string path_prefix, const mesh::Geometry &geometry)

Add Geometry xml node.

std::pair<std::vector<double>, std::array<std::size_t, 2>> read_geometry_data(MPI_Comm comm, const hid_t h5_id, const pugi::xml_node &node)

Read geometry (coordinate) data.

Returns

The coordinates of each ‘node’. The returned data is (0) an array holding the coordinates (row-major storage) and (1) the shape of the coordinate array. The shape is (num_nodes, geometric dimension).

std::pair<std::vector<std::int64_t>, std::array<std::size_t, 2>> read_topology_data(MPI_Comm comm, const hid_t h5_id, const pugi::xml_node &node)

Read topology (cell connectivity) data.

Returns

Mesh topology in DOLFINx ordering, where data row i lists the ‘nodes’ of cell i. The returned data is (0) an array holding the topology data (row-major storage) and (1) the shape of the topology array. The shape is (num_cells, num_nodes_per_cell)

namespace xdmf_meshtags

Functions

template<typename T>
void add_meshtags(MPI_Comm comm, const mesh::MeshTags<T> &meshtags, pugi::xml_node &xml_node, const hid_t h5_id, const std::string name)

Add mesh tags to XDMF file.

namespace xdmf_read

Low-level methods for reading XDMF files.

Functions

template<typename T>
std::vector<T> get_dataset(MPI_Comm comm, const pugi::xml_node &dataset_node, const hid_t h5_id, std::array<std::int64_t, 2> range = {{0, 0}})

Return data associated with a data set node.

namespace xdmf_utils

Functions

std::pair<std::string, int> get_cell_type(const pugi::xml_node &topology_node)
std::array<std::string, 2> get_hdf5_paths(const pugi::xml_node &dataitem_node)
std::filesystem::path get_hdf5_filename(const std::filesystem::path &xdmf_filename)
std::vector<std::int64_t> get_dataset_shape(const pugi::xml_node &dataset_node)

Get dimensions from an XML DataSet node.

std::int64_t get_num_cells(const pugi::xml_node &topology_node)

Get number of cells from an XML Topology node.

std::vector<double> get_point_data_values(const fem::Function<double> &u)

Get point data values for linear or quadratic mesh into flattened 2D array.

std::vector<std::complex<double>> get_point_data_values(const fem::Function<std::complex<double>> &u)
std::vector<double> get_cell_data_values(const fem::Function<double> &u)

Get cell data values as a flattened 2D array.

std::vector<std::complex<double>> get_cell_data_values(const fem::Function<std::complex<double>> &u)
std::string vtk_cell_type_str(mesh::CellType cell_type, int num_nodes)

Get the VTK string identifier.

std::pair<std::vector<std::int32_t>, std::vector<std::int32_t>> distribute_entity_data(const mesh::Mesh &mesh, int entity_dim, const std::span<const std::int64_t> &entities, const std::span<const std::int32_t> &data)

Get owned entities and associated data from input entities defined by global ‘node’ indices. The input entities and data can be supplied on any rank and this function will manage the communication.

Note

This function involves parallel distribution and must be called collectively. Global input indices for entities which are not owned by current rank could passed to this function. E.g., rank0 provides an entity with global input indices [gi0, gi1, gi2], but this identifies a triangle that is owned by rank1. It will be distributed and rank1 will receive (local) cell-vertex connectivity for this triangle.

Parameters
  • mesh[in] A mesh

  • entity_dim[in] Topological dimension of entities to extract

  • entities[in] Mesh entities defined using global input indices (‘nodes’), typically from an input mesh file, e.g. [gi0, gi1, gi2] for a triangle. Let [v0, v1, v2] be the vertex indices of some triangle (using local indexing). Each vertex has a ‘node’ (geometry dof) index, and each node has a persistent input global index, so the triangle [gi0, gi1, gi2] could be identified with [v0, v1, v2]. The data is flattened and the shape is (num_entities, nodes_per_entity).

  • data[in] Data associated with each entity in entities.

Returns

(entity-vertex connectivity of owned entities, associated data (values) with each entity)

template<typename T>
void add_data_item(pugi::xml_node &xml_node, const hid_t h5_id, const std::string &h5_path, const T &x, std::int64_t offset, const std::vector<std::int64_t> &shape, const std::string &number_type, bool use_mpi_io)

TODO: Document.