visualization.md

[Visualization](@id visualization)

There are two possible approaches to visualize results from Trixi: either directly from the REPL using Plots.jl or with ParaView/VisIt by postprocessing Trixi's output files with Trixi2Vtk.

!!! note There is also a package Trixi2Img that allows to create images from Trixi's HDF5 output files. However, it is deprecated in favor of plotting directly from the REPL with Plots.jl.

[Plots.jl [experimental]](@id Plots.jl)

By far the easiest and most convenient plotting approach is to use the powerful Plots.jl package to directly visualize Trixi's results from the REPL. In the following, you will find more information on a number of topics for how to use Plots.jl together with Trixi:

1. [Getting started](@ref getting-started-plots-jl)
2. Customizing plot results via a PlotData2D object
3. Plotting a 3D solution as a 2D plot
4. Creating a 1D plot
5. Plotting a 2D or 3D solutions as a 1D plot
6. Visualizing results during a simulation

!!! note Plotting via Plots.jl is still considered an experimental feature and might change in any future releases.

[Getting started](@id getting-started-plots-jl)

After running a simulation with Trixi in the REPL, load the Plots package with

julia> using Plots

To visualize the solution, execute

julia> plot(sol)

Here we assume that sol holds the return value of the solve(...) method (with type DiffEqBase.ODESolution), which is the default variable name when you use one of the example elixirs. This will generate a grid layout with one subplot for each solution variable, convenient for getting an overview of the current solution:

You can save the resulting file as a PNG image file by calling savefig(...) with an output file name that ends in .png, e.g.,

julia> savefig("solution-overview.png")

In Trixi, two plot types are available: 2D heatmap plots and 1D line plots. If you use plot(sol), Trixi will automatically choose the plot type that fits the dimensions of the sol input: 2D/3D data will be visualized as a heatmap, 1D data as a line plot. For more fine-grained control over what to plot, you can create such an object yourself, which can either be a PlotData2D or a PlotData1D object. For further details on both of these see below:

Customizing plot results via a PlotData2D object

For more fine-grained control over what to plot, first create a PlotData2D object by executing

julia> pd = PlotData2D(sol)

This takes the results generated by Trixi and stores them in a data format that can be understood by the Plots package, and pd holds all data relevant for plotting sol. You can pass variable names as strings to pd using a dictionary-like syntax, e.g.,

julia> plot(pd["rho"])

This will create a single 2D heatmap plot of the variable rho:

The default plot type and style can be overridden by passing any additional arguments that are understood by the Plots package. For example, to change the color scheme and add names to the axes, modify the previous command to

julia> plot(pd["rho"], seriescolor = :heat, xguide="x", yguide="y")

to yield

For more details on the various format options for plot, please consult the Plots documentation.

In addition, you can plot the mesh lines on top of the solution variables by calling the getmesh(...) function on the PlotData2D object

julia> plot!(getmesh(pd)) # here we use `plot!` with an `!` to add to the previous plot

which modifies the previous plot to

By default, PlotData2D will convert the conserved variables to primitive variables, but this can be changed by passing an appropriate conversion function in the solution_variables keyword argument, similar to the behavior of the SaveSolutionCallback:

julia> pd = PlotData2D(sol; solution_variables=cons2cons) # Plot conservative variables

There are several other keyword arguments that influence how the solution data is processed for visualization with the Plots package. A detailed explanation can be found in the docstring of the PlotData2D constructor.

Another way to change the appearance of a plot is to convert the solution to a uniformly refined mesh before plotting. This can be helpful, e.g., when trying different settings for a simulation with adaptive mesh refinement, where one would like to ignore the mesh changes when comparing solutions. This is achieved with adapt_to_mesh_level, which uses the mesh adaptation routines to adapt the solution to a uniform grid. For example, the AMR solution from above could be preprocessed with

julia> pd = PlotData2D(adapt_to_mesh_level(sol, 4)...)

When plotted together with the mesh, this will yield the following visualization:

Plotting a 3D solution as a 2D plot

It is possible to plot 2D slices from 3D simulation data with the same commands as above:

julia> plot(sol) # `sol` is from a 3D simulation

By default, plotting sol or creating a PlotData2D object from a 3D simulation will create a 2D slice of the solution in the xy-plane. You can customize this behavior by explicitly creating a PlotData2D object and passing appropriate keyword arguments:

• slice specifies the plane which is being sliced and can be :xy, :xz, or :yz (default: :xy)
• point specifies a three-dimensional point. The sliced plane is then created such that it lies on the point (default: (0.0, 0.0, 0.0)). All other attributes for PlotData2D objects apply here as well.

For example, to plot the velocity field orthogonal to the yz-plane at different x-axis locations, you can execute

julia> trixi_include(joinpath(examples_dir(), "3d", "elixir_euler_taylor_green_vortex.jl"), tspan=(0, 1))
[...]

julia> plots = []
Any[]

julia> for x in range(0, stop=pi/2, length=6)
         pd = PlotData2D(sol, slice=:yz, point=(x, 0.0, 0.0))
         push!(plots, plot(pd["v1"], clims=(-1,1), title="x = "*string(round(x, digits=2))))
       end

julia> plot(plots..., layout=(2, 3), size=(750,350))

which results in a 2x3 grid of slices of the yz-plane:

Creating a 1D plot

When plotting a 1D solution with

julia> plot(sol) # `sol` is from a 1D simulation

Trixi automatically creates a PlotData1D object and visualizes it as a line plot:

To customize your 1D plot, you can create a PlotData1D object manually as follows:

julia> pd = PlotData1D(sol)

In a very similar fashion to PlotData2D, you can customize your plot:

• plot(pd) creates the same plot as in plot(sol).
• plot(pd["rho", "p"]) only plots specific variables. In this case rho and p.
• plot!(getmesh(pd)) adds mesh lines after creating a plot.
• Any attributes from Plots can be used, e.g., plot(pd, yguide=:temperature).
• pd = PlotData1D(adapt_to_mesh_level(sol, 4)...) adapts the mesh before plotting (in this example to a mesh with refinement level 4).

You can also customize the PlotData1D object itself by passing attributes to the PlotData1D constructor:

• solution_variables specifies the variables to be plotted.
• nvisnodes sets the amount of nodes per element which the solution then is interpolated on.

Plotting a 2D or 3D solutions as a 1D plot

It is possible to extract a straight, axis-parallel line from a 2D or 3D solution and visualize it as a 1D plot. This is done by creating a PlotData1D object with a 2D/3D solution sol as input:

julia> pd = PlotData1D(sol)

The plot is then created with:

julia> plot(pd)

By default the x-axis is extracted, which can be changed with following attributes:

• slice specifies the axis which is being extracted and can be :x, :y or :z (:z is only for 3D input and default is :x)
• point specifies a two or three dimensional point. The sliced axis is then created in such a way, that it lies on the point. (default: (0.0, 0.0) or (0.0, 0.0, 0.0))

All other attributes for PlotData1D objects apply here as well.

In the following, is an example for a 2D simulation of the linear scalar advection equation. First, we have the regular 2D heatmap plot:

From this, we can extract a line plot parallel to the y-axis going through the point (1.0, 0.0) with the following commands:

julia> pd = PlotData1D(sol, slice=:y, point=(1.0, 0.0))
julia> plot(pd)

This convenient method of slicing is limited to axis-parallel slices, but for 2D/3D solutions it is also possible to create a plot along any curve you want. To do so, you first need to create a list of 2D/3D points that define your curve. Then you can create a PlotData1D with the keyword argument curve set to your list.

Let's give an example of this with the basic advection equation from above by creating a plot along the circle marked in green:

We can write a function like this, that outputs a list of points on a circle:

function circle(radius, center, n_points)
    coordinates = zeros(2, n_points)
    for i in 1:n_points
        coordinates[:,i] = radius*[cospi(2*i/n_points), sinpi(2*i/n_points)] .+ center
    end
    return coordinates
end

Then create and plot a PlotData1D object along a circle with radius one, center at (1,1), and 100 points:

pd = PlotData1D(sol, curve=circle(1.0, (1.0, 1.0), 100))
plot(pd)

This gives you the following plot:

Creating a plot like this has its downsides. For one, it is unclear what to put on the abscissa of the plot. By default, the arc length of the given curve is used. Also, with this way of plotting you lose the ability to use a mesh plot from getmesh.

Visualizing results during a simulation

To visualize solutions while a simulation is still running (also known as in-situ visualization), you can use the VisualizationCallback. It is created as a regular callback and accepts upon creation a number of keyword arguments that allow, e.g., to control the visualization interval, to specify the variables to plot, or to customize the plotting style.

During the simulation, the visualization callback creates and displays visualizations of the current solution in regular intervals. This can be useful to, e.g., monitor the validity of a long-running simulation or for illustrative purposes. An example for how to create a VisualizationCallback can be found in examples/2d/elixir_advection_amr_visualization.jl:

[...]

# Enable in-situ visualization with a new plot generated every 20 time steps
# and additional plotting options passed as keyword arguments
visualization = VisualizationCallback(interval=20; clims=(0,1))

[...]

The resulting output of the referenced elixir can be seen in the embedded video below:

  <!--
  Video creation details
  * Set up terminal size and position appropriately
  * Record video as MP4 with SimpleScreenRecorder (https://en.wikipedia.org/wiki/SimpleScreenRecorder)
  * Upload to YouTube
  * Obtain responsive code by inserting link on https://embedresponsively.com
  -->
  <style>.embed-container { position: relative; padding-bottom: 56.25%; height: 0; overflow: hidden; max-width: 100%; } .embed-container iframe, .embed-container object, .embed-container embed { position: absolute; top: 0; left: 0; width: 100%; height: 100%; }</style><div class='embed-container'><iframe src='https://www.youtube-nocookie.com/embed/UZtrqeDY1Fs' frameborder='0' allowfullscreen></iframe></div>

Trixi2Vtk

Trixi2Vtk converts Trixi's .h5 output files to VTK files, which can be read by ParaView, VisIt, and other visualization tools. It automatically interpolates solution data from the original quadrature node locations to equidistant visualization nodes at a higher resolution, to make up for the loss of accuracy from going from a high-order polynomial representation to a piecewise constant representation in VTK.

In the Julia REPL, first load the package Trixi2Vtk

julia> using Trixi2Vtk

To process an HDF5 file generated by Trixi.jl, execute

julia> trixi2vtk(joinpath("out", "solution_000000.h5"), output_directory="out")

This will create two unstructured VTK files in the out subdirectory that can be opened with ParaView or VisIt: solution_000000.vtu contains the discontinuous Galerkin solution data while solution_000000_celldata.vtu holds any cell-based values such as the current AMR indicator or the cell refinement level.

This allows you to generate VTK files for solution, restart and mesh files. By default, Trixi2Vtk generates .vtu (unstructured VTK) files for both cell/element data (e.g., cell ids, element ids) and node data (e.g., solution variables). This format visualizes each cell with the same number of nodes, independent of its size. Alternatively, you can provide format=:vti as a keyword argument to trixi2vtk, which causes Trixi2Vtk to generate .vti (image data VTK) files for the solution files, while still using .vtu files for cell-/element-based data. In .vti files, a uniform resolution is used throughout the entire domain, resulting in different number of visualization nodes for each element. This can be advantageous to create publication-quality images, but increases the file size.

If you want to convert multiple solution/restart files at once, you can just supply multiple input files as the positional arguments to trixi2vtk, e.g.,

julia> trixi2vtk("out/solution_000000.h5", "out/solution_000040.h5")

You may also use file globbing to select a range of files based on filename patterns, e.g.,

julia> trixi2vtk("out/solution_*.h5")

to convert all solution files in the out/ directory or

julia> trixi2vtk("out/restart_00[0-9]000.h5")

to convert every one-thousandth restart file (out/restart_000000.h5, out/restart_001000.h5 etc.).

When multiple solution/restart files are provided, Trixi2Vtk will also generate a .pvd file, which allows ParaView to read all .vtu/.vti files at once and which uses the time attribute in solution/restart files to inform ParaView about the solution time. A comprehensive list of all possible arguments for trixi2vtk can be found in the Trixi2Vtk.jl API.

Further information regarding the development of Trixi2Vtk can be found in the [development section](@id trixi2vtk-dev).

Trixi2Img

!!! note "Trixi2Img is deprecated" Since it is possible to visualize results from Trixi directly from the REPL using the Plots package, Trixi2Img has been deprecated. There are still some features missing when using Plots, such as postprocessing HDF5 files. Once these have been added to the Plots-based solution, Trixi2Img will be retired.

Trixi2Img can be used to directly convert Trixi's output files to image files, without having to use a third-pary visualization tool such as ParaView or VisIt. The downside of this approach is that it generally takes longer to visualize the data (especially for large files) and that it does not allow as much customization of the generated output files. Currently, PNG and PDF are supported as output formats.

In the Julia REPL, first load the package Trixi2Img

julia> using Trixi2Img

To process an HDF5 file generated by Trixi.jl, execute

julia> trixi2img(joinpath("out", "solution_000040.h5"), output_directory="out", grid_lines=true)

This will create a file solution_000040_scalar.png in the out/ subdirectory that can be opened with any image viewer:

To visualize 3D data generated by Trixi.jl, a 2D slice must be extracted. A slice can only lie in a plane orthogonal to one of the coordinate axes. The slice plane is defined by the axis to which it is orthogonal and an axis intercept.

For example, to create a 2D slice in the xy-plane at the axis intercept z = 0.5 from 3D data, execute

julia> trixi2img(joinpath("out", "solution_000000.h5"), output_directory="out", grid_lines=true, slice_axis=:z, slice_axis_intercept=0.5)

Similar to Trixi2Vtk, you can also provide multiple files to trixi2img or use file globbing, e.g.,

julia> trixi2img("out/solution_*.h5")

to convert all solution files. The default is to generate a PNG file for each variable found in the respective file. Use format=:pdf as a keyword argument to create PDF files. A comprehensive list of all possible arguments for trixi2img can be found in the Trixi2Img.jl API.