From 31016f947cfc5799b06cd5042b699f94ae995862 Mon Sep 17 00:00:00 2001
From: Dennis Ogiermann <termi-official@users.noreply.github.com>
Date: Tue, 26 Sep 2023 16:32:02 +0200
Subject: [PATCH] Introduce collections in spiral wave example.

---
 examples/spiral-wave.jl | 10 +++++++---
 src/Thunderbolt.jl      |  2 ++
 src/exports.jl          |  8 ++++++++
 src/mesh/meshes.jl      |  6 ++++++
 4 files changed, 23 insertions(+), 3 deletions(-)

diff --git a/examples/spiral-wave.jl b/examples/spiral-wave.jl
index 1e6c35e0..f1ad0242 100644
--- a/examples/spiral-wave.jl
+++ b/examples/spiral-wave.jl
@@ -516,8 +516,7 @@ function setup_initial_condition!(problem::SplitProblem, cache, initial_conditio
 end
 
 struct GalerkinDiscretization
-    # TODO interpolation collection instead of single interpolation
-    interpolations::Dict{Symbol, Interpolation}
+    interpolations::Dict{Symbol, Thunderbolt.InterpolationCollection}
 end
 
 """
@@ -525,11 +524,16 @@ end
 
 Transform a space-time model into a pure time-dependent problem.
 """
+semidiscretize
+
 function semidiscretize(split::ReactionDiffusionSplit{MonodomainModel{A,B,C,D,E}}, discretization::GalerkinDiscretization, grid::Thunderbolt.AbstractGrid) where {A,B,C,D,E}
     epmodel = split.model
 
+    ets = elementtypes(grid)
+    @assert length(ets) == 1
+
     # TODO get these from the interpolation collection in GalerkinDiscretization
-    ip = Lagrange{RefQuadrilateral, 1}()
+    ip = getinterpolation(discretization.interpolations[:φₘ], ets[1])
     dh = DofHandler(grid)
     push!(dh, :ϕₘ, ip)
     close!(dh);
diff --git a/src/Thunderbolt.jl b/src/Thunderbolt.jl
index 9680d0ff..d70ff29f 100644
--- a/src/Thunderbolt.jl
+++ b/src/Thunderbolt.jl
@@ -9,6 +9,8 @@ using JLD2
 import Ferrite: AbstractDofHandler, AbstractGrid, AbstractRefShape, AbstractCell
 import Ferrite: vertices, edges, faces, sortedge, sortface
 
+include("collections.jl")
+
 include("mesh/meshes.jl")
 include("mesh/coordinate_systems.jl")
 include("mesh/tools.jl")
diff --git a/src/exports.jl b/src/exports.jl
index a936ac9e..79e552b1 100644
--- a/src/exports.jl
+++ b/src/exports.jl
@@ -1,4 +1,11 @@
 export
+    # Collections
+    LagrangeCollection,
+    getinterpolation,
+    QuadratureRuleCollection,
+    getquadraturerule,
+    CellValueCollection,
+    FaceValueCollection,
     # Mesh generators
     generate_mesh,
     generate_ring_mesh,
@@ -53,6 +60,7 @@ export
     LoadDrivenSolver,
     # Utils
     calculate_volume_deformed_mesh,
+    elementtypes,
     value,
     # IO
     ParaViewWriter,
diff --git a/src/mesh/meshes.jl b/src/mesh/meshes.jl
index bb646a59..e862324b 100644
--- a/src/mesh/meshes.jl
+++ b/src/mesh/meshes.jl
@@ -35,6 +35,9 @@ num_faces(mgrid::SimpleMesh3D) = length(mgrid.mfaces)
 num_edges(mgrid::SimpleMesh3D) = length(mgrid.medges)
 num_vertices(mgrid::SimpleMesh3D) = length(mgrid.mvertices)
 
+elementtypes(::SimpleMesh3D{Tetrahedron}) = @SVector [Tetrahedron]
+elementtypes(::SimpleMesh3D{Hexahedron}) = @SVector [Hexahedron]
+
 function to_mesh(grid::Grid{3,C,T}) where {C, T}
     mfaces = OrderedDict{NTuple{3,Int}, Int}()
     medges = OrderedDict{NTuple{2,Int}, Int}()
@@ -124,6 +127,9 @@ num_nodes(mgrid::SimpleMesh2D) = length(mgrid.grid.nodes)
 num_faces(mgrid::SimpleMesh2D) = length(mgrid.mfaces)
 num_vertices(mgrid::SimpleMesh2D) = length(mgrid.mvertices)
 
+elementtypes(::SimpleMesh2D{Triangle}) = @SVector [Triangle]
+elementtypes(::SimpleMesh2D{Quadrilateral}) = @SVector [Quadrilateral]
+
 # Ferrite compat layer for the mesh
 @inline Ferrite.getncells(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getncells(mesh.grid)
 @inline Ferrite.getcells(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getcells(mesh.grid)