For droplets at boundaries with
%IDENT_slip% or
%IDENT_ThinFilm%, a thin film / liquid layer is computed. This is done by calculating the film height and curvature from an SPH-like approximation, based on the volumes of droplets that have accumulated around any given point at the boundary. An example is given in
LiquidLayerExample
Droplet acceleration
For droplets collecting on a surface and forming a liquid layer, the acceleration of the droplets is computed as
\begin{align} \frac{ d \mathbf{v}_{drop} }{dt} = \dot{\mathbf{v}}^{normal}_{viscous} \hspace{5pt} + \hspace{5pt} \dot{\mathbf{v}}^{tangential}_{viscous} \hspace{5pt} + \hspace{5pt} \dot{\mathbf{v}}_{\nabla p} \hspace{5pt} + \hspace{5pt} \dot{\mathbf{v}}_{wetting} \hspace{5pt} + \hspace{5pt} \hat{\mathbf{g}}\end{align}
where
-
- \( \dot{\mathbf{v}}^{normal}_{viscous}\) : is the deceleration due to viscous forces which act in the boundary layer of the film in normal direction to the wall geometry.
- \( \dot{\mathbf{v}}^{tangential}_{viscous}\) : is the deceleration due to viscous forces due to velocity gradients in tangential direction to the wall geometry.
- \( \dot{\mathbf{v}}_{\nabla p}\) : is the acceleration due to pressure gradient in the liquid film.
- \( \dot{\mathbf{v}}_{wetting}\) : is the acceleration due to additional surface tension forces (via young relation) that originate from diffrences in actual wetting angle to the defined wetting angle in equilibrium.
- \( \hat{\mathbf{g}}\) : Acceleration defined by the user. This can be simply the gravity but also include acceleration with other chambers (here air) with media like air or others. Then this could include shear forces and pressure gradient from this chamber.
The resulting acceleration is stored in
%ind_v_dot(1)% ,
%ind_v_dot(2)% and
%ind_v_dot(3)% .
Additionally, the accelerations from
DropletCollisions will be applied if defined. Also accelerations due to
DarcyConstant \( \beta > 0\) (as described in
FreeFlight ) will be applied if specified. If the darcy constant is beta>0, then accelerations from
DropletCollisions will be treated in a (semi-) implicit treatment as described in
FreeFlight .
Viscous forces in normal direction to the boundary layer
The deceleration due to viscous forces which act in the boundary layer of the film in normal direction to the wall geometry is defined as:
\begin{align} \dot{\mathbf{v}}^{normal}_{viscous} = \frac{\eta^{normal}_{drop}}{\rho_{drop}} \frac{\mathbf{v}_{geometry}-\mathbf{v}_{drop}}{H^2_{film}}\end{align}
with
-
- \( \eta^{normal}_{drop}\) : the viscosity normal to the wall, which is set by %MED_LIQUID_FILM%
- \( \rho_{drop}\) : the density of the droplet ( %ind_r% )
- \( \mathbf{v}_{geometry}\) : the velocity of the geometry
- \( \mathbf{v}_{drop}\) : the velocity of the droplet ( %ind_v(1)% , %ind_v(2)%, %ind_v(3)% )
- \( H_{film}\) : the film height of the liquid layer computed in an SPH-way
\begin{align} H_{film,i} = \sum \limits_j W_{ij} V_j n_j\end{align}
where
-
- \( V_j = \frac{\pi}{6} d_{drop,j}^3\) is the volume of the droplet (volume package represented by the MESHFREE point).
- \( n_j =\) Multiplicity %ind_mult%
- \( W_{ij} = \frac{\alpha}{\pi} \frac{1}{h_{ij}^2} \frac{1}{1-\exp(-\alpha)} \exp\left( -\alpha \frac{r_{ij}^2}{h^2} \right)\)
- \( r_{ij}^2 = x_{ij}^2+y_{ij}^2+z_{ij}^2 = (x_i-x_j)^2+(y_i-y_j)^2+(z_i-z_j)^2\)
- \( h_{ij} = \frac{h_i+h_j}{2}\)
- in general we have the approximation kernel \( \Pi u = \sum \limits_j W_{ij} u_j \frac{V_jn_j}{H_{film,j}}\) thus the layer thickness is an eigenfunction of the approximation kernel.
The film height is stored in
%ind_hwf% , in case, a boundary element shall be initialized with a thin film,
%IDENT_ThinFilm% can be used.
In order to compute the velocity contribution of the acceleration \( \dot{\mathbf{v}}^{normal}_{viscous}\) in the actual time step we apply the analytical integration of it:
\begin{align} \mathbf{v}^{viscous,n+1}_{drop,normal} = \mathbf{v}^{n}_{drop} - \left( \mathbf{v}_{geometry}-\mathbf{v}^{n}_{drop} \right) \left( 1 - e^{ - dt (2 \eta^{normal}_{drop})/(\rho_{drop} H^2_{film}) } \right)\end{align}
Viscous forces in tangential direction to the wall
The deceleration due to viscous forces due to velocity gradients in tangential direction to the wall geometry is defined as
\begin{align} \dot{\mathbf{v}}^{tangential}_{viscous} = \frac{\eta^{tangential}_{drop}}{\rho_{drop}} \Delta(\mathbf{v})\end{align}
with
-
- \( \eta^{tangential}_{drop}\) : the viscosity in tangential direction relative to the wall, which is set by %MED_LIQUID_FILM%
- \( \rho_{drop}\) : the density of the droplet ( %ind_r% )
- \( \Delta(\mathbf{v})\) : the curvature (Laplacian) of the velocity in tangential direction to the wall (stored in %ind_lap_hwf% ).
Acceleration due to the pressure gradient in the liquid film
The pressure gradient acceleration has several contributions:
\begin{align} \dot{\mathbf{v}}_{\nabla p} = \nabla \left( p_g + p_{\sigma,g} + p_{\sigma,f} + p_{dyn} \right)\end{align}
where
-
- \( p_g\) : is the hydrostatic pressure due to gravity and weight of the droplets (depth pressure)
- \( p_{\sigma,g}\) : is the hydrostatic pressure on the geometry side of the liquid film due to surface tension and curvature
- \( p_{\sigma,f}\) : is the hydrostatic pressure on the free surface side of the liquid film due to surface tension and curvature
- \( p_{dyn}\) : is the dynamic pressure due to velocity gradients
For the pressure due to gravity and weight of the droplets (depth pressure) we define:
\begin{align} p_g = - \rho_{drop}(\hat{\mathbf{g}}\cdot\mathbf{n})H_{film}\end{align}
where \( \hat{\mathbf{g}}\) is described in a separate section below.
For the pressure difference from curvature of the geometry and the curvature of the liquid film we define \( p_{\sigma,g} = \sigma \kappa_g\) and \( p_{\sigma,f} = \sigma \kappa_f\) respectively with
-
- \( \sigma\): surface tension (see sigma, %ind_SIG%)
- \( \kappa_f = \Delta H_{film}\) : surface curvature of the film at the free surface side
- \( \kappa_g = \frac{\partial\mathbf{a}}{\partial\mathbf{n}}\cdot\mathbf{a} + \frac{\partial\mathbf{b}}{\partial\mathbf{n}}\cdot\mathbf{b}\) : curvature of the film at the geometry side (see %ind_lap_geometry%)
- \( \mathbf{a}, \mathbf{b}\): Tangential vectors determined on the basis of the normal \( \mathbf{n}\)
The curvature of the liquid film can be computed by the second derivatives, as for example
\begin{align} \Delta(H_{film,i}) = \sum \limits_j W_{ij} \left( -\frac{2\alpha}{h^2} \right) V_j n_j + \sum \limits_j W_{ij} \left( -\frac{2\alpha}{h^2} \right)^2 x_{ij}^2 V_j n_j\end{align}
The dynamic pressure describes acceleration forces due to velocity gradients and is defined as:
\( p_{dyn} = - \rho_{drop}(\hat{\mathbf{g}}_{eff}\cdot\mathbf{n})H_{film}\)
with
-
- \( \hat{\mathbf{g}}_{eff} = \left(\frac{\partial\mathbf{v}_{rel}}{\partial\mathbf{n}}\cdot\mathbf{v}_{rel}\right) \mathbf{n}\) : effective gravity due to velocity gradients
- \( \mathbf{v}_{rel} = \mathbf{v}_{drop} - \mathbf{v}_{p}\) :
- \( \mathbf{v}_{p}\): Velocity of boundary movement (see %ind_v_p(1)%)
Acceleration due to differences between current wetting angle and defined wetting angle in balance
\begin{align} \dot{\mathbf{v}}_{wetting} = \frac{\sigma}{\rho_{drop}} \frac{\alpha_{balance} - \alpha}{\alpha} \frac{\nabla H}{h^2}\end{align}
where
-
- \( \alpha_{balance}\) : defined wetting angle in balance
- \( \alpha\) : current wetting angle
- \( h\) : smoothing length
The gradient of the film thickness is \begin{align} \nabla H_{film,i} = \sum \limits_j \nabla W_{ij} V_j n_j = \sum \limits_j W_{ij} \left( -\frac{2\alpha}{h^2} \right) \mathbf{x}_{ij} V_j n_j\end{align} and stored in
%ind_grad_hwf(1)% ,
%ind_grad_hwf(2)% ,
%ind_grad_hwf(3)% .
\( \nabla H\) is applied here to compare the current wetting angle \( \alpha\) against the defined wetting angle in balance \( \alpha_{balance}\) (
BC_WettingAngle ).
Acceleration defined by the user
Acceleration defined by the user. This can be simply the gravity but also include acceleration with other chambers (here air) with media like air or others. Then this could include shear forces and pressure gradient from this chamber.
If in interaction with airflow (or another medium), the effective acceleration has to be set by the user as
\begin{align} \hat{\mathbf{g}} = \frac{1}{\rho_{drop}} \frac{\tau_{W,air}}{H_{film}} - \frac{1}{\rho_{drop}} \nabla p_{air} + \mathbf{g}\end{align}
where
-
- \( \tau_{W,air}\) : is the wall shear stress of the coupled phase (e.g. LIQUID chamber) typically calculated via the approxY() functionality
- \( \nabla p_{air}\) : is the pressure gradient of the coupled phase (e.g. LIQUID chamber) typically calculated via the approxY() functionality
- \( \mathbf{g}\) : is the gravitiy defined by gravity
Remarks
-
- This approach does not solve the partial differential equations commonly referred to as shallow water equations.
- Typically, for thin film modeling, points are not forced to keep a certain distance by means of an interaction potential. Although they might therefore be very close to each other, the height computation ensures that all of their volume is still represented.
- When a droplet in FreeFlight hits a boundary element with %IDENT_slip%, it forms a thin film.
- As the points in this setting represent mainly droplets, %TOUCH_liquid% is often a good choice.
- There is a specific viscosity definition in %MED_LIQUID_FILM%.