ClimaCoreVTK.jl

ClimaCoreVTK.jl provides functionality for writing ClimaCore fields to VTK files, using the WriteVTK.jl package.

Interface

ClimaCoreVTK.writevtkFunction
writevtk(
    basename::String,
    fields;
    basis=:cell,
    latlong=false,
    vtkargs...
)

Write fields to as an unstructured mesh VTK file named basename.vtu.

fields can be either:

  • a ClimaCore Field object,
  • a FieldVector object,
  • a NamedTuple of Fields.

The basis keyword option determines the type of cells used to write.:

  • :cell (default): output values at cell centers (interpolating where necessary).
  • :point: output values at cell vertices.
  • :lagrange: output values at Lagrange nodes (valid only for spectral element spaces), using Use VTK Lagrange cells to accurately represent high-order elements.

The latlong=true keyword option will output a spherical or spherical shell domain using the Mercator projection, with longitude along the x-axis, latitude along the y-axis, and altitude along the z-axis (if applicable). Note this currently only displays correctly if the number of elements across the cubed sphere face is even.

Any additional keyword arguments are passed to WriteVTK.vtk_grid.

source
ClimaCoreVTK.writepvdFunction
writepvd(
    basename::String,
    times,
    fields;
    vtkargs...
)

Write a sequence of fields at times as a Paraview collection (.pvd) file, along with VTK files.

fields can be either be an iterable collection of fields, or a NamedTuple of collections.

source

Internal functions

WriteVTK.vtk_gridFunction
vtk_grid(basename, gridspace::ClimaCore.Spaces.AbstractSpace;
    basis=:cell, vtkargs...)

Construct a VTK grid from a ClimaCore.Spaces.AbstractSpace. If basis=:lagrange, it will construct a mesh made of Lagrange cells (valid only for spectral element spaces), otherwise it will it subdivide the space into quads, with vertices at nodal points.

source
ClimaCoreVTK.vtk_cells_lagrangeFunction

vtkcellslagrange(space)

Construct a vector of MeshCell objects representing the elements of space as an unstuctured mesh of Lagrange polynomial cells, suitable for passing to vtk_grid.

source
ClimaCoreVTK.vtk_cells_linearFunction

vtkcellslinear(space)

Construct a vector of MeshCell objects representing the elements of space as an unstuctured mesh of linear cells, suitable for passing to vtk_grid.

source
ClimaCoreVTK.vtk_grid_spaceFunction
vtk_grid_space(space::ClimaCore.Spaces.AbstractSpace)

The space for the grid used by VTK, for any field on space.

This generally does two things:

  • Modifies the horizontal space to use a ClosedUniform quadrature rule, which will use equispaced nodal points in the reference element. This is required for using VTK Lagrange elements (see 1).
  • Modifies the vertical space to be on the faces.
source
ClimaCoreVTK.vtk_cell_spaceFunction
vtk_cell_space(gridspace::ClimaCore.Spaces.AbstractSpace)

Construct a space for outputting cell data, when using outputting a grid gridspace. be stored.

This generally does two things:

  • Modifies the horizontal space to use a Uniform quadrature rule, which will use equispaced nodal points in the reference element (excluding the boundary).
  • Modifies the vertical space to be on the centers.
source
ClimaCoreVTK.addfield!Function

addfield!(vtkfile, prefix::Union{String,Nothing}, f, dataspace)

Add a field or fields f, optionally prefixing the name with prefix to the VTK file vtkfile, interpolating to dataspace.

f can be any of the following:

  • a scalar or vector field (if no prefix is provided, then the field will be named "data")
  • a composite field, which will be named accordingly
  • a NamedTuple of fields
source