A well-balanced weakly compressible SPH formulation for free-surface flows and its GPU implementation

Abstract

This study proposes a well-balanced formulation of weakly compressible smoothed particle hydrodynamics (WCSPH) for free-surface flows, which preserves hydrostatic equilibrium exactly at the discrete level--a property essential for reliable long-term simulations. Although well-balanced schemes are well established for mesh-based methods, the property remains largely unaddressed in WCSPH, where the particle approximation of the pressure gradient fails to balance the gravitational force exactly. The imbalance stems from two difficulties: the nonlinearity of the pressure-gradient-over-density term, and the approximation error of gradients evaluated by particle summation. The first is resolved by introducing an auxiliary potential variable that recasts the nonlinear term as the gradient of a single scalar, which reduces to a linear function of position under hydrostatic conditions. The second is resolved by a Riemann-based gradient approximation with kernel correction, which is first-order consistent and recovers linear fields exactly. These two ingredients ensure that the discrete potential gradient balances gravitational force exactly. Widely used techniques, including $δ-$SPH, particle shifting and tensile instability control, are readily incorporated. The formulation is further extended to three dimensions and implemented on GPU with architecture-tailored optimizations. Hydrostatic tests with rectangular, triangular and Gaussian bottom topographies show that the proposed formulation attains the well-balanced property to machine precision, reducing the spurious velocity error of conventional SPH from $10^{-3}$ to the order of $10^{-13}$. More complex benchmarks confirm its robustness, accuracy and low pressure oscillation, with simulations of up to 17.53 million particles performed on a single consumer-grade GPU.