FreeFlight

In free flight, the acceleration of the droplets is computed as \begin{align} \frac{ d \mathbf{v}_{drop} }{dt} = \hat{\mathbf{g}}\end{align} where \( \hat{\mathbf{g}}\) denotes an effective gravity specified by the user via the gravity command.

Coupling terms

If in interaction with airflow (or another medium), this effective acceleration would typically take the form \begin{align} \hat{\mathbf{g}} = \frac{\mathbf{F}_D + \mathbf{F}_p}{m_{drop}} + \mathbf{g}\end{align} with \( \mathbf{F}_D\), \( \mathbf{F}_p\), \( \mathbf{g}\) denoting the drag force, the force due to pressure gradients in the airflow and the classical gravitational acceleration, respectively. For the drag force one may use the classical drag equation \begin{align} \mathbf{F}_D = \frac{1}{2} C_D \rho_{air} A_{drop} \left( \mathbf{v}_{air} - \mathbf{v}_{drop} \right) \left\| \mathbf{v}_{air} - \mathbf{v}_{drop} \right\|\end{align} with \( C_D\), \( A_{drop}\), \( \mathbf{v}_{air}\) denoting the drag coefficient, the droplet cross-sectional area and the air velocity at the droplet position, respectively. The pressure gradient force, on the other hand, is given by \begin{align} \mathbf{F}_p = -\nabla p_{air}V_{drop}\end{align} with \( \nabla p_{air}\), \( V_{drop}\) denoting the pressure gradient in the air phase at the droplet location and the droplet volume, respectively. Now, taking into account the spherical shape of particles in DROPLETPHASE, this leads to \begin{align} \hat{\mathbf{g}} = \frac{3}{4}C_D \frac{\rho_{air}}{\rho_{drop}} \frac{1}{d_{drop}} \left( \mathbf{v}_{air} - \mathbf{v}_{drop} \right) \left\| \mathbf{v}_{air} - \mathbf{v}_{drop} \right\| - \frac{1}{\rho_{drop}} \nabla p_{air} + \mathbf{g}\end{align} We note, that \( \mathbf{v}_{air}\) and \( \nabla p_{air}\) are typically calculated via the approxY() functionality. In order to compute the velocity contribution of the acceleration due to drag in the actual time step we could apply the analytical integration of it: \begin{align} \mathbf{v}^{n+1}_{drop} = \mathbf{v}^{n}_{drop} - \left( \mathbf{v}_{air}-\mathbf{v}^{n}_{drop} \right) \left( 1 - e^{ - dt (c_D \rho_{air} || \mathbf{v}_{air}-\mathbf{v}^{n}_{drop} || ) / (2 \rho_{drop} d_{drop}) } \right)\end{align} In order to compute the velocity contribution of this acceleration in the actual time step with the analytical integration of it the user may define:

Stability and (semi-)implicit treatment of drag

Every term included into the gravity will be treated as an explicit term in the time integration of droplet velocity. Clearly, for the velocity-dependent drag force, this may become a source of instabilities within the simulation. This expresses itself in the droplet velocity oscillating around the background (air) velocity instead of exponentially decaying towards it from either above or below. To avoid these problems, there are two typical approaches:

• Keeping the explicit integration of drag, but restricting the time step (using DELT_dt_AddCond) to a fraction of the particle relaxation time \( t_r\) given by \begin{align} \frac{ d \mathbf{v}_{drop} }{dt} = \beta \left( \mathbf{v}_{air} - \mathbf{v}_{drop} \right) \ \ \Rightarrow\ \ t_r = \frac{1}{\beta}\end{align}
• Using DarcyConstant instead of gravity to specify the drag force in such a way that MESHFREE can use the analytical exponential result obtained under the assumption of a constant \( \beta\) within the time step, i.e. \begin{align} \frac{ d \mathbf{v}_{drop} }{dt} = \beta \left( \mathbf{v}_{air} - \mathbf{v}_{drop} \right) + \mathbf{g}^{\ast} \ \ \Rightarrow \ \ \mathbf{v}_{drop}^{n+1} = \mathbf{K}^n + (\mathbf{v}_{drop}^{n}-\mathbf{K}^n)e^{-\Delta t\beta^n}\quad\text{where}\quad \mathbf{K}^n=\mathbf{v}_{air}^n+\frac{\mathbf{g}^{\ast,n}}{\beta^n}\end{align} where \( \mathbf{g}^{\ast}\) subsumes all terms of \( \hat{\mathbf{g}}\) except drag (which is no longer supplied via gravity).

For an implementation of both approaches in MESHFREE input files, consult the ChannelWithFilter example.