Skip to content

Commit

Permalink
Include field data (such as time information) in parallel VTK files (#…
Browse files Browse the repository at this point in the history 
…140)

* Write field data to PVTK files as well

Field data such as the current time are now also written to PVTK files,
in addition to the VTK files associated to each "process". This seems to
be needed for ParaView to have access to the field data.

* Update tests

* Add notes to docs on including time information

* Try to fix tests on Julia 1.6
  • Loading branch information
@jipolanco
jipolanco authored Apr 4, 2024
1 parent 8813546 commit 7eb6761
Show file tree
Hide file tree
Showing 6 changed files with 121 additions and 47 deletions.
19 changes: 15 additions & 4 deletions docs/src/grids/datasets.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,8 +14,8 @@ vtk_grid("fields", x, y, z) do vtk
vtk["Pressure"] = rand(Nx - 1, Ny - 1, Nz - 1) # scalar field attached to cells
vtk["Velocity"] = rand(3, Nx, Ny, Nz) # vector field attached to points
vtk["VelocityGradients"] = rand(3, 3, Nx, Ny, Nz) # 3×3 tensor field attached to points
vtk["date"] = "31/10/2021" # metadata ("field data" in VTK)
vtk["time"] = 0.42 # metadata ("field data" in VTK)
vtk["Date"] = "31/10/2021" # metadata ("field data" in VTK)
vtk["TimeValue"] = 0.42 # metadata ("field data" in VTK)
end
```

Expand All @@ -25,6 +25,17 @@ input.
In particular, note that the `Pressure` field is attached to cells instead of points, since it has dimensions ``(N_x - 1) × (N_y - 1) × (N_z - 1)``, which is the number of cells in structured files (see [Cells in structured formats](@ref)).
For more control, see [Dataset types](@ref) below.

!!! note "Time values"

As a sidenote, the `TimeValue` field used above [is special](https://docs.vtk.org/en/latest/design_documents/IOXMLTimeInFieldData.html#field-data-as-time-meta-data-in-vtk-xml-file-formats),
as it is interpreted by VTK and ParaView as the time associated to the dataset.
This can be convenient when writing time series data (for example, results from a simulation at different timesteps).

Other alternatives for including time information are to use [ParaView collections](@ref) (`.pvd` format), or to generate
[JSON `.series` files](https://gitlab.kitware.com/paraview/paraview/-/blob/v5.5.0/Documentation/release/ParaView-5.5.0.md#json-based-new-meta-file-format-for-series-added)
(the latter can't be currently done by WriteVTK, but the format is very simple).


## Passing tuples of arrays

Note that, in the example above, the `Velocity` vector field is passed as a single ``3 × N_x × N_y × N_z`` array, which is the layout ultimately used in VTK formats.
Expand Down Expand Up @@ -124,8 +135,8 @@ vtk_grid("fields_explicit", x, y, z) do vtk
vtk["Pressure", VTKCellData()] = rand(Nx - 1, Ny - 1, Nz - 1)
vtk["Velocity", VTKPointData()] = rand(3, Nx, Ny, Nz)
vtk["VelocityGradients", VTKPointData()] = rand(3, 3, Nx, Ny, Nz)
vtk["date", VTKFieldData()] = "31/10/2021"
vtk["time", VTKFieldData()] = 0.42
vtk["Date", VTKFieldData()] = "31/10/2021"
vtk["TimeValue", VTKFieldData()] = 0.42
end
```

Expand Down
10 changes: 10 additions & 0 deletions src/WriteVTK.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -107,6 +107,16 @@ end
DatasetFile(dtype, xdoc::XMLDocument, fname::AbstractString, args...; kwargs...) =
DatasetFile(xdoc, add_extension(fname, dtype), xml_name(dtype), args...; kwargs...)

function data_format(vtk::DatasetFile)
if vtk.appended
:appended
elseif vtk.ascii
:ascii
else
:inline
end
end

function show(io::IO, vtk::DatasetFile)
open_str = isopen(vtk) ? "open" : "closed"
print(io, "VTK file '$(vtk.path)' ($(vtk.grid_type) file, $open_str)")
Expand Down
42 changes: 40 additions & 2 deletions src/gridtypes/pvtk_grid.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,11 @@ struct PVTKFile <: VTKFile
path::String
end

# This is just to make a PVTKFile work like a DatasetFile.
# Only used when writing VTKFieldData to a PVTK file.
# Note that data is always written as text (ASCII).
data_format(::PVTKFile) = :ascii

# Returns true if the arguments do *not* contain any cell vectors.
_pvtk_is_structured(x::AbstractVector{<:AbstractMeshCell}, args...) = Val(false)
_pvtk_is_structured(x, args...) = _pvtk_is_structured(args...)
Expand Down Expand Up @@ -130,11 +135,44 @@ end

# Add point and cell data as usual
function Base.setindex!(
pvtk::PVTKFile, data, name::AbstractString, args...,
pvtk::PVTKFile, data, name::AbstractString, loc::AbstractFieldData,
)
pvtk.vtk[name, loc] = data
end

# In the case of field data, also add it to the pvtk file.
# Otherwise field data (typically the time / "TimeValue") is not seen by ParaView.
function Base.setindex!(
pvtk::PVTKFile, data, name::AbstractString, loc::VTKFieldData,
)
pvtk.vtk[name, loc] = data
if pvtk.pvtkargs.ismain
add_field_data(pvtk, data, name, loc)
end
data
end

# Used in add_field_data.
# We need to find the PUnstructuredGrid / PStructuredGrid / ... node.
function find_base_xml_node_to_add_field(pvtk::PVTKFile, loc)
xroot = root(pvtk.xdoc)
pgrid_type = "P" * pvtk.vtk.grid_type
find_element(xroot, pgrid_type)
end

# If `loc` was not passed, try to guess which kind of data was passed.
function Base.setindex!(
pvtk::PVTKFile, data, name::AbstractString,
)
pvtk.vtk[name, args...] = data
loc = guess_data_location(data, pvtk.vtk) :: AbstractFieldData
pvtk[name, loc] = data
end

# Write XML attribute.
# Example:
#
# pvtk[VTKCellData()] = Dict("HigherOrderDegrees" => "HigherOrderDegreesDataset")
#
function Base.setindex!(
pvtk::PVTKFile, attributes, loc::AbstractFieldData,
)
Expand Down
47 changes: 29 additions & 18 deletions src/write_data.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -110,17 +110,23 @@ function write_array(io, data::Tuple)
n
end

function add_data_ascii(xml, x::Union{Number,String})
function add_data_ascii(xml, x::Union{Number,AbstractString})
add_text(xml, " ")
add_text(xml, string(x))
end

add_data_ascii(xml, x::AbstractArray) = map(v -> add_data_ascii(xml, v), x)
function add_data_ascii(xml, x::AbstractArray)
for v ∈ x
add_data_ascii(xml, v)
end
nothing
end

function add_data_ascii(xml, data::Tuple)
for i in eachindex(data...), x in data
add_data_ascii(xml, x[i])
end
nothing
end

function set_num_components(xDA, vtk, data, loc)
Expand Down Expand Up @@ -173,9 +179,10 @@ function data_to_xml(vtk, xParent, data, name,
set_attribute(xDA, "ComponentName$(i - 1)", n) # 0-based
end
end
if vtk.appended
fmt = data_format(vtk)
if fmt === :appended
data_to_xml_appended(vtk, xDA, data)
elseif vtk.ascii
elseif fmt === :ascii
data_to_xml_ascii(vtk, xDA, data)
else
data_to_xml_inline(vtk, xDA, data)
Expand Down Expand Up @@ -209,7 +216,7 @@ containing the size of the data array in bytes.
"""
function data_to_xml_appended(vtk::DatasetFile, xDA::XMLElement, data)
@assert vtk.appended
@assert data_format(vtk) === :appended

buf = vtk.buf # append buffer
compress = vtk.compression_level > 0
Expand Down Expand Up @@ -258,7 +265,7 @@ end
Add inline, base64-encoded data to VTK XML file.
"""
function data_to_xml_inline(vtk::DatasetFile, xDA::XMLElement, data)
@assert !vtk.appended && !vtk.ascii
@assert data_format(vtk) === :inline
compress = vtk.compression_level > 0

# DataArray node
Expand Down Expand Up @@ -305,8 +312,8 @@ end
Add inline data to VTK XML file in ASCII format.
"""
function data_to_xml_ascii(vtk::DatasetFile, xDA::XMLElement, data)
@assert !vtk.appended && vtk.ascii
function data_to_xml_ascii(vtk::VTKFile, xDA::XMLElement, data)
@assert data_format(vtk) === :ascii
set_attribute(xDA, "format", "ascii")
add_text(xDA, "\n")
add_data_ascii(xDA, data)
Expand All @@ -320,18 +327,10 @@ end
Add either point or cell data to VTK file.
"""
function add_field_data(vtk::DatasetFile, data, name::AbstractString,
function add_field_data(vtk::VTKFile, data, name::AbstractString,
loc::AbstractFieldData;
component_names::Union{AbstractVector, Nothing}=nothing)
# Find Piece node.
xroot = root(vtk.xdoc)
xGrid = find_element(xroot, vtk.grid_type)

xbase = if loc === VTKFieldData()
xGrid
else
find_element(xGrid, "Piece")
end
xbase = find_base_xml_node_to_add_field(vtk, loc)

# Find or create "nodetype" (PointData, CellData or FieldData) node.
nodetype = node_type(loc)
Expand All @@ -344,6 +343,18 @@ function add_field_data(vtk::DatasetFile, data, name::AbstractString,
xDA
end

function find_base_xml_node_to_add_field(vtk::DatasetFile, loc)
# Find Piece node.
xroot = root(vtk.xdoc)
xGrid = find_element(xroot, vtk.grid_type)
xbase = if loc === VTKFieldData()
xGrid
else
find_element(xGrid, "Piece")
end
xbase
end

vtk_point_data(args...; kwargs...) = add_field_data(args..., VTKPointData(); kwargs...)
vtk_cell_data(args...; kwargs...) = add_field_data(args..., VTKCellData(); kwargs...)
vtk_field_data(args...; kwargs...) = add_field_data(args..., VTKFieldData(); kwargs...)
Expand Down
46 changes: 23 additions & 23 deletions test/checksums.sha1
View file
Original file line number Diff line number Diff line change
Expand Up @@ -55,28 +55,28 @@ b20bf14b2531aaf323b891e090d0615d0e943ef6 collection_03.vtr
8005b27b071ed50765d0d859a9e4e84a52c3963c arraydata_2D.vti
d5215c84b3c1630bedd98994963560ff14200dc5 arraydata_3D.vti
aa49366e66d8930fa639a8361d153a95baca361e arrays.vti
2626bef24152fcfe28c810ca68adf473c8fad90c simulation.pvtu
4a0be2eaf5c7ffdb60970f0396a11ef2ab3affca simulation/simulation_1.vtu
3edeb95b6d2f955c8ee5539cb0a64be1f53f0bf1 simulation.pvtu
1348c734f1b79881f6b74461aaee059fa2f01c09 simulation/simulation_1.vtu
ee3ea76e543c6e3a9c1c28d80596b68c5b22d7ac higherorderdegrees.pvtu
9a524fb558d6e8e2a42c6bf704da88c2a835b296 higherorderdegrees/higherorderdegrees_1.vtu
e23c8f0db331e1db635aa5d98b9f1375132e7ffd pimage.pvti
e581d4824d67655584f573743f42788fb26e162c pimage/pimage_1.vti
150ce17cb8ddbb974344954a3193a3373ee1c3fc pimage/pimage_2.vti
e389b7948af0bd9cd9d178dcc6d0fb86493738f3 pimage/pimage_3.vti
669fdb06ab07c289cb5c01dd511da380fb06bdfd pimage/pimage_4.vti
1ea223999c053e38ab84726c0a4bfdd0468435ef pimage/pimage_5.vti
21dd9d15e51720f9d5822bf663c718974094d779 pimage/pimage_6.vti
d35ceff83c6c43b43eefe07007cbba8752cead15 prectilinear.pvtr
b7edbc7bb8907b472551d456f60cbc7f0c915b53 prectilinear/prectilinear_1.vtr
d26dd7cf0570b1fbad8714e16e784f09bb6affed prectilinear/prectilinear_2.vtr
e36a533248051994e7da804bff1b168d2d866bba prectilinear/prectilinear_3.vtr
4b73039f56d1a92fe958a98c86dde64b3fb86b63 prectilinear/prectilinear_4.vtr
87a07f8911949b002501da5bdd5daf497a65d09a prectilinear/prectilinear_5.vtr
4fd1d4de3867230a0636bef0260417ed17b7c17f prectilinear/prectilinear_6.vtr
fc3a26b061585c093520bbb720714d180ca22335 pstructured.pvts
584520780f94c58a9919cb216d116733f8466e0f pstructured/pstructured_1.vts
0ff33138288a96f8f6f4377f2e061748aa702d97 pstructured/pstructured_2.vts
e365dc765a17a288640e04f726cf8a7421559f39 pstructured/pstructured_3.vts
0daf5cccfaec20cdcfae94294513ef211954c9e1 pstructured/pstructured_4.vts
a2b4e0e63dda63eee19466b08d5d4ae8d34dcfd5 pstructured/pstructured_5.vts
0adbea4b9bebd2b95569398666028c19f3d8cb4a pstructured/pstructured_6.vts
03f54866860e984b2537a118def1c857f855d6ab pimage.pvti
3badee9575f3feb593109644d6124416731085da pimage/pimage_1.vti
a898a771d8e0962628118d73a216a3e4660a5295 pimage/pimage_2.vti
436cb1c2024e5795664a1355c8d19b7f336c2184 pimage/pimage_3.vti
4badae7ec804773110e662e105531e5e3d4fd288 pimage/pimage_4.vti
99843bc0f87c436d9ed385b40707ee1e0327f709 pimage/pimage_5.vti
4e37ea3120610664425712254b88d88f78d4e45d pimage/pimage_6.vti
9c0d3e7ee95c6e7bfc2f0554972d31fb9d832172 prectilinear.pvtr
e7ff55d7c8bb2c1739a2dc73c0c1fdeff2f7d536 prectilinear/prectilinear_1.vtr
6b0b2f8e5d381ec4f9644905677be8f29f432af2 prectilinear/prectilinear_2.vtr
68f065b66a0ff9ba80bc18776912f248ebc83390 prectilinear/prectilinear_3.vtr
b8b456b99f45f45e2e0fcc0de56963185b3f24ec prectilinear/prectilinear_4.vtr
a6fd8a0367efc31e7616dd7046f048c424aa8ae0 prectilinear/prectilinear_5.vtr
b1ef0b3a9e3c25781c31e0d2ceb552c78dab23f8 prectilinear/prectilinear_6.vtr
d9c0b2367bf7ba42534a134b3ef4dace753560a9 pstructured.pvts
92ff18a7c1c631800de3adc3b0bd4631e16be845 pstructured/pstructured_1.vts
9bb94bd6b2374156e60df995340e3b395f2ae6b9 pstructured/pstructured_2.vts
ad4b841fdf277c6138bbd763f4db2a359c46e400 pstructured/pstructured_3.vts
651905ba5698b503f58d40fc1051267fefa4b49e pstructured/pstructured_4.vts
bba1504744a03b3f9594fcfe63435b7311e030fe pstructured/pstructured_5.vts
2e8acea4c25a7ef64c25c05ff3f65b3890a52e72 pstructured/pstructured_6.vts
4 changes: 4 additions & 0 deletions test/pvtk_grid.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@ function pvtk_unstructured()
pvtk["Processor"] = [1,2]
pvtk["Temperature", VTKPointData()] = y
pvtk["Id", VTKCellData()] = [2, 1]
pvtk["TimeValue"] = 4.2
end
@test isfile(vtufile)
@test vtufile ∈ outfiles
Expand Down Expand Up @@ -73,6 +74,7 @@ function pvtk_imagedata()
) do vtk
vtk["point_data"] = point_data
vtk["process_id"] = processid
vtk["TimeValue"] = 1.2
end
end
end
Expand All @@ -95,6 +97,7 @@ function pvtk_rectilinear()
) do vtk
vtk["point_data"] = point_data
vtk["process_id"] = processid
vtk["TimeValue"] = 3.4
end
end
end
Expand Down Expand Up @@ -124,6 +127,7 @@ function pvtk_structured()
) do vtk
vtk["point_data"] = point_data
vtk["process_id"] = processid
vtk["TimeValue"] = 3.5
end
end
end
Expand Down

0 comments on commit 7eb6761

Please sign in to comment.