8.6. Contact line models
In the previous example, we had a static interface and free flow at the side walls "bottom"/"top" towards the liquid-gas interface. For e.g. an evaporating droplet, as already discussed in Section 7.6, the consideration of contact line dynamics is required. Contact lines can either be pinned or moving and of course this also can change during the drying of the droplet, causing a stick-slip behavior of the contact line. Moreover, when the surface tensions of the three interfaces in contact is changing at the contact line, usually also the equilibrium contact angle in case of a moving contact angle will change.
For multi-component flow, pyoomph has a single equation class, the DynamicContactLineEquations, which can be added to the contact line domain to handle all considerable cases. The particular dynamics is implemented in versatile contact angle models, which are passed to the DynamicContactLineEquations.
Both the versatile contact line models and the DynamicContactLineEquations are defined in the file pyoomph.equations.contact_angle.
8.6.1. PinnedContactLine
As we have seen in Section 7.6, a pinned contact line can be realized by enforcing partial_t(var("mesh"))\(=0\) on the test space \(\vec{v}\) of the velocity. Hence, if an instance of the PinnedContactLine model is passed to the DynamicContactLineEquations, it will add the following weak contribution
with the position-enforcing Lagrange multiplier \(\lambda\) with test function \(\mu\). \(\vec{X}\) is the position of the contact line and \(\vec{t}_\text{w}\) is here the tangent of the wall (i.e. the substrate in the case of an evaporating droplet). The wall tangent can be set by the wall_tangent keyword argument of the DynamicContactLineEquations class. It should be tangential to the wall, having a magnitude of unity and pointing inward (i.e. towards the droplet domain in the evaporating droplet case). Since we must have the possibility to add a contribution to the velocity test space in direction of the wall tangent, it is necessary to allow for slip by using e.g. a slip length boundary condition on the substrate.
8.6.2. UnpinnedContactLine
When the contact angle freely moves, we have seen how to impose an equilibrium contact angle \(\theta_\text{eq}\) in Section 6.6. More general, when the contact angle is measured with respect to the wall tangent \(\vec{t}_\text{w}\), we add
where \(\vec{n}_\text{w}\) is the wall normal, pointing into the droplet domain. With these vectors, the contact angle \(\theta\) approaches zero when the droplet is very flat, \(\theta=90\:\mathrm{^\circ}\) holds for a hemi-spherical droplet and for droplets on hydrophobic substrates, \(\theta>90\:\mathrm{^\circ}\) can be observed.
However, as also discussed in Section 6.6, the speed of the contact line motion by what \(\theta\) approaches the equilibrium \(\theta_\text{eq}\) is just controlled by the slip length. This is hard to control and to adjust with e.g. experimental observations. Instead, the UnpinnedContactLine let you control the motion velocity by a relation
with some typical velocity scale \(U_\text{cl}^0\). In the spirit of the PinnedContactLine, which is just given by \(U_\text{cl}^0=0\), the UnpinnedContactLine class therefore adds the weak contribution
The wall normal \(\vec{n}_\text{w}\) can be passed by the wall_normal argument to the DynamicContactLineEquations class, whereas the speed scale \(U_\text{cl}^0\) and the equilibrium contact angle can be passed via cl_speed_scale (defaults to \(10^{-5}\:\mathrm{m}/\mathrm{s}\)) and theta_eq (defaults to the initial contact angle) of the UnpinnedContactLine. If you want to modify the relation (8.7), you can do so by inheriting a custom contact line model from the UnpinnedContactLine and override the method get_unpinned_motion_velocity_expression() according to your demands.
To adjust the contact angle at infinite speed, i.e. instantaneously to the equilibrium contact angle, you can set cl_speed_scale=None. Thereby, the droplet will always be at \(\theta=\theta_\text{eq}\). This can be used e.g. to impose exactly some experimental dynamics, if you set theta_eq to a function of time that gives the experimental data.
8.6.3. StickSlipContactLine
We have seen similarities between the pinned and the unpinned contact line, where the former is just the latter, but with \(U_{\text{cl}}=0\). The additional Neumann term (8.6) imposing the equilibrium contact angle in (8.8) will be compensated anyhow by the Lagrange multiplier \(\lambda\) in the PinnedContactLine. Therefore, a stick slip motion can be realized by toggling between both contact angle modes with an indicator function \(\psi\), which is \(0\) if the contact line is pinned and \(1\) if the contact line moves according to (8.7). Hence, the relaxation velocity is just augmented by the factor \(\psi\) in the weak form:
The factor \(\psi\) can be controlled by the pin() and unpin() methods of the StickSlipContactLine. Initially, it is \(1\), i.e. the contact line could freely move, unless you call pin() before running the simulation. The speed of the contact line and the equilibrium contact angle can be set analogous to the UnpinnedContactLine.
However, it is cumbersome to switch the contact line dynamics by hand during the simulation. In reality, the contact line of an evaporating droplet starts usually pinned and switches to an unpinned motion once the contact angle \(\theta\) falls below some receding unpinning contact angle \(\theta^\text{rec}_\text{unpin}\). It can then re-pin, when the contact line recedes so that the contact angle has risen above some value \(\theta^\text{rec}_\text{pin}\), where \(\theta^\text{rec}_\text{unpin}<\theta^\text{rec}_\text{pin}<\theta_\text{eq}\) must hold for an alternating stick-slip motion during evaporation. If the droplet growth, e.g. due to condensation, we can have similar effects. Here, the contact line will unpin once the contact angle \(\theta\) raises above \(\theta^\text{adv}_\text{unpin}\) and pins again once the contact angle has fallen again below \(\theta^\text{adv}_\text{pin}\) due to the outward motion of the contact line. For stick-slip behavior, therefore \(\theta_\text{eq}<\theta^\text{adv}_\text{pin}<\theta^\text{adv}_\text{unpin}\) must hold. To set these contact angles, you can use the following methods:
Desired event |
relevant quantity |
method to use |
|---|---|---|
unpin if below and contact angle decreases |
\(\theta^\text{rec}_\text{unpin}\) |
|
pin if above and contact angle increases |
\(\theta^\text{rec}_\text{pin}\) |
|
pin if above and contact angle decreases |
\(\theta^\text{adv}_\text{pin}\) |
|
unpin if above and contact angle increases |
\(\theta^\text{adv}_\text{unpin}\) |
The first argument angle passed to these methods is the contact angle threshold, where this transition should happen. Note that the angles must be in increasing order following the table to to bottom and the equilibrium contact angle must be in between the receding angles and the advancing angles. If one passes by_factor=True, the first argument angle is not interpreted as angle, but as factor and the resulting threshold angle is this factor times the equilibrium angle. This means, that e.g. set_receding_unpin_below_angle(0.9,by_factor=True) sets \(\theta^\text{rec}_\text{unpin}=0.9\theta_\text{eq}\). Pass angle=None to remove a previously set angle. Further arguments are only_if_decaying/only_if_growing (default True) to indeed check whether the actual contact angle grows or shrinks when crossing the theshold. explicit (defaults to True) tells pyoomph to evaluate the actual contact angle \(\theta\) from the previous time step. If set to False, the contact line dynamics are fully implicitly considered. However, since the transitions from pinned and unpinned are discontinuous, it is likely not converging well in the Newton solver. One can improve it a bit by smearing out the transition angles with the heaviside_smoothing argument. Thereby, the transitions, which are implemented by heaviside step functions \(\Theta\), e.g. \(\Theta(\theta^\text{rec}_\text{unpin}-\theta)\), will be smoothed out by \(\operatorname{atan}((\theta^\text{rec}_\text{unpin}-\theta)/S)/\pi+1/2\) with \(S\) given by heaviside_smoothing. If explicit=True, heaviside_smoothing does not improve the convergence.
Once such contact angle threshold are set, the methods pin() and unpin() will not fix or freely move the contact line, unless it is allowed by the contact angle thresholds (i.e. in the hysteretic regions, e.g. if \(\theta\) is between \(\theta^\text{rec}_\text{unpin}\) and \(\theta^\text{rec}_\text{pin}\)). If we are outside these ranges, unpin() and pin() will have no effect. If we want to override this, i.e. force the contact line to be either pinned or unpinned, irrespectively of the set threshold angles, we must pass the forced=True to pin() or unpin(). The contact line will then remain is this state until one calls pin() or unpin() again. If this time the forced argument is False or omitted, the contact line dynamics based on the set angles will take over again.
8.6.4. YoungDupreContactLine
So far, we had to set the equilibrium contact angle \(\theta_\text{eq}\) by hand or it will default to the initial contact angle. However, when the surface tension of the liquid-gas interface changes, e.g. due to preferential evaporation of a multi-component droplet, the equilibrium contact angle will change as well. According to the Young-Dupré equation, the equilibrium contact angle is given by
where \(\sigma\), \(\sigma_\text{sg}\) and \(\sigma_\text{sl}\) are the surface tensions of the liquid-gas, solid-gas and liquid-solid interfaces at the contact line. The YoungDupreContactLine inherits from the StickSlipContactLine, i.e. all the stick-slip dynamics works as well. In the constructor, one can pass (besides the cl_speed_scale) the surface tensions sigma_sg and sigma_sl. If one do not pass these, its difference in the numerator of (8.10) will be determined by the initial contact angle and initial surface tension, so that the contact angle is initially at equilibrium.
8.6.5. KwokNeumannContactLine
The Kwok-Neumann contact line model was developed semi-empirical by measuring the contact angle of versatile liquid-substrate combinations [33]. The equilibrium contact angle reads
where the coefficient \(\beta\) (set via beta in the constructor) defaults to the value of Kwok and Neumann, \(\beta=124.7 \:\mathrm{m^4}/\mathrm{J^2}\). The solid-gas surface tension \(\sigma_\text{sg}^0\) can be set by sigma_sg_0 via the constructor. If not set, it will be calculated so that the initial contact angle is in equilibrium given the initial surface tension. Since the KwokNeumannContactLine inherits from the StickSlipContactLine, also the cl_speed_scale can be set and the stick-slip methods are available.