Merge dev into main for major release 0.4 (#832)

Author: efaulhaber

NEWS.md

@@ -4,6 +4,59 @@ TrixiParticles.jl follows the interpretation of
 [semantic versioning (semver)](https://julialang.github.io/Pkg.jl/dev/compatibility/#Version-specifier-format-1)
 used in the Julia ecosystem. Notable changes will be documented in this file for human readability.
 
+## Version 0.4
+
+### API Changes
+
+- Renamed `BoundarySPHSystem` to `WallBoundarySystem` and the keyword argument
+  `movement` to `prescribed_motion`.
+
+- Renamed `OpenBoundarySPHSystem` to `OpenBoundarySystem`.
+
+- Renamed `BoundaryMovement` to `PrescribedMotion`. The `movement_function` must now be
+  a function of `(x, t)` returning the *new position* instead of an offset.
+  For example, `movement_function(t) = SVector(t, 0.0)` now needs to be
+  `movement_function(x, t) = x + SVector(t, 0.0)`.
+
+- Renamed directory `solid` to `structure` in the examples file tree.
+  VTK files for the `TotalLagrangianSPHSystem` are now also called `structure_*`.
+
+- Renamed keyword argument `n_fixed_particles` of the `TotalLagrangianSPHSystem`
+  to `n_clamped_particles`.
+
+- API for `OpenBoundarySystem` and `BoundaryZone` changed.
+  It is now possible to pass multiple `BoundaryZone`s to a single `OpenBoundarySystem`.
+  Reference values are now assigned individually to each `BoundaryZone`. (#866)
+
+- Rename keyword arguments `plane` and `plane_normal` for `BoundaryZone` to `boundary_face` and `face_normal` (#597).
+
+- The argument of `TransportVelocityAdami` is now a keyword argument.
+  `TransportVelocityAdami(1000.0)` now becomes
+  `TransportVelocityAdami(background_pressure=1000.0)` (#884).
+
+- Combined transport velocity formulation (TVF) and particle shifting technique (PST) into
+  one unified framework.
+  The keyword argument `transport_velocity` now changed to `shifting_technique`.
+  The `ParticleShiftingCallback` has been removed. To use PST, use the `UpdateCallback`
+  instead, and pass `shifting_technique=ParticleShiftingTechniqueSun2017()` to the system (#884).
+
+- Renamed the keyword argument `tlsph` to `place_on_shell` for `ParticlePackingSystem`,
+  `sample_boundary`, `extrude_geometry`, `RectangularShape`, and `SphereShape` (#814).
+
+- Custom quantity functions passed to `SolutionSavingCallback` or `PostprocessCallback`
+  that were not using the documented API but were functions of
+  `(system, v_ode, u_ode, semi, t)` now need to be functions
+  of `(system, dv_ode, du_ode, v_ode, u_ode, semi, t)` (#879).
+
+- Renamed `each_moving_particle` to `each_integrated_particle`,
+  `n_moving_particles` to `n_integrated_particles`
+  and `active_particles` to `each_active_particle`.
+
+### Features
+
+- Added consistent particle shifting by Sun et al. (2019) as `ConsistentShiftingSun2019` (#888).
+
+
 ## Version 0.3.2
 ### Features
 
@@ -19,18 +72,17 @@ used in the Julia ecosystem. Notable changes will be documented in this file for
 
 - Fix the coordinates used for TLSPH in Adami extrapolation (#853)
 - Fix PST for small smoothing length factors (#834)
-- The TVF model has been improved to integrate correctly with time stepping (#864) 
- 
+- The TVF model has been improved to integrate correctly with time stepping (#864)
 
 
 ## Version 0.3.1
 
 ### Features
 
-- **Simplified SGS Viscosity Models**: Added ViscosityMorrisSGS and ViscosityAdamiSGS, 
+- **Simplified SGS Viscosity Models**: Added ViscosityMorrisSGS and ViscosityAdamiSGS,
   which implement a simplified Smagorinsky-type sub-grid-scale viscosity. (#753)
 
-- **Multithreaded Integration Array**: Introduced a new array type for CPU backends 
+- **Multithreaded Integration Array**: Introduced a new array type for CPU backends
   that enables multithreaded broadcasting, delivering speed-ups of up to 5× on systems
   with many threads when combined with thread pinning. (#722)
 
@@ -40,17 +92,17 @@ used in the Julia ecosystem. Notable changes will be documented in this file for
 - **DXF file format support**: Import complex geometries using the DXF file format. (#821)
 
 - **Improved Plane interpolation**: Massively improved interpolation performance for planes (#763).
-  
+
 ### GPU
 
  - Make PST GPU-compatible (#813).
-   
+
  - Make open boundaries GPU-compatible (#773).
-  
+
  - Make interpolation GPU-compatible (#812).
 
 ### Important Bugfixes
- 
+
  - Fix validation setups (#801).
 
  - Calculate interpolated density instead of computed density when using interpolation (#808).
@@ -230,4 +282,3 @@ Features:
 #### TLSPH
 
 An implementation of TLSPH (Total Lagrangian Smoothed Particle Hydrodynamics) for solid bodies enabling FSI (Fluid Structure Interactions).
-

Project.toml

@@ -1,7 +1,7 @@
 name = "TrixiParticles"
 uuid = "66699cd8-9c01-4e9d-a059-b96c86d16b3a"
 authors = ["erik.faulhaber <44124897+efaulhaber@users.noreply.github.com>"]
-version = "0.3.2"
+version = "0.4"
 
 [deps]
 Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e"
@@ -27,6 +27,7 @@ RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
+Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 StrideArraysCore = "7792a7ef-975c-4747-a70f-980b88e8d1da"
 TimerOutputs = "a759f4b9-e2f1-59dc-863e-4aeb61b1ea8f"
 TrixiBase = "9a0f1c46-06d5-4909-a5a3-ce25d3fa3284"

docs/src/gpu.md

@@ -2,7 +2,7 @@
 
 GPU support is still an experimental feature that is actively being worked on.
 Currently, the [`WeaklyCompressibleSPHSystem`](@ref), [`TotalLagrangianSPHSystem`](@ref)
-and [`BoundarySPHSystem`](@ref) support GPU execution.
+and [`WallBoundarySystem`](@ref) support GPU execution.
 We have tested GPU support on Nvidia, AMD and Apple GPUs.
 Note that most Apple GPUs do not support `Float64`.
 See [below on how to run single precision simulations](@ref single_precision).

docs/src/refs.bib

@@ -707,6 +707,17 @@ @article{Sun2018
 	doi       = {10.1016/j.cpc.2017.11.016},
 }
 
+@article{Sun2019,
+	author    = {Sun, P. N. and Colagrossi, A. and Marrone, S. and Antuono, M. and Zhang, A. -M.},
+	title     = {A consistent approach to particle shifting in the δPlus-{SPH} model},
+	journal   = {Computer Methods in Applied Mechanics and Engineering},
+	volume    = {348},
+	year      = {2019},
+	issn      = {0045-7825},
+	pages     = {912--934},
+	doi       = {10.1016/j.cma.2019.01.045},
+}
+
 @Article{Tafuni2018,
   author    = {A. Tafuni and J.M. Dom{\'{\i}}nguez and R. Vacondio and A.J.C. Crespo},
   journal   = {Computer Methods in Applied Mechanics and Engineering},

docs/src/systems/boundary.md

@@ -1,15 +1,15 @@
 # Boundary System
 
 ```@docs
-    BoundarySPHSystem
+    WallBoundarySystem
 ```
 
 ```@docs
     BoundaryDEMSystem
 ```
 
 ```@docs
-    BoundaryMovement
+    PrescribedMotion
 ```
 
 
@@ -61,7 +61,7 @@ We provide six options to compute the boundary density and pressure, determined
    This option usually yields the best results of the options listed here.
 2. (Only relevant for FSI) With [`BernoulliPressureExtrapolation`](@ref), the pressure is extrapolated from the
    pressure similar to the [`AdamiPressureExtrapolation`](@ref), but a relative velocity-dependent pressure part
-   is calculated between moving solids and fluids, which increases the boundary pressure in areas prone to
+   is calculated between moving bodies and fluids, which increases the boundary pressure in areas prone to
    penetrations.
 3. With [`SummationDensity`](@ref), the density is calculated by summation over the neighboring particles,
    and the pressure is computed from the density with the state equation.
@@ -103,7 +103,7 @@ Identical to the pressure ``p_b `` calculated via [`AdamiPressureExtrapolation`]
 p_b = \frac{\sum_f (p_f + \frac{1}{2} \, \rho_{\text{neighbor}} \left( \frac{ (\mathbf{v}_f - \mathbf{v}_{\text{body}}) \cdot (\mathbf{x}_f - \mathbf{x}_{\text{neighbor}}) }{ \left\| \mathbf{x}_f - \mathbf{x}_{\text{neighbor}} \right\| } \right)^2 \times \text{factor} +\rho_f (\bm{g} - \bm{a}_b) \cdot \bm{r}_{bf}) W(\Vert r_{bf} \Vert, h)}{\sum_f W(\Vert r_{bf} \Vert, h)}
 ```
 where ``\mathbf{v}_f`` is the velocity of the fluid and ``\mathbf{v}_{\text{body}}`` is the velocity of the body.
-This adjustment provides a higher boundary pressure for solid bodies moving with a relative velocity to the fluid to prevent penetration.
+This adjustment provides a higher boundary pressure for bodies moving with a relative velocity to the fluid to prevent penetration.
 This modification is original and not derived from any literature source.
 
 ```@docs
@@ -217,9 +217,8 @@ a no-slip condition is imposed. When omitting the viscous interaction
 !!! warning
     The no-slip conditions for `BoundaryModelMonaghanKajtar` have not been verified yet.
 
-```@autodocs
-Modules = [TrixiParticles]
-Pages = [joinpath("schemes", "boundary", "monaghan_kajtar", "monaghan_kajtar.jl")]
+```@docs
+    BoundaryModelMonaghanKajtar
 ```
 
 # [Open Boundaries](@id open_boundary)
@@ -234,6 +233,11 @@ Modules = [TrixiParticles]
 Pages = [joinpath("schemes", "boundary", "open_boundary", "boundary_zones.jl")]
 ```
 
+```@autodocs
+Modules = [TrixiParticles]
+Filter = t -> typeof(t) === typeof(TrixiParticles.planar_geometry_to_face)
+```
+
 # [Open Boundary Models](@id open_boundary_models)
 We offer two models for open boundaries, with the choice depending on the specific problem and flow characteristics near the boundary:
 1. [**Method of characteristics**](@ref method_of_characteristics): The method of characteristics is typically used in problems where tracking of wave propagation

docs/src/systems/dem.md

@@ -19,12 +19,12 @@ and moments acting upon them.
 
 ```@autodocs
 Modules = [TrixiParticles]
-Pages = [joinpath("schemes", "solid", "discrete_element_method", "system.jl")]
+Pages = [joinpath("schemes", "structure", "discrete_element_method", "system.jl")]
 ```
 
 ### Contact Models
 
 ```@autodocs
 Modules = [TrixiParticles]
-Pages = [joinpath("schemes", "solid", "discrete_element_method", "contact_models.jl")]
+Pages = [joinpath("schemes", "structure", "discrete_element_method", "contact_models.jl")]
 ```

docs/src/systems/entropically_damped_sph.md

@@ -7,7 +7,8 @@ this scheme uses a pressure evolution equation to calculate the pressure
 ```
 which is derived by [Clausen (2013)](@cite Clausen2013). This equation is similar to the continuity equation (first term, see
 [`ContinuityDensity`](@ref)), but also contains a pressure damping term (second term, similar to density diffusion
-see [`DensityDiffusion`](@ref)), which reduces acoustic pressure waves through an entropy-generation mechanism.
+see [`AbstractDensityDiffusion`](@ref TrixiParticles.AbstractDensityDiffusion)),
+which reduces acoustic pressure waves through an entropy-generation mechanism.
 
 The pressure evolution is discretized with the SPH method by [Ramachandran (2019)](@cite Ramachandran2019) as following:
 
@@ -45,64 +46,3 @@ is a good choice for a wide range of Reynolds numbers (0.0125 to 10000).
 Modules = [TrixiParticles]
 Pages = [joinpath("schemes", "fluid", "entropically_damped_sph", "system.jl")]
 ```
-
-## [Transport Velocity Formulation (TVF)](@id transport_velocity_formulation)
-Standard SPH suffers from problems like tensile instability or the creation of void regions in the flow.
-To address these problems, [Adami (2013)](@cite Adami2013) modified the advection velocity and added an extra term to the momentum equation.
-The authors introduced the so-called Transport Velocity Formulation (TVF) for WCSPH. [Ramachandran (2019)](@cite Ramachandran2019) applied the TVF
-also for the [EDAC](@ref edac) scheme.
-
-The transport velocity ``\tilde{v}_a`` of particle ``a`` is used to evolve the position of the particle ``r_a`` from one time step to the next by
-
-```math
-\frac{\mathrm{d} r_a}{\mathrm{d}t} = \tilde{v}_a
-```
-
-and is obtained at every time-step ``\Delta t`` from
-
-```math
-\tilde{v}_a (t + \Delta t) = v_a (t) + \Delta t \left(\frac{\tilde{\mathrm{d}} v_a}{\mathrm{d}t} - \frac{1}{\rho_a} \nabla p_{\text{background}} \right),
-```
-
-where ``\rho_a`` is the density of particle ``a`` and ``p_{\text{background}}`` is a constant background pressure field.
-The tilde in the second term of the right hand side indicates that the material derivative has an advection part.
-
-The discretized form of the last term is
-
-```math
- -\frac{1}{\rho_a} \nabla p_{\text{background}} \approx  -\frac{p_{\text{background}}}{m_a} \sum_b \left(V_a^2 + V_b^2 \right) \nabla_a W_{ab},
-```
-
-where ``V_a``, ``V_b`` denote the volume of particles ``a`` and ``b`` respectively.
-Note that although in the continuous case ``\nabla p_{\text{background}} = 0``, the discretization is not 0th-order consistent for **non**-uniform particle distribution,
-which means that there is a non-vanishing contribution only when particles are disordered.
-That also means that ``p_{\text{background}}`` occurs as prefactor to correct the trajectory of a particle resulting in uniform pressure distributions.
-Suggested is a background pressure which is in the order of the reference pressure but can be chosen arbitrarily large when the time-step criterion is adjusted.
-
-The inviscid momentum equation with an additional convection term for a particle moving with ``\tilde{v}`` is
-
-```math
-\frac{\tilde{\mathrm{d}} \left( \rho v \right)}{\mathrm{d}t} = -\nabla p +  \nabla \cdot \bm{A},
-```
-
- where the tensor ``\bm{A} = \rho v\left(\tilde{v}-v\right)^T`` is a consequence of the modified
- advection velocity and can be interpreted as the convection of momentum with the relative velocity ``\tilde{v}-v``.
-
-The discretized form of the momentum equation for a particle ``a`` reads as
-
-```math
-\frac{\tilde{\mathrm{d}} v_a}{\mathrm{d}t} = \frac{1}{m_a} \sum_b \left(V_a^2 + V_b^2 \right) \left[ -\tilde{p}_{ab} \nabla_a W_{ab} + \frac{1}{2} \left(\bm{A}_a + \bm{A}_b \right) \cdot \nabla_a W_{ab} \right].
-```
-
-Here, ``\tilde{p}_{ab}`` is the density-weighted pressure
-
-```math
-\tilde{p}_{ab} = \frac{\rho_b p_a + \rho_a p_b}{\rho_a + \rho_b},
-```
-
-with the density  ``\rho_a``,  ``\rho_b`` and the pressure  ``p_a``,  ``p_b`` of particles ``a`` and ``b`` respectively. ``\bm{A}_a`` and ``\bm{A}_b`` are the convection tensors for particle ``a`` and ``b`` respectively and are given, e.g. for particle ``a``, as ``\bm{A}_a = \rho v_a\left(\tilde{v}_a-v_a\right)^T``.
-
-```@autodocs
-Modules = [TrixiParticles]
-Pages = [joinpath("schemes", "fluid", "transport_velocity.jl")]
-```

docs/src/systems/total_lagrangian_sph.md

@@ -59,7 +59,7 @@ The term $\bm{f}_a^{PF}$ is an optional penalty force. See e.g. [`PenaltyForceGa
 
 ```@autodocs
 Modules = [TrixiParticles]
-Pages = [joinpath("schemes", "solid", "total_lagrangian_sph", "system.jl")]
+Pages = [joinpath("schemes", "structure", "total_lagrangian_sph", "system.jl")]
 ```
 
 ## Penalty Force
@@ -100,7 +100,7 @@ where the error vector is defined as
 
 ```@autodocs
 Modules = [TrixiParticles]
-Pages = [joinpath("schemes", "solid", "total_lagrangian_sph", "penalty_force.jl")]
+Pages = [joinpath("schemes", "structure", "total_lagrangian_sph", "penalty_force.jl")]
 ```
 
 ## Viscosity