plot broken scalar FE Types in broken plots

Author: pjaap

CHANGELOG.md

@@ -1,5 +1,11 @@
 # CHANGES
 
+## v1.3.0
+
+### Added
+  - `plot` for discontinuous the FE Spaces `L2P0` and `L2P1` are now discontinuous to avoid averaging of values at the grid nodes.
+     This can be disabled by passing the kwarg `average_broken_plots = true` to the `plot` call.
+
 ## v1.2.0 June 4, 2025
 
 ### Changed

Project.toml

@@ -1,7 +1,7 @@
 name = "ExtendableFEM"
 uuid = "a722555e-65e0-4074-a036-ca7ce79a4aed"
+version = "1.3.0"
 authors = ["Christian Merdon <merdon@wias-berlin.de>", "Patrick Jaap <patrick.jaap@wias-berlin.de>"]
-version = "1.2"
 
 [deps]
 CommonSolve = "38540f10-b2f7-11e9-35d8-d573e4eb0ff2"
@@ -18,6 +18,7 @@ Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
 SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 SparseDiffTools = "47a9eef4-7e08-11e9-0b38-333d64bd3804"
+StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
 TimerOutputs = "a759f4b9-e2f1-59dc-863e-4aeb61b1ea8f"
 UnicodePlots = "b8865327-cd53-5732-bb35-84acbb429228"
@@ -63,10 +64,9 @@ Metis = "2679e427-3c69-5b7f-982b-ece356f1e94b"
 OrdinaryDiffEqRosenbrock = "43230ef6-c299-4910-a778-202eb28ce4ce"
 OrdinaryDiffEqSDIRK = "2d112036-d095-4a1e-ab9a-08536f3ecdbf"
 SimplexGridFactory = "57bfcd06-606e-45d6-baf4-4ba06da0efd5"
-StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 TetGen = "c5d3f3f7-f850-59f6-8a2e-ffc6dc1317ea"
 Triangulate = "f7e6ffb2-c36d-4f8f-a77e-16e897189344"
 
 [targets]
-test = ["Aqua", "ExampleJuggler", "ExplicitImports", "IncompleteLU", "Metis", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqSDIRK", "SimplexGridFactory", "StaticArrays", "Test", "TetGen", "Triangulate"]
+test = ["Aqua", "ExampleJuggler", "ExplicitImports", "IncompleteLU", "Metis", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqSDIRK", "SimplexGridFactory", "Test", "TetGen", "Triangulate"]

src/ExtendableFEM.jl

@@ -49,7 +49,7 @@ using ExtendableGrids: ExtendableGrids, AT_NODES, AbstractElementGeometry,
     ON_BEDGES, ON_BFACES, ON_CELLS, ON_EDGES, ON_FACES,
     ON_IFACES, SerialVariableTargetAdjacency,
     UniqueBFaceGeometries, UniqueCellGeometries,
-    UniqueFaceGeometries, append!, dim_element, eval_trafo!,
+    UniqueFaceGeometries, append!, dim_element, dim_space, eval_trafo!,
     facetype_of_cellface, interpolate!,
     max_num_targets_per_source, num_cells, num_faces,
     num_nodes, num_sources, num_targets, simplexgrid,
@@ -74,6 +74,7 @@ using Symbolics: Symbolics
 using SciMLBase: SciMLBase
 using TimerOutputs: TimerOutput, print_timer, @timeit
 using UnicodePlots: UnicodePlots
+using StaticArrays: SVector
 
 ## reexport stuff from ExtendableFEMBase and ExtendableGrids
 export FESpace, FEMatrix, FEVector

src/plots.jl

@@ -19,7 +19,20 @@ function plot!(p::GridVisualizer, ops, sol; kwargs...)
 Plots the operator evaluations ops of blocks in sol into the GridVisualizer.
 
 """
-function plot!(p::GridVisualizer, ops, sol; rasterpoints = 10, linewidth = 1, keep = [], ncols = size(p.subplots, 2), do_abs = true, do_vector_plots = true, title_add = "", kwargs...)
+function plot!(
+        p::GridVisualizer,
+        ops,
+        sol;
+        rasterpoints = 10,
+        linewidth = 1,
+        keep = [],
+        ncols = size(p.subplots, 2),
+        do_abs = true,
+        do_vector_plots = true,
+        title_add = "",
+        average_broken_plots = false,
+        kwargs...
+    )
     col, row, id = 0, 1, 0
     for op in ops
         col += 1
@@ -48,7 +61,11 @@ function plot!(p::GridVisualizer, ops, sol; rasterpoints = 10, linewidth = 1, ke
                 title = op[2] == Identity ? "$(sol[op[1]].name)" : "$(op[2])($(sol[op[1]].name))"
             end
             if resultdim == 1
-                GridVisualize.scalarplot!(p[row, col], sol[op[1]].FES.dofgrid, view(nodevalues(sol[op[1]], op[2]; abs = false), 1, :), title = title * title_add; kwargs...)
+                if !average_broken_plots && ExtendableFEMBase.broken(sol[op[1]].FES)
+                    broken_scalarplot!(p[row, col], sol[op[1]]; title = title * title_add, kwargs...)
+                else
+                    GridVisualize.scalarplot!(p[row, col], sol[op[1]].FES.dofgrid, view(nodevalues(sol[op[1]], op[2]; abs = false), 1, :), title = title * title_add; kwargs...)
+                end
             elseif do_abs == true
                 GridVisualize.scalarplot!(p[row, col], sol[op[1]].FES.dofgrid, view(nodevalues(sol[op[1]], op[2]; abs = true), 1, :), title = "|" * title * "|" * title_add; kwargs...)
             else
@@ -72,6 +89,31 @@ function plot!(p::GridVisualizer, ops, sol; rasterpoints = 10, linewidth = 1, ke
 end
 
 
+"""
+    broken_scalarplot!(vis, feVectorBlock; kwargs...)
+
+A "broken" scalarplot of a broken finite element vector.
+Instead of averaging the discontinuous values on the grid nodes, each grid cell is plotted
+independently. Thus, a discontinuous plot is generated.
+
+All kwargs of the calling method are transferred to the scalarplot in this method.
+"""
+function broken_scalarplot!(vis, feVectorBlock::FEVectorBlock; kwargs...)
+
+    dofgrid = feVectorBlock.FES.dofgrid
+    cell_nodes = dofgrid[CellNodes]
+    coords = dofgrid[Coordinates]
+
+    all_values = nodevalues(feVectorBlock; cellwise = true) # cellwise evaluation of the FE
+    all_coords = @views coords[:, cell_nodes[:]]
+    all_cells = reshape(1:length(all_values), size(all_values))
+
+    GridVisualize.scalarplot!(vis, all_coords, all_cells, view(all_values, :); kwargs...)
+
+    return nothing
+end
+
+
 """
 ````
 function plot!(p::GridVisualizer, ops, sol; Plotter = nothing, kwargs...)