/** ![Water droplets may rest on a hydrophobic surface. Image via [the university of Wisconsin](https://www.chem.wisc.edu/content/hydrophobic-interactions).](https://www.chem.wisc.edu/deptfiles/Water.png) # Equilibrium shapes of a droplet on a hydrophobic surface An atmosphere with a density $\rho_2$ hosts a denser, watery droplet with a density $\rho_1 = \rho_2 + \Delta\rho$. The droplet is affected by the acceleration due to gravity ($g$) and rests on a perfect hydrophobic surface. A characteristic lengthscale of the droplet ($R$) may be introduced by considering it's volume; $V=4/3 \pi R^3$. We assume that there exists a single equilibrium state of the system that is governed by the so-called Bond number: $$Bo = \frac{\Delta \rho g R^2 }{\sigma}$$ On this page we aim to find the shape of the droplet for a range of $Bo$ values. Because of symmetry reasons, we choose to use the axi-symmetric solver to save some computational effort. */ #include "navier-stokes/centered.h" #include "two-phase.h" #include "tension.h" #include "axi.h" #include "reduced.h" #include "view.h" int maxlevel = 8; double R = 1.; f[left] = 0.; u.t[left] = dirichlet (0.); int sw = 1; char fname[99]; int it = 0; double AT = 0.5*10E-4; /** ## set-up We choose to use normalized values for all parameters that define the Bond number except for the value of the surface tension ($\sigma$). This parameter is used to control the aforementioned dimensionless ratio $Bo$. */ int main(){ G.x = -1; mu1 = 1.; mu2 = 1.; rho1 = 2.; rho2 = 1.; L0 = 4*R; init_grid (64); run(); } /** For $Bo \rightarrow 0$, a spherical droplet shape is expected. We therefore start with an small value for the Bond number and initialize the droplet as a sphere. */ event init (t=0){ TOLERANCE = 10E-4; f.sigma = 32.; refine (sq(x-R) + sq(y) < 1.1*R && level < maxlevel); fraction (f, R - sq(x-R) - sq(y)); } /** ## Finding equilibrium solutions Due to the effect of gravity, the droplet will deform. Deformation is characterized by velocities and we assume that equilibrium solutions have no kinetic energy in the system. We apply an arbritrary threshold (AT) that sets a lower limit for the kinetic energy. When the kinetic energy is lower than this value, the solution is considered to be in equilibrium. At this point we output the shape of the bubble to a file. Next, the value of the Bond number is increased by a factor of two and a new equilibrium is calculated, starting from the previous solution. Since this is not really a physically consistent set-up we also monitor the volume of the droplet to diagnose possible changes in the $R$ parameter. */ double e = 0.; event check_equi (i+=25;t<100){ double en=0.; foreach(reduction(+:en)){ foreach_dimension() en+=sq(u.x[])*sq(Delta); } if (sw == 0){ it++; if (it > 10) sw=1; } if (en