Dynamic contact angle does not limit to advancing or receding contact angles at high speed
Summary
When using dynamicAlphaContactAngle, I enter 3 values - an equalibrium CA, plus advancing and receding CA. From everything I've read on the topic, at high speeds the contact angle should limit to the advancing/receding angles, while at rest it will relax to the equilibrium CA.
Steps to reproduce
Use a VOF solver (e.g. interFoam) and set a boundary to dynamicAlphaContactAngle. Create a geometry/flow condition that results in a moving contact line along this boundary.
What is the current bug behaviour?
Contact angle limits to values that are not the advancing or receding contact angles specified in the boundary condition.
What is the expected correct behavior?
Contact angle should limit to advancing and receding contact angles specified in the boundary condition.
Relevant logs and/or images
This plot shows a simple case of contact angle hysteresis - equalibrium of 90deg and advancing and receding of 100 and 80 respectively. The current calculation limits to 110 and 70, while my code (below) correctly limits to the values specified.
In this case the angles are set asymetrically; my code again limits to the specified values, while the existing code does not. (It does not produce asymmetry at all, as it only really has two parameters: theta0 and (thetaA - thetaR)).
Environment information
- OpenFOAM version : v1912
- Operating system : Ubuntu 18.04
Possible fixes
Line 145 of dynamicAlphaContactAngleFvPatchScalarField is where the calculation is being done. I suggest changing it to:
return theta0_ * (1 - tanh(uwall/uTheta_) * sign(uwall/uTheta_) ) - thetaA_ * neg(uwall) * tanh(uwall/uTheta_) + thetaR_ * pos(uwall) * tanh(uwall/uTheta_);
It is probably worthwhile to do some sanity checking on thetaA>=theta0>=thetaR as if this is violated then it will probably give non-physical results.