Stress balance for a viscous flow with a single rolling particle

One of the most important aspects in hydraulic engineering is to describe flows over mobile porous media in a continuum sense to derive models for sediment transport. This remains a challenging task due to the complex coupling of the particle and the fluid phase. Computational Fluid Dynamics can provide the data needed to understand the coupling of the two phases. To this end, we carry out grain-resolving Direct Numerical Simulations of multiphase flow. The particle phase is introduced by the Immersed Boundary Method and the particle-particle interaction is described by a sophisticated Discrete Element Method. We derive the stress budgets of the fluid and the particle phases separately through a rigorous analysis of the governing equations using the Double Averaging Methodology and the Coarse-Graining Method. As a next step, we perform a simple simulation of a heavy particle exposed to a Poiseuille flow rolling along a wall to understand the physical implications of the fluid-particle coupling. All terms of the stress balances can be computed in a straightforward manner allowing to close the budgets for the two phases separately. However, we encounter problems when attempting to combine the fluid-resolved local stresses with the coarse-grained particle stresses into a single balance for the fluid-particle mixture.


Introduction
Shallow flows over fixed or mobile granular beds in open channels are of vital interest to hydraulic engineering.Despite its importance, it remains a challenging task to scale physical processes such as erosion and bed-load transport from the grain scale up to the field scale.Recently, this topic has seen a big push through the formulation of the Double Averaging Methodology [DAM, 1,2].In particular, applying this concept to datasets generated by particle-resolving Direct Numerical Simulations (DNS) employing the Immersed Boundary Method (IBM) has provided great insights into the fundamental physical mechanisms governing of complex multiphase flows [e.g 3, 4].These studies have revealed the importance of the nonlinear terms of the double-averaged momentum budget.However, it was not possible to close the momentum budget entirely due to the difficulties arising from the computation of the interfacial stress terms for pressure and viscous drag.
The present study addresses this issue.We use DAM principles to write the stress balance of the fluid phase for the IBM framework by directly accounting for the interfacial stresses through a rigorous coupling of the fluid and the particle phase.We then employ the Coarse-Graining Method [CGM,5] to derive the analogous stress balance for the particle phase.The two budgets are then computed for a simple test case of a heavy particle exposed to a laminar Poiseuille flow, which is rolling along the wall in the streamwise direction.Understanding and closing these budgets will help us to understand the importance of the different parts of the stress budget to build higher order closures within the DAM framework [6] or to employ constitutive relations for the granular phase, such as the µ(I)-rheology [7].

Fully-resolved simulations of particle-laden flows
We solve the unsteady Navier-Stokes equations for an incompressible Newtonian fluid everywhere in the domain, including within the particles: as well as the continuity equation, ∇ • u = 0, on a uniform rectangular grid with grid cell size ∆x = ∆y = ∆z = h.Here, u = (u, v, w) T designates the fluid velocity vector in Cartesian components, ρ f the fluid density, p the pressure, µ f the dynamic viscosity, t the time, and f IBM an artificial volume force introduced by the IBM [8,9].This volume force, which acts on the right-hand side of (1) in the vicinity of the inter-phase boundaries, connects the motion of the particles to the fluid phase.The volume force f b drives the flow and can be interpreted as a pressure gradient.Within the framework of the IBM, we calculate the motion of each individual spherical particle by solving an ordinary differential equation for its translational velocity and its angular velocity Here, m p is the particle mass, Ω − is the part of the fluid domain inside the particle and L is the shell volume surrounding the particle surface as shown in Figure 1a.Furthermore, ρ p is the particle density, V p the particle volume, g the gravitational acceleration, I p = 8πρ p R 5 p /15 the moment of inertia, R p the particle radius, r = x − x p is the position vector of the surface point with respect to the center of mass x p of a particle, and F c,p and T c,p are the force and torque due to particle collisions, respectively.For brevity, we define F r,p and T r,p to be the rigid body force and torque, F IBM,p and T IBM,p to be the IBM force and torque, and F g,p to be the buoyancy force acting on particle p.
The computational approach for modeling the cohesionless particle-particle interaction has been described in detail in [10] and [11] including a validation of normal and oblique binary collisions as well as for the collective motion of a sediment bed sheared by a Poiseuille flow.The particle-particle interaction is comprised of short-range hydrodynamic effects due to lubrication forces F ln and F lt acting in normal and tangential direction as well as contact forces acting in the normal and tangential direction, denoted as F n and F t , respectively.The resulting collision force on particle p is the sum off all these effects: where the subscripts pq and pw indicate interactions with particle q or a wall, respectively.For all the simulations presented in the following, we have chosen a common parameterization of the different parts of the collision model to reflect the properties of silicate materials [e.g. 3, 4, 10, 12].

Simulation scenario
We consider a sphere rolling along the bottom wall of a channel in a Poiseuille flow.The relative particle density was set to ρ p /ρ f = 2.1, the spatial domain was chosen to be The particle was discretized by 24 grid cells per diameter.
In the absence of the particle, there would be a laminar flow with a bulk velocity of u ref = −L 2 y /(12µ f ) f b,x , a channel height of L y , and a Reynolds number of Re ref = ρ f u ref L y /µ f = 10.The flow was initialized with the analytical solution of a laminar, undisturbed Poiseuille flow.The particle velocity was initially set to the spatially averaged analytical velocity field of the fluid within the particle.The presence of the particle, however, changes the bulk velocity and Reynolds number.At this low Reynolds number, the flow remains laminar, but takes some time to develop as the presence of the particle constricts and slows the flow (Figure 1b).We

Momentum balance for the fluid phase
For our stress analysis, we are interested in analyzing a control volume that extends from the top wall to an arbitrary height y in the vertical dimension, and encompasses the entire domain in the streamwise and spanwise directions.We can write the integral form of (1) as where CV comprises all surfaces of the control volume shown in Figure 2a and n + is the normal vector pointing outwards from Ω + CV .In addition, we recast the pressure and viscous terms using the fluid stress tensor, τ + = −pI + µ f (∇u + (∇u) T ), where I is the identity matrix.Here we neglect the effect of the immersed boundary force, assuming that it is implicitly handled by the fluid stress at the fluid/particle interface.We do, however, include the forcing term f b , which represents the background pressure gradient used to drive the flow.In a steady-state configuration, the time-dependent term evaluates to zero, allowing us to balance the stresses by looking at the advective terms and stresses at the boundaries.
All of these terms except for the fluid stress at the particle surface are straightforward to calculate.However, we can evaluate the fluid stress indirectly using the IBM force, as was done to obtain the particle equations of motion (2) and (3).That is, the IBM force acts as a jump in stress: where we are careful to distinguish between n + , the outward surface normal for the volume Ω + CV , and n − , the outward surface normal for the volume Ω − CV , which point in opposite directions.
To evaluate the fluid stress τ − • n − at the particle interior, we can perform a stress balance on Ω − CV .The integral form of the Navier-Stokes equations together with divergence theorem give us Using ( 7) together with ( 6) and ( 5), we obtain where Using the definition of the horizontal, superficial average [2], we can rewrite (8) for the x-velocity component as where γ is an indicator function for the fluid volume fraction (γ = 1 outside the particle and γ = 0 inside the particle) and φ is an indicator function for the particle volume fraction (φ = 1−γ).We have also used the fact that µ f , ρ f , and f b,x are constant throughout the domain.
Note that this form is equivalent to the double-averaged momentum balance analyzed in [4] since there are no fluctuations or unsteady effects for the present scenario of a laminar flow, but the interfacial terms of the DAM framework had been replaced by the particle stress.

Momentum balance for the particle phase
Although the fluid stress analysis accounts for the effects of the particles through the IBM force, we can perform an analysis on the particle phase by itself.We can apply the CGM based on the works of [13] and [5] to (2).The CGM conserves quantities of interest by smoothing out discrete variables of the particle-phase in a continuum field.This method can represent smooth fields even for instantaneous data to obtain: where E3S Web of Conferences 40, 04003 (2018) https://doi.org/10.1051/e3sconf/20184004003River Flow 2018 is the coarse-grained local particle acceleration, and are the coarse-grained forces acting on the particle center of mass through (11).The coarsegraining function W(r) spreads the particle-centered quantities onto a Eulerian mesh, allowing us to treat them similar to a continuum field.The coarse-graining function W plays a very similar role to that of the delta functions used in the IBM: smoothly spreading a quantity from one mesh to another.The main property is that R 3 W(r) dr = 1.We implemented a coarse-graining function based on the Dirac Delta function of [14].We could coarse-grain the collision forces in a manner similar to that of [5], which allows us to evaluate the entire stress tensor for collisions.However, the information on hydrodynamic forces is given at the particle centers.Hence, for the sake of consistency, we applied the CGM to all quantities at the particle center.We have also found the kinetic stress, σ k , mentioned by [5] to be negligible in these simulations.Similar to the fluid stresses, we can analyze the coarse-grained particle forces within a control volume spanning the entire domain in the streamwise and spanwise directions and extending from the top wall to an arbitrary height y.Integrating (11) over this volume, we obtain We can again apply the averaging operator (9) to recast (14) as a line integral in the wallnormal direction: If the particles are in a steady state, then the acceleration term would be zero.

Results
We now investigate the momentum balance relationships (10) and (15) for the single rolling sphere, which are shown in Figure 3.The reference stress for these plots is the wall shear stress for the reference case, σ ref = µ f du/dy| y=0 .Thus, in the absence of the particle, the external stress in Figure 3a would extend from σ x /σ ref = −1 at y/L y = 1 to σ x /σ ref = 1 at y/L y = 0.The presence of the particle causes this curve to shift to the right, decreasing the magnitude of the stress at the top wall and increasing the stress at the bottom wall, most of which is taken up by the particle through collision forces.
We can see in Figure 3a, which represents instantaneous data at t = 2.5t ref , that the stresses arising from the fluid and particle phases, which act as a momentum sink, are in perfect equilibrium with the external stress, which drives the flow.Similarly, the particle stress balance in Figure 3b is also in equilibrium, the hydrodynamic forces driving the particle in the positive x-direction balanced by the collision forces acting in the negative x-direction, indicate that the particle is not accelerating.As expected, Figures 3c and 3d show us that the convective terms in (10) are negligible.Thus, the viscous term accounts for the fluid stress, and the viscous and IBM terms account for the particle stress.
The particle stress in Figure 3a and the hydrodynamic stress in Figure 3b should represent the same quantity, the former representing the effect of the particle acting on the fluid and the which is plotted in Figure 4.While the total stress does match the external stress in the clear fluid layer above the particle and at the lower wall, they do not match within the particle region.This result is due to the fact that the coarse-grained collision stress represents the net stress acting on the particle (at its center).As an illustration, the fluid phase particle stress, which represents the local hydrodynamic interactions that occur along the particle surface, has been included in the plot.The collision stress, which perfectly balances the particle stress at the lower wall, once the entire domain is considered, is simply smeared over the volume by the coarse-graining method and cannot recover local variations.Thus, using hybrid fluid/coarse-grained quantities can lead to local gaps in the stress balance, but will provide consistent information once the entire particle space is considered.

Conclusion
Two approaches to compute the stress balances for the test case of a single particle rolling along a wall in a Poiseuille flow were presented.We have successfully computed all terms of both balances.Since the out-of-balance terms were shown to be zero, the proposed set

2 E3SFigure 1 .
Figure 1.a) Sketch of the computational setup together with an indication of the different parts of the fluid domain as they enter Equations (2) and (3).b) Simulation result of a contour plot cutting through the particle center.

Figure 2 .
Figure 2. Control volumes for fluid (a) surrounding the particle and (b) within the particle.

E3S−Figure 3 .
Figure 3. Stress balance in the x-direction for the fluid/particle mixture, given by (16).The total stress is the sum of the fluid and collision stresses.