JuliaVTK / WriteVTK.jl

Julia package for writing VTK XML files
Other
151 stars 32 forks source link

Velocity vectors aren't viewable as such in paraview or VisIt #101

Closed maxrudolph closed 2 years ago

maxrudolph commented 2 years ago

I used writevtk to save a velocity field to a .vtr file and I'm able to open the result in paraview, but for some reason I cannot display the velocity as arrows. I suspect that this is because in the XML format, WriteVTK does not declare Vectors="..." and Scalars="..." Although the velocity is written (correctly) as multicomponent pointdata, it may be that because it's not flagged as a vector quantity in the PointData field, Paraview and VisIt do not recognize it as a vector field.

From the VTK manual: https://vtk.org/Wiki/VTK/Tutorials/ScalarsVsVectors

From page 14 of the VTK file format manual here: https://vtk.org/wp-content/uploads/2015/04/file-formats.pdf

<PointData Scalars=”Temperature” Vectors=”Velocity”>
<DataArray Name=”Velocity” .../>
<DataArray Name=”Temperature” .../>
<DataArray Name=”Pressure” .../>
</PointData>

test.vtr.zip

This is the Julia code that was used to write the vtr file:

function visualization(grid::CartesianGrid,rho::Matrix,eta::Matrix,vn::Array{Float64} ; filename="test.vtr")
    vtk_grid(filename, grid.x, grid.y) do vtk
        vtk["rho"] = transpose(rho)
        vtk["viscosity"] = transpose(eta)
        vtk["Velocity"] = vn
    end
end

I've attached the .vtr file as well.

jipolanco commented 2 years ago

From the VTK manual: https://vtk.org/Wiki/VTK/Tutorials/ScalarsVsVectors

Thanks, I wasn't aware that VTK made the distinction between physical vectors and "multi-component data".

I have the impression that the issue is actually because VTK expects vectors to be 3D, and your example is in two dimensions. I say this because I tried manually adding Vectors="Velocity" to your example file (in a text editor):

      <PointData Vectors="Velocity">
        <DataArray type="Float64" Name="rho" NumberOfComponents="1" format="appended" offset="337"/>
        <DataArray type="Float64" Name="viscosity" NumberOfComponents="1" format="appended" offset="889"/>
        <DataArray type="Float64" Name="Velocity" NumberOfComponents="2" format="appended" offset="1483"/>
      </PointData>

and ParaView complains with:

(   3.593s) [paraview        ]vtkDataSetAttributes.cx:1294  WARN| vtkPointData (0x1986c9a0): Can not set attribute Vectors. Incorrect number of components.

And if I also set NumberOfComponents="3" for the "Velocity" field, things fail for a different reason.

In 3D, one can visualise vector fields without needing the Vectors attribute. See for instance the minimal example below. In ParaView, one just needs to select "velocity" as the orientation array when creating glyphs:

using WriteVTK
xs = (1:10, 1:4, 1:4)
Ns = length.(xs)
v⃗s = randn(3, Ns...)
vtk_grid("vectors", xs...; compress = false) do vtk
    vtk["velocity"] = v⃗s
end

image

So, in 2D, I guess a possible workaround is to add a fake third dimension (say z = 1:1).

It may be worth it to add optional support for "Vectors", "Scalars" and other attributes to WriteVTK, but for now it's not clear to me how ParaView (or Visit) uses this information. OK, I see this is explained in the VTK wiki page you linked:

In VTK, as in physics, the number of components has nothing to do with whether a component is a scalar. The only thing that makes a scalar a "scalar" is that it is unchanged when a geometric transformation is applied to the data. [...]

fredrikekre commented 2 years ago

In Ferrite.jl we pad with zeroes for this reason: https://github.com/Ferrite-FEM/Ferrite.jl/blob/9bfc3ad57a9a31728a9853f5f63625da6e43cff2/src/Export/VTK.jl#L54 but that is not necessary for the grid points I believe, just the data points.

maxrudolph commented 2 years ago

Thank you for helping with this. I verified that if I add a third dimension to the 2D velocity, it is interpreted as velocity by paraview and I can construct stream tracers. I was not aware that even for a 2D dataset, VTK keeps information about the third dimension.