diff --git a/NEWS.md b/NEWS.md
index 40a7a09b694..00ee67c2aaf 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -14,7 +14,15 @@ for human readability.
 #### Deprecated
 
 #### Removed
-- The shallow-water equation types `ShallowWaterEquations1D`, `ShallowWaterEquations2D`, and `ShallowWaterEquationsQuasi1D` have been removed from Trixi.jl and are now available via [TrixiShallowWater.jl](https://github.com/trixi-framework/TrixiShallowWater.jl/). This also affects the related functions `hydrostatic_reconstruction_audusse_etal`, `flux_nonconservative_audusse_etal`, and `FluxHydrostaticReconstruction`. ([#2379])
+
+- The shallow-water equation types `ShallowWaterEquations1D`, `ShallowWaterEquations2D`, and 
+  `ShallowWaterEquationsQuasi1D` have been removed from Trixi.jl and are now available via 
+  [TrixiShallowWater.jl](https://github.com/trixi-framework/TrixiShallowWater.jl/). 
+  This also affects the related functions `hydrostatic_reconstruction_audusse_etal`, 
+  `flux_nonconservative_audusse_etal`, and `FluxHydrostaticReconstruction`. ([#2379])
+- The additional `ìnitial_cache` entries in the caches of `SemidiscretizationHyperbolic`
+  and `SemidiscretizationHyperbolicParabolic`, and the corresponding keyword arguments of
+  their constructors have been removed. ([#2399])
 
 ## Changes in the v0.11 lifecycle
 
diff --git a/docs/literate/src/files/behind_the_scenes_simulation_setup_plots/src/SemidiscretizationHyperbolic_structure_figure.jl b/docs/literate/src/files/behind_the_scenes_simulation_setup_plots/src/SemidiscretizationHyperbolic_structure_figure.jl
index 3066a0a515a..4b2afbb273f 100644
--- a/docs/literate/src/files/behind_the_scenes_simulation_setup_plots/src/SemidiscretizationHyperbolic_structure_figure.jl
+++ b/docs/literate/src/files/behind_the_scenes_simulation_setup_plots/src/SemidiscretizationHyperbolic_structure_figure.jl
@@ -4,7 +4,7 @@ plot(Shape([(-2.3, 4.5), (2.35, 4.5), (2.35, 2.5), (-2.3, 2.5)]), linecolor = "b
      grid = false, xlim = (-2.4, 2.8), ylim = (-25, 5.5))
 annotate!(2.3, 3.5,
           ("SemidiscretizationHyperbolic(mesh, equations, initial_conditions, solver; source_terms,
-boundary_conditions, RealT, uEltype, initial_cache)          ", 10, :black, :right))
+boundary_conditions, RealT, uEltype)                               ", 10, :black, :right))
 annotate!(-2.3, 1.5,
           ("creates and returns SemidiscretizationHyperbolic object, initialized using a mesh, equations,
 initial_conditions, boundary_conditions, source_terms, solver and cache", 9, :black, :left))
diff --git a/src/semidiscretization/semidiscretization_hyperbolic.jl b/src/semidiscretization/semidiscretization_hyperbolic.jl
index c909196b5db..7496a345661 100644
--- a/src/semidiscretization/semidiscretization_hyperbolic.jl
+++ b/src/semidiscretization/semidiscretization_hyperbolic.jl
@@ -53,8 +53,7 @@ end
                                  source_terms=nothing,
                                  boundary_conditions=boundary_condition_periodic,
                                  RealT=real(solver),
-                                 uEltype=RealT,
-                                 initial_cache=NamedTuple())
+                                 uEltype=RealT)
 
 Construct a semidiscretization of a hyperbolic PDE.
 """
@@ -63,12 +62,10 @@ function SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver
                                       boundary_conditions = boundary_condition_periodic,
                                       # `RealT` is used as real type for node locations etc.
                                       # while `uEltype` is used as element type of solutions etc.
-                                      RealT = real(solver), uEltype = RealT,
-                                      initial_cache = NamedTuple())
+                                      RealT = real(solver), uEltype = RealT)
     @assert ndims(mesh) == ndims(equations)
 
-    cache = (; create_cache(mesh, equations, solver, RealT, uEltype)...,
-             initial_cache...)
+    cache = create_cache(mesh, equations, solver, RealT, uEltype)
     _boundary_conditions = digest_boundary_conditions(boundary_conditions, mesh, solver,
                                                       cache)
 
diff --git a/src/semidiscretization/semidiscretization_hyperbolic_parabolic.jl b/src/semidiscretization/semidiscretization_hyperbolic_parabolic.jl
index 1af65995e79..72e79dde853 100644
--- a/src/semidiscretization/semidiscretization_hyperbolic_parabolic.jl
+++ b/src/semidiscretization/semidiscretization_hyperbolic_parabolic.jl
@@ -87,8 +87,7 @@ end
                                           source_terms=nothing,
                                           both_boundary_conditions=(boundary_condition_periodic, boundary_condition_periodic),
                                           RealT=real(solver),
-                                          uEltype=RealT,
-                                          both_initial_caches=(NamedTuple(), NamedTuple()))
+                                          uEltype=RealT)
 
 Construct a semidiscretization of a hyperbolic-parabolic PDE.
 """
@@ -100,12 +99,9 @@ function SemidiscretizationHyperbolicParabolic(mesh, equations::Tuple,
                                                                       boundary_condition_periodic),
                                                # `RealT` is used as real type for node locations etc.
                                                # while `uEltype` is used as element type of solutions etc.
-                                               RealT = real(solver), uEltype = RealT,
-                                               initial_caches = (NamedTuple(),
-                                                                 NamedTuple()))
+                                               RealT = real(solver), uEltype = RealT)
     equations_hyperbolic, equations_parabolic = equations
     boundary_conditions_hyperbolic, boundary_conditions_parabolic = boundary_conditions
-    initial_hyperbolic_cache, initial_cache_parabolic = initial_caches
 
     return SemidiscretizationHyperbolicParabolic(mesh, equations_hyperbolic,
                                                  equations_parabolic,
@@ -113,9 +109,7 @@ function SemidiscretizationHyperbolicParabolic(mesh, equations::Tuple,
                                                  solver_parabolic, source_terms,
                                                  boundary_conditions = boundary_conditions_hyperbolic,
                                                  boundary_conditions_parabolic = boundary_conditions_parabolic,
-                                                 RealT, uEltype,
-                                                 initial_cache = initial_hyperbolic_cache,
-                                                 initial_cache_parabolic = initial_cache_parabolic)
+                                                 RealT, uEltype)
 end
 
 function SemidiscretizationHyperbolicParabolic(mesh, equations, equations_parabolic,
@@ -126,11 +120,8 @@ function SemidiscretizationHyperbolicParabolic(mesh, equations, equations_parabo
                                                boundary_conditions_parabolic = boundary_condition_periodic,
                                                # `RealT` is used as real type for node locations etc.
                                                # while `uEltype` is used as element type of solutions etc.
-                                               RealT = real(solver), uEltype = RealT,
-                                               initial_cache = NamedTuple(),
-                                               initial_cache_parabolic = NamedTuple())
-    cache = (; create_cache(mesh, equations, solver, RealT, uEltype)...,
-             initial_cache...)
+                                               RealT = real(solver), uEltype = RealT)
+    cache = create_cache(mesh, equations, solver, RealT, uEltype)
     _boundary_conditions = digest_boundary_conditions(boundary_conditions, mesh, solver,
                                                       cache)
     _boundary_conditions_parabolic = digest_boundary_conditions(boundary_conditions_parabolic,
@@ -138,11 +129,9 @@ function SemidiscretizationHyperbolicParabolic(mesh, equations, equations_parabo
 
     check_periodicity_mesh_boundary_conditions(mesh, _boundary_conditions)
 
-    cache_parabolic = (;
-                       create_cache_parabolic(mesh, equations, equations_parabolic,
-                                              solver, solver_parabolic, RealT,
-                                              uEltype)...,
-                       initial_cache_parabolic...)
+    cache_parabolic = create_cache_parabolic(mesh, equations, equations_parabolic,
+                                             solver, solver_parabolic, RealT,
+                                             uEltype)
 
     SemidiscretizationHyperbolicParabolic{typeof(mesh), typeof(equations),
                                           typeof(equations_parabolic),