Using high-order Lagrange elements in VTK

[mattfalcone1997]

This question relates more strictly to the implementation of the vtkNek5000Reader in VTK rather than nekRS itself, but I want to make sure that my reasoning is sound before I propose any additions to VTK. Currently, the vtkNek5000Reader breaks each spectral element into many smaller vtkHexahedron elements; hence, the solution is not correctly represented in VTK or Paraview. This also applies to gradients calculated in Paraview. Recently, VTK introduced functionality to use high-order Lagrange elements (vtkLagrangeHexahedron). If I am not mistaken, nekRS uses Lagrange basis functions based on the GLL quadrature points. Therefore, it should be possible to have the option in VTK to create the mesh using these elements and get a better representation of the solution in VTK and Paraview. Is this reasoning correct, and would it be helpful to implement?

I have created a test Python script that converts the current output into one using vtkLagrangeHexahedron cells. While not thoroughly tested, it gives a basis for how it would work and appears to yield reasonable results.

[fischer1]

Dear Matthew, Thanks for the note. Our approach to this has been that the cell-based approach is generally sufficient for viz—one of the more recent features of Nek is that the user can (from inside Nek) interpolate the data onto a finer mesh in each element using either the std. GLL distribution or one that is uniform in the reference element. For example, an 8x8x8 GLL distribution could be mapped to a 10x10x10 uniform distribution in each element. (The map GLL—> uniform is stable; the map uniform—>GLL is, however, unstable, so these are not good files for restart.) Generally, we compute gradients inside Nek and then dump the corresponding vector fields—a significant reason for this is that, indeed, we get the spectral accuracy of the high-order representation, which is what you're describing as being possible, but also, we get the full 64-bit working precision. Right now, the field files are 32 bit, unless the user asks for 64 bit output. Again, 32 bit seems to be OK for viz and for most restart situations. I think a more pressing matter might be whether we can have zf/fz or other compression techniques that can be used both for reading and writing on the Nek side and for reading on the viz side. All that being said, an accurate in-viz gradient feature would allow computation of derivative-based quantities such a vorticity and Lambda 2 from base data, which indeed does reduce I/O traffic since one wouldn't need to write the derived fields. Best, Paul

…

________________________________ From: Matthew Falcone ***@***.***> Sent: Thursday, March 27, 2025 5:25 AM To: Nek5000/nekRS ***@***.***> Cc: Subscribed ***@***.***> Subject: [Nek5000/nekRS] Using high-order Lagrange elements in VTK (Discussion #632) This question relates more strictly to the implementation of the vtkNek5000Reader in VTK rather than nekRS itself, but I want to make sure that my reasoning is sound before I propose any additions to VTK. Currently, the vtkNek5000Reader breaks each spectral element into many smaller vtkHexahedron elements; hence, the solution is not correctly represented in VTK or Paraview. This also applies to gradients calculated in Paraview. Recently, VTK introduced functionality to use high-order Lagrange elements (vtkLagrangeHexahedron). If I am not mistaken, nekRS uses Lagrange basis functions based on the GLL quadrature points. Therefore, it should be possible to have the option in VTK to create the mesh using these elements and get a better representation of the solution in VTK and Paraview. Is this reasoning correct, and would it be helpful to implement? I have created a test Python script that converts the current output into one using vtkLagrangeHexahedron cells. While not thoroughly tested, it gives a basis for how it would work and appears to yield reasonable results. — Reply to this email directly, view it on GitHub<https://urldefense.com/v3/__https://github.com/Nek5000/nekRS/discussions/632__;!!DZ3fjg!-9RdeK7BBXL4zvjGlcIhS22ipWBszKbsvUgtEmdpiGem9pKQRj5UKGUl9Y2TWbQ-KlWXadASO0cQWQVZ_Bdt6W9b8hU$>, or unsubscribe<https://urldefense.com/v3/__https://github.com/notifications/unsubscribe-auth/ACV35EMNLDPXOZ3OR6CUF6T2WO7XVAVCNFSM6AAAAABZ4SRO42VHI2DSMVQWIX3LMV43ERDJONRXK43TNFXW4OZYGEZTKNJQHE__;!!DZ3fjg!-9RdeK7BBXL4zvjGlcIhS22ipWBszKbsvUgtEmdpiGem9pKQRj5UKGUl9Y2TWbQ-KlWXadASO0cQWQVZ_BdtEIhc3yw$>. You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[mattfalcone1997]

The main motivation for my suggestion is to improve fidelity of calculation based on VTK filters such gradients, Q criterion or even sampling filters where the shape functions can be evaluated directly. From an implementation perspective, this would not replace the existing cell based approach but be an option if desired by the user, which can be activated via a method like SetUseHighOrder in VTK (I assume in Paraview this could be turned into a checkbox).

[mattfalcone1997]

In terms of compression, I am less knowledgeable. If zf/fz is referring to this, they are apparently implemented in VTK-m, albeit only for structured grids. So it may be possible to decompress within VTK on a per spectral element basis, although I am not sure.

[yslan]

Hi Matthew,

Here are some related discussions in the past:

• https://groups.google.com/g/nek5000/c/bSxxr3Kb2Yw/m/OJsXN4G6AAAJ
• https://gitlab.kitware.com/paraview/paraview/-/issues/20196#note_828117
• Support high order Cell type for Nek5000 reader visit-dav/visit#5043

There are 4 stages for rendering the data.
a. Simulation, GLL
b. Nek's f00001 files,
c. Paraview or Visit readers
d. Pixels

And, three workflows...

1. Default workflow:
   b. f00001 file: GLL, 32 bit
   c. Paraview or Visit readers: linear cells
   d. Pixels: linear interpolation

2. Our typical workaround for visualization:
   b. f00001 file: finer uniform grid, 32 bit
   c. Paraview or Visit readers: linear cells
   d. Pixels: linear interpolation

   For b. you can do it in (next branch) via

   void nrs_t::writeCheckpoint(double t, int step, bool enforceOutXYZ, bool enforceFP64, int N_, bool uniform)
   

   or call it in Nek5000's userchk inside udf.

   ifreguo = .true.
   call outpost(vx,vy,vz,pr,t,'aaa')
   

3. What your are proposing
   b. f00001 file: GLL, 64 bit
   c. Paraview or Visit readers: VTK_HIGHER_ORDER_HEXAHEDRON
   d. Pixels: high order interpolation?

   I put a question mark because I doubt all of the actual values on each pixel are rendered via high-order interpolation. Especially, things like isosurface and transparent volume rendering from high-order cells are likely done using some low-order operations. Otherwise, I have more questions about its robustness and performance.

For accuracy, I'll be interested if you found a case that workflow 3 shows visible improvement.
For performance, from a to b is fast. As we increase the problem size, we are willing to preprocess more on our side (including various QOI) to ease the cost (both memory and work) inside c and d.

[mattfalcone1997]

So far, I have only checked how to do the conversion rather than perform detailed comparisons, so cannot yet comment on performance or potential visualisation improvements. My main consideration is having VTK filters better reflect the solution that nekRS produces, particularly in relation to sampling, gradients and integration. I also do not think this should replace the existing setup but be an option within VTK/Paraview which users can select if desired.