Add capability to handle with node variables

src/convert.jl

@@ -124,7 +124,7 @@ function trixi2vtk(filename::AbstractString...;
     if is_datafile
       verbose && println("| Reading data file...")
       @timeit "read data" (labels, data, n_elements, n_nodes,
-                           element_variables, time) = read_datafile(filename)
+                           element_variables, node_variables, time) = read_datafile(filename)
 
       assert_cells_elements(n_elements, mesh, filename, meshfile)
 
@@ -202,6 +202,26 @@ function trixi2vtk(filename::AbstractString...;
             verbose && println("| | Element variable: $label...")
             @timeit label vtk_celldata[label] = variable
           end
+
+          # Add node variables
+          for (label, variable) in node_variables
+            verbose && println("| | Node variable: $label...")
+            if reinterpolate
+              if index == 1
+                println("WARNING: The limiting coefficients are no continuous field but happens " *
+                "to be represented by a piecewise-constant approximation. Thus, reinterpolation " *
+                "does not give a meaningful representation.")
+              end
+              @timeit "interpolate data" interpolated_cell_data = interpolate_data(Val(format),
+                                                                    reshape(variable, size(variable)..., 1),
+                                                                    mesh, n_visnodes, verbose)
+            else
+              @timeit "interpolate data" interpolated_cell_data = reshape(variable,
+                                                                          n_visnodes^ndims_ * n_elements)
+            end
+            # Add the "interpolated" cell_data to celldata, not node_data
+            @timeit label vtk_nodedata[label] = interpolated_cell_data
+          end
         end
       end
     end

src/io.jl

@@ -76,6 +76,21 @@ function read_datafile(filename::String)
       index +=1
     end
 
-    return labels, data, n_elements, n_nodes, element_variables, time
+    # Extract node variable arrays
+    node_variables = Dict{String, Union{Array{Float64}, Array{Int}}}()
+    index = 1
+    while haskey(file, "node_variables_$index")
+      varname = read(attributes(file["node_variables_$index"])["name"])
+      nodedata = read(file["node_variables_$index"])
+      if ndims_ == 2
+        node_variables[varname] = Array{Float64}(undef, n_nodes, n_nodes, n_elements)
+        @views node_variables[varname][:, :, :] .= nodedata
+      else
+        error("Unsupported number of spatial dimensions: ", ndims_)
+      end
+      index +=1
+    end
+
+    return labels, data, n_elements, n_nodes, element_variables, node_variables, time
   end
 end

test/test_2d.jl

@@ -395,6 +395,43 @@ end
       end
     end
   end
+
+  @testset "Subcell limiting coefficients" begin
+    isdir(outdir) && rm(outdir, recursive=true)
+    run_trixi(joinpath(examples_dir(), "tree_2d_dgsem", "elixir_euler_sedov_blast_wave_sc_subcell.jl"), maxiters=10)
+
+    @timed_testset "do not reinterpolate" begin
+      # Create and test output without reinterpolation
+      @test_nowarn trixi2vtk(joinpath(outdir, "solution_000010.h5"), output_directory=outdir, reinterpolate=false)
+      outfilename = "solution_000010.vtu"
+      out_file = joinpath(outdir, outfilename)
+
+      # save output file to `artifacts` to facilitate debugging of failing tests
+      testname = "2d-treemesh-shockcapturing-subcell-no-reinterp"
+      cp(out_file, joinpath(artifacts_dir, testname * "-" * outfilename), force=true)
+
+      # remote file path is actually a URL so it always has the same path structure
+      remote_filename = "2d/treemesh/dgsem_sedov_subcell_no_interp_10.vtu"
+      ref_file = get_test_reference_file("dgsem_sedov_subcell_no_interp_10.vtu", remote_filename)
+      compare_point_data(out_file, ref_file)
+    end
+
+    @timed_testset "do reinterpolate" begin
+      # Create and test output without reinterpolation
+      @test_nowarn trixi2vtk(joinpath(outdir, "solution_000010.h5"), output_directory=outdir, reinterpolate=true)
+      outfilename = "solution_000010.vtu"
+      out_file = joinpath(outdir, outfilename)
+
+      # save output file to `artifacts` to facilitate debugging of failing tests
+      testname = "2d-treemesh-shockcapturing-subcell-reinterp"
+      cp(out_file, joinpath(artifacts_dir, testname * "-" * outfilename), force=true)
+
+      # remote file path is actually a URL so it always has the same path structure
+      remote_filename = "2d/treemesh/dgsem_sedov_subcell_interp_10.vtu"
+      ref_file = get_test_reference_file("dgsem_sedov_subcell_interp_10.vtu", remote_filename)
+      compare_cell_data(out_file, ref_file)
+    end
+  end
 end
 
 # Clean up afterwards: delete Trixi output directory and reference file directory