IO (dolfinx::io)

namespace io

Support for file IO.

IO to files for checkpointing and visualisation.

Enums

enum class FidesMeshPolicy

Mesh reuse policy.

Values:

enumerator update

Re-write the mesh to file upon every write of a fem::Function.

enumerator reuse

Write the mesh to file only the first time a fem::Function is written to file

enum class VTXMeshPolicy

Mesh reuse policy.

Values:

enumerator update

Re-write the mesh to file upon every write of a fem::Function.

enumerator reuse

Write the mesh to file only the first time a fem::Function is written to file

Functions

template<typename U, typename T>
VTXWriter(MPI_Comm comm, U filename, T mesh) -> VTXWriter<typename std::remove_cvref<typename T::element_type>::type::geometry_type::value_type>

Type deduction.

template<typename U, typename T>
FidesWriter(MPI_Comm comm, U filename, T mesh) -> FidesWriter<typename std::remove_cvref<typename T::element_type>::type::geometry_type::value_type>

Type deduction.

template<typename T>
std::tuple<std::vector<T>, 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<T> &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(MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan<const std::int32_t, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents<std::size_t, 2>> dofmap_x, mesh::CellType cell_type)

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:
  • dofmap_x[in] Geometry dofmap.

  • cell_type[in] Cell type.

Returns:

Cell topology in VTK ordering and in term of the DOLFINx geometry ‘nodes’.

class ADIOS2Writer
#include <ADIOS2Writers.h>

Base class for ADIOS2-based writers.

Public Functions

void close()

Close the file.

template<std::floating_point T>
class FidesWriter : public ADIOS2Writer

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

Public Functions

inline FidesWriter(MPI_Comm comm, const std::filesystem::path &filename, std::shared_ptr<const mesh::Mesh<T>> mesh, std::string engine = "BPFile")

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:
inline FidesWriter(MPI_Comm comm, const std::filesystem::path &filename, const typename adios2_writer::U<T> &u, std::string engine, const FidesMeshPolicy mesh_policy = FidesMeshPolicy::update)

Create Fides writer for list of functions.

Note

All functions in u must share the same Mesh. The element family and degree, and degree-of-freedom map (up to the blocksize) must be the same for all functions.

Parameters:
  • comm[in] The MPI communicator

  • filename[in] Name of output file

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

  • engine[in] ADIOS2 engine type. See https://adios2.readthedocs.io/en/latest/engines/engines.html.

  • mesh_policy[in] Controls if the mesh is written to file at the first time step only or is re-written (updated) at each time step.

inline FidesWriter(MPI_Comm comm, const std::filesystem::path &filename, const typename adios2_writer::U<T> &u, const FidesMeshPolicy mesh_policy = FidesMeshPolicy::update)

Create Fides writer for list of functions using default ADIOS2 engine,.

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.

  • mesh_policy[in] Controls if the mesh is written to file at the first time step only or is re-written (updated) at each time step.

FidesWriter(FidesWriter &&file) = default

Move constructor.

~FidesWriter() = default

Destructor.

FidesWriter &operator=(FidesWriter&&) = default

Move assignment.

inline void write(double t)

Write data with a given time.

Parameters:

t[in] The time step

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 for 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.

template<std::floating_point U>
void write(const mesh::Mesh<U> &mesh, double time = 0.0)

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

Parameters:
  • mesh[in] Mesh to write to file.

  • time[in] Time parameter to associate with mesh.

template<dolfinx::scalar T, std::floating_point U = scalar_value_type_t<T>>
void write(const std::vector<std::reference_wrapper<const fem::Function<T, U>>> &u, double t)

Write finite elements function with an associated time step.

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

  • t[in] Time parameter to associate with u

Pre:

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

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. Interpolate fem::Function before output if required.

Pre:

All Functions in u must share the same mesh

template<std::floating_point T>
class VTXWriter : public ADIOS2Writer

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

inline VTXWriter(MPI_Comm comm, const std::filesystem::path &filename, std::shared_ptr<const mesh::Mesh<T>> mesh, std::string engine = "BPFile")

Create a VTX writer for a mesh.

This format supports arbitrary degree meshes.

Note

This format supports 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] MPI communicator to open the file on.

  • filename[in] Name of output file.

  • mesh[in] Mesh to write.

  • engine[in] ADIOS2 engine type.

inline VTXWriter(MPI_Comm comm, const std::filesystem::path &filename, const typename adios2_writer::U<T> &u, std::string engine, VTXMeshPolicy mesh_policy = VTXMeshPolicy::update)

Create a VTX writer for a list of fem::Functions.

This format supports arbitrary degree meshes.

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, and degree-of-freedom map (up to the blocksize) must be the same for all functions.

  • engine[in] ADIOS2 engine type.

  • mesh_policy[in] Controls if the mesh is written to file at the first time step only or is re-written (updated) at each time step.

inline VTXWriter(MPI_Comm comm, const std::filesystem::path &filename, const typename adios2_writer::U<T> &u, VTXMeshPolicy mesh_policy = VTXMeshPolicy::update)

Create a VTX writer for a list of fem::Functions using the default ADIOS2 engine.

This format supports arbitrary degree meshes.

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.

  • mesh_policy[in] Controls if the mesh is written to file at the first time step only or is re-written (updated) at each time step.

VTXWriter(VTXWriter &&file) = default

Move constructor.

~VTXWriter() = default

Destructor.

VTXWriter &operator=(VTXWriter&&) = default

Move assignment.

inline void write(double t)

Write data with a given time stamp.

Parameters:

t[in] Time stamp to associate with 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, std::string file_mode, Encoding encoding = default_encoding)

Constructor.

XDMFFile(XDMFFile&&) = default

Move 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.

template<std::floating_point U>
void write_mesh(const mesh::Mesh<U> &mesh, 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<double> &geometry, std::string name, 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<double> read_mesh(const fem::CoordinateElement<double> &element, mesh::GhostMode mode, std::string name, 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(std::string name, 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::variant<std::vector<float>, std::vector<double>>, std::array<std::size_t, 2>> read_geometry_data(std::string name, 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(std::string grid_name, 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

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

Write a fem::Function to file.

Note

User interpolation to a suitable Lagrange space may be required to satisfy the precondition on u. The VTX output (io::VTXWriter) format is recommended over XDMF for discontinuous and/or high-order spaces.

Parameters:
  • u[in] Function to write to file.

  • t[in] Time stamp to associate with u.

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

Pre:

The fem::Function u must be (i) a lowest-order (P0) discontinuous Lagrange element or (ii) a continuous Lagrange element where the element ‘nodes’ are the same as the nodes of its mesh::Mesh. Otherwise an exception is raised.

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

Write MeshTags

Parameters:
  • meshtags[in]

  • x[in] Mesh geometry

  • 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 mesh::Mesh<double> &mesh, std::string name, 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(std::string name, std::string value, 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(std::string name, 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 adios2_writer

Typedefs

template<std::floating_point T>
using U = std::vector<std::variant<std::shared_ptr<const fem::Function<float, T>>, std::shared_ptr<const fem::Function<double, T>>, std::shared_ptr<const fem::Function<std::complex<float>, T>>, std::shared_ptr<const fem::Function<std::complex<double>, T>>>>
namespace cells

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

Basix ordering is used for the geometry nodes, and is shown below for a range of cell types.

Triangle:               Triangle6:          Triangle10:
v
^
|
2                       2                    2
|`\                     |`\                  | \
|  `\                   |  `\                6   4
|    `\                 4    `3              |     \
|      `\               |      `\            5   9   3
|        `\             |        `\          |         \
0----------1 --> u      0-----5----1         0---7---8---1


Quadrilateral:         Quadrilateral9:         Quadrilateral16:
v
^
|
2-----------3          2-----7-----3           2--10--11---3
|           |          |           |           |           |
|           |          |           |           7  14  15   9
|           |          5     8     6           |           |
|           |          |           |           6  12  13   8
|           |          |           |           |           |
0-----------1 --> u    0-----4-----1           0---4---5---1


Tetrahedron:                    Tetrahedron10:               Tetrahedron20
            v
            /
          2                            2                         2
        ,/|`\                        ,/|`\                     ,/|`\
      ,/  |  `\                    ,/  |  `\                 13  |  `9
      ,/    '.   `\                ,8    '.   `6           ,/     4   `\
    ,/       |     `\            ,/       4     `\       12    19 |     `8
  ,/         |       `\        ,/         |       `\   ,/         |       `\
0-----------'.--------1 -> u 0--------9--'.--------1 0-----14----'.--15----1
  `\.         |      ,/        `\.         |      ,/   `\.  17     |  16 ,/
    `\.      |    ,/             `\.      |    ,5       10.   18 5    ,6
        `\.   '. ,/                  `7.   '. ,/            `\.   '.  7
          `\. |/                       `\. |/                 11. |/
              `3                            `3                    `3
                `\.
                    w


Hexahedron:          Hexahedron27:
        w
    6----------7               6----19----7
    /|   ^   v /|              /|         /|
  / |   |  / / |            17 |  25    18|
  /  |   | / /  |            / 14    24 /  15
4----------5   |           4----16----5   |
|   |   +--|---|--> u      |22 |  26  | 23|
|   2------+---3           |   2----13+---3
|  /       |  /           10  / 21   12  /
| /        | /             | 9    20  | 11
|/         |/              |/         |/
0----------1               0-----8----1


Prism:                      Prism15:
            w
            ^
            |
            3                       3
          ,/|`\                   ,/|`\
        ,/  |  `\               12  |  13
      ,/    |    `\           ,/    |    `\
    4------+------5         4------14-----5
    |      |      |         |      8      |
    |    ,/|`\    |         |      |      |
    |  ,/  |  `\  |         |      |      |
    |,/    0    `\|         10     0      11
  ./|    ,/ `\    |\        |    ,/ `\    |
  /  |  ,/     `\  | `\      |  ,6     `7  |
u    |,/         `\|   v     |,/         `\|
    1-------------2         1------9------2


Pyramid:                     Pyramid13:
                4                            4
              ,/|\                         ,/|\
            ,/ .'|\                      ,/ .'|\
  v      ,/   | | \                   ,/   | | \
    \.  ,/    .' | `.                ,/    .' | `.
      \.      |  '.  \             11      |  12  \
    ,/  \.  .'  w |   \          ,/       .'   |   \
  ,/      \. |  ^ |    \       ,/         7    |    9
2----------\'--|-3    `.     2-------10-.'----3     `.
  `\       |  \.|  `\    \      `\        |      `\    \
    `\     .'   +----`\ - \ -> u  `6     .'         8   \
      `\   |           `\  \        `\   |           `\  \
        `\.'               `\         `\.'               `\
          0-----------------1           0--------5--------1

Functions

std::vector<std::uint16_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::uint16_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::uint16_t> transpose(std::span<const std::uint16_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(std::span<const std::int64_t> cells, std::array<std::size_t, 2> shape, std::span<const std::uint16_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 hdf5

Functions

template<typename T>
hid_t hdf5_type()

C++ type to HDF5 data type.

hid_t open_file(MPI_Comm comm, const std::filesystem::path &filename, const std::string &mode, 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

void close_file(hid_t handle)

Close HDF5 file

Parameters:

handle[in] HDF5 file handle

void flush_file(hid_t handle)

Flush data to file to improve data integrity after interruption

Parameters:

handle[in] HDF5 file handle

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

Get filename

Parameters:

handle[in] HDF5 file handle return The filename

bool has_dataset(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

hid_t open_dataset(hid_t handle, const std::string &path)

Open dataset

Parameters:
  • handle[in] HDF5 file handle.

  • path[in] Data set path.

Returns:

Data set handle. Should be closed by caller using H5Dclose.

std::vector<std::int64_t> get_dataset_shape(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)

void set_mpi_atomicity(hid_t handle, 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.

bool get_mpi_atomicity(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

void add_group(hid_t handle, const std::string &dataset_path)

Add group to HDF5 file

Parameters:
  • handle[in] HDF5 file handle

  • dataset_path[in] Data set path to add

template<typename T>
void write_dataset(hid_t file_handle, const std::string &dataset_path, const T *data, std::array<std::int64_t, 2> range, const std::vector<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:
  • file_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>
std::vector<T> read_dataset(hid_t dset_id, std::array<std::int64_t, 2> range, bool allow_cast)

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

Template Parameters:

T – The data type to read into.

Parameters:
  • dset_id[in] HDF5 file handle.

  • range[in] The local range on this processor.

  • allow_cast[in] If true, allow casting from HDF5 type to type T.

Returns:

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

namespace impl

Functions

template<typename T>
std::tuple<std::vector<T>, std::array<std::size_t, 2>, std::vector<std::int64_t>, std::vector<std::uint8_t>> tabulate_lagrange_dof_coordinates(const fem::FunctionSpace<T> &V)

Tabulate the coordinate for every ‘node’ in a Lagrange function space.

Parameters:

V[in] The function space

Pre:

V must be (discontinuous) Lagrange and must not be a subspace

Returns:

Mesh coordinate data

  1. Node coordinates (shape={num_dofs, 3}) where the ith row corresponds to the coordinate of the ith dof in V (local to process)

  2. Node coordinates shape

  3. Unique global index for each node

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

namespace impl_adios2

Functions

template<class T>
adios2::Attribute<T> define_attribute(adios2::IO &io, std::string name, const T &value, std::string var_name = "", std::string separator = "/")

Safe definition of an attribute. First check if it has already been defined and return it. If not defined create new attribute.

template<class T>
adios2::Variable<T> define_variable(adios2::IO &io, std::string name, const adios2::Dims &shape = adios2::Dims(), const adios2::Dims &start = adios2::Dims(), const adios2::Dims &count = adios2::Dims())

Safe definition of a variable. First check if it has already been defined and return it. If not defined create new variable.

template<std::floating_point T>
std::shared_ptr<const mesh::Mesh<T>> extract_common_mesh(const typename adios2_writer::U<T> &u)

Extract common mesh from list of Functions.

Variables

constexpr std::array field_ext = {"_real", "_imag"}

String suffix for real and complex components of a vector-valued field

namespace impl_fides

Functions

void initialize_mesh_attributes(adios2::IO &io, mesh::CellType type)

Initialize mesh related attributes for the ADIOS2 file used in Fides.

Parameters:
  • io[in] ADIOS2 IO object

  • type[in] Cell type

template<std::floating_point T>
void initialize_function_attributes(adios2::IO &io, const typename adios2_writer::U<T> &u)

Initialize function related attributes for the ADIOS2 file used in Fides

Parameters:
  • io[in] The ADIOS2 IO

  • u[in] List of functions

template<typename T, std::floating_point U>
std::vector<T> pack_function_data(const fem::Function<T, U> &u)

Pack Function data at vertices. The mesh and the function must both be ‘P1’.

template<typename T, std::floating_point U>
void write_data(adios2::IO &io, adios2::Engine &engine, const fem::Function<T, U> &u)

Write a first order Lagrange function (real or complex) using ADIOS2 in Fides format. Data is padded to be three dimensional if vector and 9 dimensional if tensor.

Parameters:
  • io[in] The ADIOS2 io object

  • engine[in] The ADIOS2 engine object

  • u[in] The function to write

template<std::floating_point T>
void write_mesh(adios2::IO &io, adios2::Engine &engine, const mesh::Mesh<T> &mesh)

Write mesh geometry and connectivity (topology) for Fides.

Parameters:
  • io[in] The ADIOS2 IO

  • engine[in] The ADIOS2 engine

  • mesh[in] The mesh

namespace impl_vtx

Functions

std::stringstream create_vtk_schema(const std::vector<std::string> &point_data, const std::vector<std::string> &cell_data)

Create VTK xml scheme to be interpreted by the VTX reader https://adios2.readthedocs.io/en/latest/ecosystem/visualization.html#saving-the-vtk-xml-data-model

template<std::floating_point T>
std::vector<std::string> extract_function_names(const typename adios2_writer::U<T> &u)

Extract name of functions and split into real and imaginary component.

template<typename T, std::floating_point X>
void vtx_write_data(adios2::IO &io, adios2::Engine &engine, const fem::Function<T, X> &u)

Given a Function, write the coefficient to file using ADIOS2.

Note

Only supports (discontinuous) Lagrange functions.

Note

For a complex function, the coefficient is split into a real and imaginary function.

Note

Data is padded to be three dimensional if vector and 9 dimensional if tensor.

Parameters:
  • io[in] ADIOS2 io object.

  • engine[in] ADIOS2 engine object.

  • u[in] Function to write.

template<std::floating_point T>
void vtx_write_mesh(adios2::IO &io, adios2::Engine &engine, const mesh::Mesh<T> &mesh)

Write mesh to file using VTX format

Parameters:
  • io[in] The ADIOS2 io object

  • engine[in] The ADIOS2 engine object

  • mesh[in] The mesh

template<std::floating_point T>
std::pair<std::vector<std::int64_t>, std::vector<std::uint8_t>> vtx_write_mesh_from_space(adios2::IO &io, adios2::Engine &engine, const fem::FunctionSpace<T> &V)

Given a FunctionSpace, create a topology and geometry based on the function space dof coordinates. Writes the topology and geometry using ADIOS2 in VTX format.

Note

Only supports (discontinuous) Lagrange functions.

Parameters:
  • io[in] The ADIOS2 io object

  • engine[in] The ADIOS2 engine object

  • V[in] The function space

namespace xdmf_function

Low-level methods for reading/writing XDMF files.

Functions

template<dolfinx::scalar T, std::floating_point U>
void add_function(MPI_Comm comm, const fem::Function<T, U> &u, double t, pugi::xml_node &xml_node, const hid_t h5_id)

Write a fem::Function to XDMF.

namespace xdmf_mesh

Low-level methods for reading XDMF files.

Functions

template<std::floating_point U>
void add_mesh(MPI_Comm comm, pugi::xml_node &xml_node, hid_t h5_id, const mesh::Mesh<U> &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.

template<std::floating_point U>
void add_topology_data(MPI_Comm comm, pugi::xml_node &xml_node, hid_t h5_id, std::string path_prefix, const mesh::Topology &topology, const mesh::Geometry<U> &geometry, int cell_dim, 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.

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

Add Geometry xml node.

std::pair<std::variant<std::vector<float>, std::vector<double>>, std::array<std::size_t, 2>> read_geometry_data(MPI_Comm comm, 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, 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)

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

Add mesh tags to XDMF file.

namespace xdmf_utils

Functions

std::pair<std::string, int> get_cell_type(const pugi::xml_node &topology_node)

Get DOLFINx cell type string from XML topology node

Returns:

DOLFINx cell type and polynomial degree

std::array<std::string, 2> get_hdf5_paths(const pugi::xml_node &dataitem_node)

Return (0) HDF5 filename and (1) path in HDF5 file from a 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::string vtk_cell_type_str(mesh::CellType cell_type, int num_nodes)

Get the VTK string identifier.

template<typename T>
std::pair<std::vector<std::int32_t>, std::vector<T>> distribute_entity_data(const mesh::Topology &topology, std::span<const std::int64_t> nodes_g, std::int64_t num_nodes_g, const fem::ElementDofLayout &cmap_dof_layout, MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan<const std::int32_t, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents<std::size_t, 2>> xdofmap, int entity_dim, MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan<const std::int64_t, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents<std::size_t, 2>> entities, std::span<const 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 be 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 the (local) cell-vertex connectivity for this triangle.

Parameters:
  • topology[in] A mesh topology.

  • nodes_g[in] Global ‘input’ indices for the mesh, as returned by Geometry::input_global_indices.

  • num_nodes_g[in] Global number of geometry nodes, as returned by Geometry::index_map()->size_global().

  • cmap_dof_layout[in] Coordinate element dof layout, computed using Geometry::cmap().create_dof_layout().

  • xdofmap[in] Dofmap for the mesh geometry (Geometry::dofmap).

  • 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, hid_t h5_id, const std::string &h5_path, std::span<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.

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

Get data associated with a data set node.

Warning

Data will be silently cast to type T if requested type and storage type differ.

Template Parameters:

T – Data type to read into.