Add argument t to source_terms (#582)

* add time to source terms

* remove example file

examples/fluid/oscillating_drop_2d.jl

@@ -19,7 +19,7 @@ sigma = 0.5
 # Make this a constant because global variables in the source terms are slow
 const OMEGA = 1.0
 
-source_terms = (coords, velocity, density, pressure) -> -OMEGA^2 * coords
+source_terms = (coords, velocity, density, pressure, t) -> -OMEGA^2 * coords
 
 # 1 period in the exact solution as computed below (but integrated with a small timestep)
 period = 4.567375

src/general/semidiscretization.jl

@@ -451,7 +451,8 @@ function kick!(dv_ode, v_ode, u_ode, semi, t)
         @trixi_timeit timer() "system interaction" system_interaction!(dv_ode, v_ode, u_ode,
                                                                        semi)
 
-        @trixi_timeit timer() "source terms" add_source_terms!(dv_ode, v_ode, u_ode, semi)
+        @trixi_timeit timer() "source terms" add_source_terms!(dv_ode, v_ode, u_ode,
+                                                               semi, t)
     end
 
     return dv_ode
@@ -513,7 +514,7 @@ function update_nhs!(semi, u_ode)
     end
 end
 
-function add_source_terms!(dv_ode, v_ode, u_ode, semi)
+function add_source_terms!(dv_ode, v_ode, u_ode, semi, t)
     foreach_system(semi) do system
         dv = wrap_v(dv_ode, system, semi)
         v = wrap_v(v_ode, system, semi)
@@ -522,7 +523,7 @@ function add_source_terms!(dv_ode, v_ode, u_ode, semi)
         @threaded system for particle in each_moving_particle(system)
             # Dispatch by system type to exclude boundary systems
             add_acceleration!(dv, particle, system)
-            add_source_terms_inner!(dv, v, u, particle, system, source_terms(system))
+            add_source_terms_inner!(dv, v, u, particle, system, source_terms(system), t)
         end
     end
 
@@ -544,13 +545,13 @@ end
     return dv
 end
 
-@inline function add_source_terms_inner!(dv, v, u, particle, system, source_terms_)
+@inline function add_source_terms_inner!(dv, v, u, particle, system, source_terms_, t)
     coords = current_coords(u, system, particle)
     velocity = current_velocity(v, system, particle)
     density = particle_density(v, system, particle)
     pressure = particle_pressure(v, system, particle)
 
-    source = source_terms_(coords, velocity, density, pressure)
+    source = source_terms_(coords, velocity, density, pressure, t)
 
     # Loop over `eachindex(source)`, so that users could also pass source terms for
     # the density when using `ContinuityDensity`.
@@ -561,7 +562,7 @@ end
     return dv
 end
 
-@inline add_source_terms_inner!(dv, v, u, particle, system, source_terms_::Nothing) = dv
+@inline add_source_terms_inner!(dv, v, u, particle, system, source_terms_::Nothing, t) = dv
 
 @doc raw"""
     SourceTermDamping(; damping_coefficient)
@@ -592,7 +593,7 @@ struct SourceTermDamping{ELTYPE}
     end
 end
 
-@inline function (source_term::SourceTermDamping)(coords, velocity, density, pressure)
+@inline function (source_term::SourceTermDamping)(coords, velocity, density, pressure, t)
     (; damping_coefficient) = source_term
 
     return -damping_coefficient * velocity

src/schemes/boundary/open_boundary/system.jl

@@ -165,8 +165,6 @@ end
 
 update_callback_used!(system::OpenBoundarySPHSystem) = system.update_callback_used[] = true
 
-@inline source_terms(system::OpenBoundarySPHSystem) = nothing
-
 @inline hydrodynamic_mass(system::OpenBoundarySPHSystem, particle) = system.mass[particle]
 
 @inline function particle_density(v, system::OpenBoundarySPHSystem, particle)

src/schemes/fluid/entropically_damped_sph/system.jl

@@ -31,7 +31,7 @@ See [Entropically Damped Artificial Compressibility for SPH](@ref edac) for more
 - `buffer_size`:    Number of buffer particles.
                     This is needed when simulating with [`OpenBoundarySPHSystem`](@ref).
 - `source_terms`:   Additional source terms for this system. Has to be either `nothing`
-                    (by default), or a function of `(coords, velocity, density, pressure)`
+                    (by default), or a function of `(coords, velocity, density, pressure, t)`
                     (which are the quantities of a single particle), returning a `Tuple`
                     or `SVector` that is to be added to the acceleration of that particle.
                     See, for example, [`SourceTermDamping`](@ref).

src/schemes/fluid/weakly_compressible_sph/system.jl

@@ -31,7 +31,7 @@ See [Weakly Compressible SPH](@ref wcsph) for more details on the method.
                     This is needed when simulating with [`OpenBoundarySPHSystem`](@ref).
 - `correction`:     Correction method used for this system. (default: no correction, see [Corrections](@ref corrections))
 - `source_terms`:   Additional source terms for this system. Has to be either `nothing`
-                    (by default), or a function of `(coords, velocity, density, pressure)`
+                    (by default), or a function of `(coords, velocity, density, pressure, t)`
                     (which are the quantities of a single particle), returning a `Tuple`
                     or `SVector` that is to be added to the acceleration of that particle.
                     See, for example, [`SourceTermDamping`](@ref).

test/general/semidiscretization.jl

@@ -144,7 +144,7 @@
         # Avoid `SystemBuffer` barrier
         TrixiParticles.each_moving_particle(system::Union{System1, System2}) = TrixiParticles.eachparticle(system)
 
-        TrixiParticles.add_source_terms!(dv_ode, v_ode, u_ode, semi)
+        TrixiParticles.add_source_terms!(dv_ode, v_ode, u_ode, semi, 0.0)
 
         dv1 = TrixiParticles.wrap_v(dv_ode, system1, semi)
         @test dv1 == -0.1 * v1