Usage Reference for pyvisfile.vtk
¶
Constants¶
Vector formats¶
- pyvisfile.vtk.VF_LIST_OF_COMPONENTS¶
[[x0, y0, z0], [x1, y1, z1]]
- pyvisfile.vtk.VF_LIST_OF_VECTORS¶
[[x0, x1], [y0, y1], [z0, z1]]
Element types¶
- pyvisfile.vtk.VTK_VERTEX¶
- pyvisfile.vtk.VTK_POLY_VERTEX¶
- pyvisfile.vtk.VTK_LINE¶
- pyvisfile.vtk.VTK_POLY_LINE¶
- pyvisfile.vtk.VTK_TRIANGLE¶
- pyvisfile.vtk.VTK_TRIANGLE_STRIP¶
- pyvisfile.vtk.VTK_POLYGON¶
- pyvisfile.vtk.VTK_PIXEL¶
- pyvisfile.vtk.VTK_QUAD¶
- pyvisfile.vtk.VTK_TETRA¶
- pyvisfile.vtk.VTK_VOXEL¶
- pyvisfile.vtk.VTK_HEXAHEDRON¶
- pyvisfile.vtk.VTK_WEDGE¶
- pyvisfile.vtk.VTK_PYRAMID¶
- pyvisfile.vtk.VTK_LAGRANGE_CURVE¶
- pyvisfile.vtk.VTK_LAGRANGE_TRIANGLE¶
- pyvisfile.vtk.VTK_LAGRANGE_QUADRILATERAL¶
- pyvisfile.vtk.VTK_LAGRANGE_TETRAHEDRON¶
- pyvisfile.vtk.VTK_LAGRANGE_HEXAHEDRON¶
- pyvisfile.vtk.VTK_LAGRANGE_WEDGE¶
XML elements¶
- class pyvisfile.vtk.XMLElementBase[source]¶
Base type for XML elements.
- children: List[str | XMLElement]¶
- add_child(child: str | XMLElement) None [source]¶
Append a new child to the current element.
- class pyvisfile.vtk.XMLElement(tag: str, **attributes: Any)[source]¶
Bases:
XMLElementBase
- copy(children: Sequence[str | XMLElement] | None = None) XMLElement [source]¶
Make a copy of the element with new children.
- class pyvisfile.vtk.XMLRoot(child: str | XMLElement | None = None)[source]¶
Bases:
XMLElementBase
Binary encoders¶
- class pyvisfile.vtk.EncodedBuffer[source]¶
An interface for binary buffers for XML data (inline and appended).
- abstract raw_buffer() ByteString [source]¶
The raw buffer object that was used to construct this encoded buffer.
- abstract add_to_xml_element(xml_element: XMLElement) int [source]¶
Add encoded buffer to the given xml_element.
- Returns:
total size of encoded buffer in bytes.
- class pyvisfile.vtk.BinaryEncodedBuffer(buffer: ByteString)[source]¶
Bases:
EncodedBuffer
An encoded buffer that uses raw uncompressed binary data.
- __init__(buffer: ByteString) None [source]¶
- class pyvisfile.vtk.Base64EncodedBuffer(buffer: memoryview)[source]¶
Bases:
EncodedBuffer
An encoded buffer that uses
base64
data.- __init__(buffer: memoryview) None [source]¶
- class pyvisfile.vtk.Base64ZLibEncodedBuffer(buffer: ByteString)[source]¶
Bases:
EncodedBuffer
An encoded buffer that uses
base64
andzlib
compression.- __init__(buffer: ByteString) None [source]¶
Building blocks¶
- class pyvisfile.vtk.Visitable[source]¶
A generic class for objects that can be mapped to XML elements.
- generator_method: ClassVar[str]¶
Name of the method called in
invoke_visitor()
.
- invoke_visitor(visitor: XMLGenerator) XMLElement [source]¶
Visit the current object with the given visitor and generate the corresponding XML element.
- class pyvisfile.vtk.DataArray(name: str, container: Any, vector_padding: int = 3, vector_format: int = 0, components: int | None = None)[source]¶
Bases:
Visitable
A representation of a generic VTK DataArray.
The storage format (inline or appended) is determined by the
XMLGenerator
at writing time.- __init__(name: str, container: Any, vector_padding: int = 3, vector_format: int = 0, components: int | None = None) None [source]¶
- Parameters:
name – name of the data array.
container – a
numpy.ndarray
or anotherDataArray
.vector_padding – pad any
ndarray
with additional zeros given by this variable.vector_format –
VF_LIST_OF_COMPONENTS
orVF_LIST_OF_VECTORS
.components – number of components in the container (not used).
- get_encoded_buffer(encoder: str, compressor: str | None = None) EncodedBuffer [source]¶
Re-encode the underlying buffer of the current
DataArray
.- Parameters:
encoder – new encoder name.
compressor – new compressor name.
- encode(compressor: str | None, xml_element: XMLElement) int [source]¶
Encode the underlying buffer with the given compressor and add it to the xml_element.
The re-encoding is performed using
get_encoded_buffer()
and the result is added to the element usingEncodedBuffer.add_to_xml_element()
. The encoding is always done inbase64
.
- class pyvisfile.vtk.UnstructuredGrid(points: tuple[int, DataArray], cells: ndarray[Any, dtype[Any]] | tuple[int, DataArray, DataArray], cell_types: ndarray[Any, dtype[Any]] | DataArray)[source]¶
Bases:
Visitable
- __init__(points: tuple[int, DataArray], cells: ndarray[Any, dtype[Any]] | tuple[int, DataArray, DataArray], cell_types: ndarray[Any, dtype[Any]] | DataArray) None [source]¶
- Parameters:
- class pyvisfile.vtk.StructuredGrid(mesh: ndarray[Any, dtype[Any]])[source]¶
Bases:
Visitable
XML generators¶
- class pyvisfile.vtk.XMLGenerator(compressor: str | None = None, vtk_file_version: str | None = None)[source]¶
- __init__(compressor: str | None = None, vtk_file_version: str | None = None) None [source]¶
- Parameters:
vtk_file_version –
a string
"x.y"
with the desired VTK XML file format version. Relevant versions are as follows:"0.1"
is the original version."1.0"
added support for 64-bit indices and offsets, as described here."2.0"
added support for ghost array data, as described here."2.1"
: added support for writing additional information attached to aDataArray
using information keys."2.2"
: changed the node numbering of the hexahedron, as described here.
- class pyvisfile.vtk.InlineXMLGenerator(compressor: str | None = None, vtk_file_version: str | None = None)[source]¶
Bases:
XMLGenerator
An XML generator that uses inline
DataArray
entries.
- class pyvisfile.vtk.AppendedDataXMLGenerator(compressor: str | None = None, vtk_file_version: str | None = None)[source]¶
Bases:
InlineXMLGenerator
An XML generator that uses appended data for
DataArray
entries.This creates a special element called
AppendedData
and each data array will index into it. Additional compression can be added to the appended data.
Convenience functions¶
- pyvisfile.vtk.write_structured_grid(file_name: str | Path, mesh: ndarray[Any, dtype[Any]], cell_data: Sequence[tuple[str, ndarray[Any, dtype[Any]]]] | None = None, point_data: Sequence[tuple[str, ndarray[Any, dtype[Any]]]] | None = None, overwrite: bool = False) None [source]¶
Write a structure grid to filename.
This constructs a
StructuredGrid
and adds the relevant point and cell data, as necessary. The data is all flattened to one dimensional arrays.- Parameters:
overwrite – if True, existing files are overwritten, otherwise an exception is raised.
VTK High-Order Lagrange Elements¶
The high-order elements are described in this blog post. The ordering of the element nodes is as follows:
the vertices (in an order that matches the linear elements, e.g.
VTK_TRIANGLE
).the interior edge (or face in 2D) nodes, i.e. without the endpoints
the interior face (3D only) nodes, i.e. without the edge nodes.
the remaining interior nodes.
For simplices, the interior nodes are defined recursively by using the same rules. However, for box elements the interior nodes are just listed in order, with the last coordinate moving slowest.
To a large extent, the VTK ordering matches the ordering used by gmsh
and
described here.
- pyvisfile.vtk.vtk_ordering.vtk_lagrange_simplex_node_tuples(dims: int, order: int, vtk_version: tuple[int, int] = (2, 1)) Sequence[tuple[int, ...]] [source]¶
- Parameters:
dims – dimension of the simplex, i.e. 1 corresponds to a curve, 2 to a triangle, etc.
order – order of the polynomial representation, which also defines the number of nodes on the simplex.
vtk_version – a
tuple
of two elements containing the version of the VTK XML file format in use. The ordering of some of the high-order elements changed between versions 2.1 and 2.2.
- Returns:
a
list
ofdims
-dimensional tuples of integers up toorder
in the ordering expected by VTK.
- pyvisfile.vtk.vtk_ordering.vtk_lagrange_simplex_node_tuples_to_permutation(node_tuples: Sequence[tuple[int, ...]]) Sequence[int] [source]¶
Construct a permutation from the simplex node ordering of VTK to that of
modepy
.- Returns:
a
list
of indices in[0, len(node_tuples)]
.
- pyvisfile.vtk.vtk_ordering.vtk_lagrange_quad_node_tuples(dims: int, order: int, vtk_version: tuple[int, int] = (2, 1)) Sequence[tuple[int, ...]] [source]¶
- Parameters:
dims – dimension of the box, i.e. 1 corresponds to a curve, 2 to a quadrilateral, and 3 to a hexahedron.
order – order of the polynomial representation, which also defines the number of nodes on the box.
vtk_version – a
tuple
of two elements containing the version of the VTK XML file format in use. The ordering of some of the high-order elements changed between versions 2.1 and 2.2.
- Returns:
a
list
ofdims
-dimensional tuples of integers up toorder
in the ordering expected by VTK.
Examples¶
Writing a structured mesh¶
# contributed by Luke Olson
import numpy as np
from pyvisfile.vtk import write_structured_grid
n = 50
x, y = np.meshgrid(np.linspace(-1, 1, n), np.linspace(-1, 1, n))
u = np.exp(-50 * (x**2 + y**2))
mesh = np.rollaxis(np.dstack((x, y)), 2)
write_structured_grid(
"test.vts",
mesh,
point_data=[("u", u[np.newaxis, :, :])])
(You can find this example as
examples/vtk-structured-2d-plain.py
in the PyVisfile
source distribution.)
Writing a collection of points¶
Note
Observe that this is written as a ‘unstructured grid’, even though there is not much grid here. However, by supplying connectivity data, it is possible to generalize from this to actual unstructured meshes.
import pathlib
import numpy as np
from pyvisfile.vtk import (
VF_LIST_OF_COMPONENTS,
VF_LIST_OF_VECTORS,
VTK_VERTEX,
AppendedDataXMLGenerator,
DataArray,
UnstructuredGrid,
)
rng = np.random.default_rng(seed=42)
n = 5000
points = rng.normal(size=(n, 3))
data = [
("pressure", rng.normal(size=n)),
("velocity", rng.normal(size=(3, n)))]
grid = UnstructuredGrid(
(n, DataArray("points", points, vector_format=VF_LIST_OF_VECTORS)),
cells=np.arange(n, dtype=np.uint32),
cell_types=np.array([VTK_VERTEX] * n, dtype=np.uint8))
for name, field in data:
grid.add_pointdata(
DataArray(name, field, vector_format=VF_LIST_OF_COMPONENTS))
file_name = pathlib.Path("points.vtu")
compressor = None
if file_name.exists():
raise FileExistsError(f"Output file '{file_name}' already exists")
with open(file_name, "w") as outf:
AppendedDataXMLGenerator(compressor)(grid).write(outf)
(You can find this example as
examples/vtk-unstructured-points.py
in the PyVisfile
source distribution.)