Additional work on parabolic terms

src/basic_types.jl

@@ -70,5 +70,19 @@ const boundary_condition_periodic = BoundaryConditionPeriodic()
 
 Base.show(io::IO, ::BoundaryConditionPeriodic) = print(io, "boundary_condition_periodic")
 
+"""
+    boundary_condition_do_nothing = BoundaryConditionDoNothing()
+
+Imposing no boundary condition just evaluates the flux at the inner state.
+"""
+struct BoundaryConditionDoNothing end
+
+@inline function (boundary_condition::BoundaryConditionDoNothing)(inner_flux_or_state, other_args...)
+  return inner_flux_or_state
+end
+
+const boundary_condition_do_nothing = BoundaryConditionDoNothing()
+
+Base.show(io::IO, ::BoundaryConditionDoNothing) = print(io, "boundary_condition_do_nothing")
 
 end # @muladd

src/equations/laplace_diffusion_2d.jl

@@ -12,6 +12,9 @@ end
 LaplaceDiffusion2D(diffusivity, equations) =
   LaplaceDiffusion2D{typeof(equations), nvariables(equations), typeof(diffusivity)}(diffusivity, equations)
 
+varnames(variable_mapping, equations_parabolic::LaplaceDiffusion2D) =
+  varnames(variable_mapping, equations_parabolic.equations)
+
 # no orientation specified since the flux is vector-valued
 function flux(u, grad_u, equations::LaplaceDiffusion2D)
   dudx, dudy = grad_u

src/solvers/dgmulti/dg.jl

@@ -31,9 +31,8 @@ mul_by_accum!(A::UniformScaling) = MulByAccumUniformScaling()
 # solution storage formats.
 @inline apply_to_each_field(f::MulByUniformScaling, out, x, args...) = copy!(out, x)
 @inline function apply_to_each_field(f::MulByAccumUniformScaling, out, x, args...)
-  # TODO: DGMulti speed up using threads
-  for (i, x_i) in enumerate(x)
-    out[i] = out[i] + x_i
+  @threaded for i in eachindex(x)
+    out[i] = out[i] + x[i]
   end
 end
 

test/test_parabolic_2d.jl

@@ -30,6 +30,7 @@ isdir(outdir) && rm(outdir, recursive=true)
     semi = SemidiscretizationHyperbolicParabolic(mesh, equations, equations_parabolic, initial_condition, dg)
     @test_nowarn_debug show(stdout, semi)
     @test_nowarn_debug show(stdout, MIME"text/plain"(), semi)
+    @test_nowarn_debug show(stdout, boundary_condition_do_nothing)
 
     @test nvariables(semi)==nvariables(equations)
     @test Base.ndims(semi)==Base.ndims(mesh)
@@ -40,15 +41,17 @@ isdir(outdir) && rm(outdir, recursive=true)
     Trixi.compute_coefficients!(u0, 0.0, semi)
     @test u0 ≈ ode.u0
 
+    # test "do nothing" BC just returns first argument
+    @test boundary_condition_do_nothing(u0, nothing) == u0
+
     @unpack cache, cache_parabolic, equations_parabolic = semi
     @unpack u_grad = cache_parabolic
     for dim in eachindex(u_grad)
       fill!(u_grad[dim], zero(eltype(u_grad[dim])))
     end
 
     t = 0.0
-    # pass in `boundary_condition_periodic` to fake "do-nothing"
-    # TODO: DGMulti. Make `boundary_condition_do_nothing` a callable struct like BoundaryConditionPeriodic
+    # pass in `boundary_condition_periodic` to skip boundary flux/integral evaluation
     Trixi.calc_gradient!(u_grad, ode.u0, t, mesh, equations_parabolic,
                          boundary_condition_periodic, dg, cache, cache_parabolic)
     @unpack x, y = mesh.md