# Droplet dynamics in capillary parallel plate channel

Damien Thomas, Stéphanie Lacour, Stéphane Zaleski

► To cite this version:

Damien Thomas, Stéphanie Lacour, Stéphane Zaleski. Droplet dynamics in capillary parallel plate channel. Computers and Fluids, 2025, 300, pp.106728. 10.1016/j.compfluid.2025.106728. hal-05154291

HAL Id: hal-05154291

https://hal.inrae.fr/hal-05154291v1

Submitted on 9 Jul 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Droplet dynamics in capillary parallel plate channel

Damien Thomas<sup>a,b</sup>,*, Stéphanie Lacour<sup>b</sup>, Stéphane Zaleski<sup>a</sup>

## ARTICLE INFO

**Keywords:**

Breakup

Contact lines

Droplet

Clogging

VOF

Spray cooling

## ABSTRACT

Sprays are more effective in dissipating heat compared to forced air convection. However, in heat exchangers, water films can block airflow through the fins during spray injection. This obstruction leads to a reduction in heat transfer and an increase in pressure drop. To address this issue, enhanced fins are essential to maintain a high heat transfer rate during spray cooling. Consequently, the study investigates the clogging processes related to the water films generated by sprays.

Numerical methods using the Volume Of Fluid (VOF) approach were employed to model the water-air interface of a droplet crossing a heat exchanger. The heat exchanger is modeled as a channel with two parallel plates as an embedded boundary. A variable-size droplet, representing an aggregate of spray droplets at the channel entry, is introduced and slides down the plates.

The outcome of the droplet penetration depends on factors such as the plate gap width, the droplet size, and the contact angle. In this context, overcoming an energy barrier is crucial for droplet penetration. The energy required, influenced by capillarity, makes hydrophobic surfaces challenging to penetrate. However, hydrophilic surfaces complicate droplet exit. Large droplets tend to break down into smaller ones during penetration. The crossing time, of similar magnitude for both hydrophobic and hydrophilic surfaces, depends mainly on the aspect ratio between the droplet diameter and the gap width. The drop break-up and the crossing time are thoroughly analyzed to identify the delicate balance of parameters essential for preventing channel clogging.

## 1. Introduction

Effective heat transfer is critical in thermal processes such as refrigeration, where improving cooling rates is essential. Spray cooling, a well-established method, utilizes a cloud of water droplets that impact the heat exchanger, boosting heat transfer through evaporation [1]. In the context of refrigeration, spray cooling has been extensively studied [2–13]. These studies have shown that the coefficient of performance can be enhanced by up to 60%, depending on air flow velocity and spray characteristics. When the heat exchanger temperature is below 373 K, the rate of evaporation is regulated mainly by air humidity. Air flow plays a crucial role by removing the forming water vapor, thus promoting continuous evaporation. Therefore, ensuring the unobstructed passage of air is crucial. The air obstruction, caused by droplets that block the heat exchanger, results in a decrease in heat exchange efficiency and an increase in pressure drop for air flow [14–22]. Vende's experiments showed that a surface clogging of 15% results in a 45% increase in head loss and a heat transfer reduction of 20% to 25% [23]. Consequently, understanding clogging is crucial for maintaining high heat transfer rates during spray cooling. This

research, therefore, investigates the processes that lead to the clogging of the heat exchanger surface by water films.

Clogging is largely governed by the interaction between the liquid and the surface, especially the contact angle, which determines how well a droplet spreads on or sticks to the surface. It is also influenced by several factors, including gravity, droplet velocity, and surface texture. However, the primary factor behind clogging is the wetting behavior, the object of numerous studies since the seminal work of de Gennes [24]. Clogging phenomena are also related to the study of droplet impacts, with comprehensive reviews provided by Yarin [25] and Josserand [26]. Papers on droplet impact on various textured, granular surfaces or pores [27–31] have reported that droplets simultaneously spread across the surface and penetrate into the pores.

Hydrophobic surfaces were thought to always prevent liquids from entering cylindrical capillaries and microfluidic channels. Marmur [32] developed a theory based on the balance of the Laplace pressure, which depends on the surface tension $\gamma$ and the droplet curvature, with the hydrostatic pressure. The contact angle $\theta$ affects this balance. His theory showed that a droplet fully penetrates a hydrophilic capillary,

<footnote><sup>a</sup> Sorbonne University and CNRS, Institut Jean Le Rond d'Alembert, France</footnote>

<footnote><sup>b</sup> University Paris-Saclay, INRAE, FRISE, 92761 Antony, France</footnote>

<footnote>* Corresponding author at: Sorbonne University and CNRS, Institut Jean Le Rond d'Alembert, France.</footnote>

<footnote>E-mail address: damien.mg.thomas@gmail.com (D. Thomas).</footnote>
## Nomenclature

### Interface & boundary quantities

$\delta_i$ Dirac function at the fluid interface

$\mathbf{n}$ Boundary normal vector

$\mathbf{n}_i$ Fluid interface normal vector

$\mathbf{t}$ Boundary tangential vector

### Characteristic numbers

$\text{Bo}$ Bond number

$\text{Mo}$ Morton number

$t_\gamma$ Characteristic time based on capillarity

$t_\mu$ Characteristic time based on viscosity

$t_c$ Crossing time

$t_g$ Characteristic time based on gravity

### Dynamic variables

$\alpha$ Angle between the tangent to the wall and the vertical direction

$\hat{\theta}$ Angle between the fluid interface and the vertical direction

$\kappa$ Droplet curvature

$\tau_{xy}$ Shear stress

$\mathbf{D}$ Deformation tensor

$\mathbf{f}_\gamma$ Surface tension force

$\mathbf{u}$ Velocity field

$\mathbf{x}$ Position vector

$c$ Fluid volume fraction

$c_s$ Solid volume fraction

$p$ Pressure field

$t$ Time

$V_f$ Final droplet volume

$V_i$ Initial droplet volume

$y_c$ Droplet mass center along y-axis

### Geometric parameters

$\Delta x$ Cell size

$\mathbf{D}$ Droplet diameter

$\mathbf{H}$ Channel height

$\mathbf{L}$ Channel gap width

$r$ Channel corner radius

$\mathbf{W}$ Channel width

### Indices & exponents

$b$ Back

$f$ Front

$g$ Gas

$l$ Liquid

$\sim$ Ejected mass incorporated

### Physical & numerical properties

$\gamma$ Surface tension

$\mu$ Dynamic viscosity

$\rho$ Density

$\theta$ Contact angle

$g$ Gravity

as expected. However, Marmur also demonstrated that a droplet with a contact angle up to 114° can fully penetrate a capillary when gravity is absent. For higher contact angles, only a part of the droplet penetrates

the capillary. In fact, penetration occurs when the droplet pressure exceeds the meniscus pressure plus the gravitational head.

Building on these insights, several experimental studies and simulations have further expanded the understanding of droplet behavior. Experimental research [33–37] has examined how droplets penetrate capillaries and porous media, while simulations using molecular dynamics [38–41], phase-field [42,43], and Volume of Fluid [44] (VOF) methods have provided deeper insights. These results have helped explain how surface tension and droplet interactions with porous surfaces influence pore penetration [45–49].

Denesuk [47] examined droplet behavior on an array of glass capillaries and found that droplets depleted faster when their contact area remained constant due to a pinned contact line, rather than when the contact line moved while maintaining a constant contact angle. Ding and Theofanous [50] observed that air bubbles can get trapped during droplet spreading and receding. Furthermore, inertia has been shown to drive droplets into capillaries, regardless of the capillary’s hydrophobic or hydrophilic nature [34].

Bouchard further investigated droplet impact between two parallel plates, revealing that at high velocities, droplet cleavage reduces the amount of liquid entering the channel [31]. This cleavage was also observed in cylindrical tubes [34,50] and orifice plates [51], indicating a common effect in different geometries. Bordoloi [51] found that higher Bond numbers lead to more fragmentation, especially on sharp-edged surfaces.

Marmur’s research on radial capillaries [52] and the work of Piroird et al. on capillary tubes [53] contributed to the development of an analytical model to extract droplets from capillaries using Laplace pressure. Researchers, such as Schebarchov and Hendy [54,55] used a model based on surface energies to determine penetration into a cylindrical capillary. Bormashenko [56–58] also explored how surface energy impacts wetting transitions between Cassie–Baxter and Wenzel states. Recently, Hagg [59] produced a theory based on Helmholtz free energy, predicting the likelihood of penetration of a drop within the gap between two hydrophobic parallel plates. His theory sheds light on the energy barrier that a droplet needs to overcome in order to penetrate a channel successfully. Droplet sliding occurs only when full penetration leads to a lower energy state than remaining outside the channel.

The studies conducted so far have provided valuable insight into the impact process of droplets and their penetration into channels. However, a notable gap exists in the literature with respect to the departure of droplets from channels. Our experimental observations indicate that clogging in the heat exchangers is not solely linked to droplet penetration at the channel entrance; it also occurs during the exit phase. Hagg’s analysis [59] is expanded to the entire channel in Fig. 1. At entry and exit, capillary forces at the contact lines in the hydrophilic cases point in the same directions. In contrast, the capillary forces in the hydrophobic cases have opposing directions. The directional difference between the top and bottom contact lines of the droplet arises from the vertical and horizontal orientations of the channel walls. Therefore, capillary forces can act as driving forces (Figs. 1(a) and 1(d)) or resistive forces (Figs. 1(b) and 1(c)). This suggests that clogging may occur at the bottom of a hydrophilic channel and at the top of a hydrophobic one.

Therefore, this study investigates droplet penetration and exit dynamics in hydrophilic and hydrophobic channels through 2D simulations between parallel plates. While the 2D simulations are well-suited to describe the geometry of a heat exchanger fin, they do not fully capture the detailed droplet behavior but are useful for understanding the dynamics influenced by capillarity. A single droplet — representing an aggregation of spray droplets — descends into the channel with a static contact angle imposed. The simulations are verified for convergence and the results are compared with Hagg’s prior observations. The focus is on characterizing the crossing time (the duration for a droplet to traverse the channel), establishing criteria that distinguish free-flowing from clogged states, and examining how droplet fragmentation influences clogging.
<figure><img src="image_3.png" /><figcaption>Fig. 1. Possible clogging scenarios for a droplet can occur either at the entrance or exit of a channel. The penetration process of a droplet in a hydrophilic channel is illustrated in (a), (b) depicts a droplet exiting a hydrophilic channel. On the other hand, (c) shows a droplet penetrating a hydrophobic channel, and (d) illustrates the exit of a droplet from a hydrophobic channel.</figcaption></figure>

## 2. Numerical method

### 2.1. Governing equations and flow modeling

A two-phase flow simulation is performed to model a water droplet moving within an airflow by solving the Navier-Stokes equations. Our approach utilizes a VOF solver implemented on quadtree meshes, employing the free software platform Basilisk. Detailed discussions on the Navier-Stokes solver can be found in [60–62]. This section outlines the structure and numerical methods applied.

In the context of multiphase flow, the fluid is divided into two phases: phase l (liquid, representing water) and phase g (gas, representing air). The phases are distinguished by a color function, denoted $c(x, t)$. The VOF method treats $c$ as a volume fraction, with values ranging from 0 (pure phase g) to 1 (pure phase l), and intermediate values indicating interface cells. The advection equation governs the evolution of this color function over time. The fluid volume fraction $c$ is passively advected by the velocity field $u(x, t)$, ensuring the tracking of the interface dynamics. Its evolution is governed by the advection equation

$$ \frac{\partial c}{\partial t} + \mathbf{u} \cdot \nabla c = 0. \qquad (1) $$

For incompressible and isothermal flow, the conservation equations for mass and momentum are expressed as

$$ \nabla \cdot \mathbf{u} = 0, \qquad (2) $$

$$ \frac{\partial \rho \mathbf{u}}{\partial t} + \nabla \cdot (\rho \mathbf{u} \otimes \mathbf{u}) = -\nabla p + \nabla \cdot (2\mu \mathbf{D}) + \rho \mathbf{g} + \mathbf{f}_{\gamma}, \qquad (3) $$

where $\rho(x, t)$ and $\mathbf{D} = \frac{1}{2}(\nabla \mathbf{u} + \nabla \mathbf{u}^T)$ are, respectively, the pressure field and the deformation tensor with $x$ the position and $t$ the time. The properties of the flow are denoted $\rho$ for the density and $\mu$ for the dynamic viscosity. In this formulation, a single set of velocity and pressure fields is defined throughout the domain, while the fluid properties $\rho$ and $\mu$ vary according to the local phase composition as determined by the volume fraction $c$. The vector $\rho \mathbf{g}$ is the gravitational force. The last term $\mathbf{f}_{\gamma}$ represents the surface tension force

$$ \mathbf{f}_{\gamma} = \gamma \kappa \mathbf{n}_i \delta_i, \qquad (4) $$

where $\gamma$ is the surface tension. The force $\mathbf{f}_{\gamma}$ operates specifically on the free surface, which is represented by the Dirac function $\delta_i$ supported at the fluid interface. This force depends on both the curvature $\kappa$ and the unit normal vector $\mathbf{n}_i$ at the fluid interface.

Momentum advection is performed by a combination of the Bell–Collela–Glaz scheme typical of Basilisk and a momentum-consistent scheme [63], see also [64]. The density averaging is tuned to the momentum-consistent method. Viscosity averaging in mixed cells is performed by the harmonic mean. As explained in [65], the harmonic mean is more accurate for flows parallel to the interface, but leads

to less stable computations than the arithmetic mean. For that reason the arithmetic mean used to be more frequently used. However, with the expanded use of momentum-consistent schemes, stability issues for water–air systems at moderate Reynolds and Weber numbers such as the current one are not present, which allows the use of the more accurate harmonic mean.

### 2.2. Embedded boundary and contact line

The solid channel is modeled using an embedded boundary approach, as described in [66,67]. A volume fraction field $c_s$ represents the solid phase, and the contact angle $\theta$ is imposed at the contact line. To apply the contact angle at the contact line, the interface is extended into the embedded boundary $c_s$ and adjusted to match the specified contact angle. This adjustment involves two main steps [68]:

1. A piecewise linear reconstruction is assumed for both the fluid–fluid and solid–fluids interfaces to determine the geometry near the contact line.

2. The geometric information is used to extend the volume fractions and the fluid interface into the solid phase. This in turn allows to advect the volume fraction despite the presence of the solid.

In the current work this contact line treatment is implemented within the existing VOF framework of Basilisk with no modifications or new developments.

This method was tested with a sessile drop on flat and inclined surfaces. For the flat surface, the relative error in the contact angle was within 3% for angles ranging from 15° to 165° [69]. However, on an inclined surface, the accuracy varied. The relative error increased to about 6% for contact angles between 40° and 150° [70].

For dynamic contact line modeling, the embedded boundary approach currently lacks a well-tested method in the quadtree VOF framework [71]. Although the method employed has been validated for static cases, its performance in dynamic situations remains unverified in the literature. Clogging depends primarily on energy barriers, which are dictated by the static contact angle rather than by dynamic contact line motion. The absence of a dynamic contact line model is not expected to significantly affect our results. In fact, contact line motion on embedded boundaries has been known to cause numerical instabilities or spurious currents, particularly when solid surfaces are misaligned with the grid. No such instabilities were observed, probably because the solid surfaces remain mostly aligned with the computational grid. Moreover, the convergence study (detailed in the following) confirms that dynamic contact line singularities do not influence convergence of the wetting dynamics. Therefore, any errors introduced by the dynamic contact line are expected to be minimal.

Finally, spatial discretization used a collocated quadtree grid, which is a tree grid data structure with AMR (Adaptive Mesh Refinement). The computational efficiency is improved by employing AMR, which
enables the mesh resolution to dynamically adjust based on estimated errors in the color function of fluid c, volume fraction c_s and the velocity field u. Furthermore, the Navier-Stokes equations are solved using a fully implicit scheme for the viscous term and a second-order projection method, the Bell-Collela-Glaz advection scheme. To handle sharp interfaces, the Piecewise Linear Interface Capturing Volume Of Fluid (PLIC-VOF) approach was used as described in [72,73], for interface reconstruction.

### 3. Simulation set-up

#### 3.1. Set-up

Simulations involve a water droplet of diameter D within air, the droplet center is initially located at (0, D/2). The droplet descends without an initial velocity between two parallel plates separated by a length L and a height H with a corner radius r. An additional constraint is imposed to ensure that the channel width W is sufficiently large, preventing the droplet from escaping laterally. Consequently, the plates extend from the top at coordinates (±L/2, 0) to (±max(2D, 4L), -H). The ratios L/H and r/L are maintained respectively at a constant value of 8.23 × 10⁻² and 0.1, values representative of a fin-and-tube heat exchanger, where the fin height is H = 15.8 mm and the gap width L = 1.3 mm. This setup approximates the pore geometry encountered in real condenser fins.

The boundary conditions are illustrated in Fig. 2, the left, right, and bottom boundaries are set as outflow conditions,

$$ \frac{\partial u_n}{\partial n} = 0, \quad p = 0 \quad \text{and} \quad \frac{\partial c}{\partial n} = 0, \qquad (5) $$

with n the normal vector to the domain boundary. The top boundary is set as a low-intensity inflow, with a speed equal to 0.1 times the free-fall velocity of the droplet falling from height H,

$$ u_n = -0.1 \sqrt{2 g H}, \quad \frac{\partial p}{\partial n} = 0 \quad \text{and} \quad c = 0. \qquad (6) $$

The embedded boundary applied on the solid, uses impermeable, no-slip conditions to stabilize the contact line relaxation,

$$ u_n = u_r = 0, \qquad (7) $$

with n and r respectively denoting the normal and tangential components to the interface.

The no-slip boundary condition causes a force and dissipation singularity at a moving contact line: viscous dissipation and stresses at the contact line become infinite [74]. In numerical simulations, methods such as VOF introduce a small numerical slip length [75] of the order Δx/2 that allows fluid motion, with Δx the cell size. In large viscosity fluids and well-resolved simulations, the effects of the singularity become more visible as the mesh is refined. Non-convergence of the solution has often been observed [75]. Nevertheless, for water and length scales of centimeters, the Laplace number (equivalent to the Reynolds number but based on the capillary wave speed) is of order 10⁶ and viscous effects including the singularity are small. Therefore, a degree of convergence is achieved [76].

#### 3.2. Scaling

In this study, the flow is driven by gravity, with capillary forces that may oppose it. A key parameter to characterize the competition between these forces is the Bond number Bo defined as the ratio of gravitational to capillary forces. Several geometric parameters influence the flow outcome, such as the gap width L, the droplet diameter D, the channel height H or the channel curvature radius r. To better understand the role of these parameters and to identify which ones should be used to define the Bond number, a simplified model is considered. Specifically, the case of a droplet penetrating a hydrophobic channel is examined (see Fig. 1(c)), with a small portion already inside.

<figure><img src="image_4.png" /><figcaption>Fig. 2. Visualization of the boundary conditions and initial position of the droplet. The graph's origin is denoted by O, with its x- and y-axes.</figcaption></figure>

The droplet is assumed to be in static equilibrium, and gravitational effects on its curvature are considered negligible. In this configuration, two types of pressure act on the droplet, the hydrostatic pressure ρ_l g h where h ≈ D is the droplet height and the capillary pressure γ(κ_f + κ_b) where κ_b ≈ 2/D is the curvature at the back interface of the droplet and κ_f = 2 cos θ/L is the curvature at the front interface inside the channel. Setting the total pressure difference to zero for equilibrium gives the balance ρ_l g D + γ(2 cos θ/L + 2/D) = 0. This can be rearranged into the dimensionless form

$$ \mathrm{Bo} = \frac{\rho_l g D^2}{\gamma} = -2 \left( 1 + \frac{\cos \theta}{L/D} \right). \qquad (8) $$

Thus, the behavior of the system is governed by two key dimensionless parameters the Bond number Bo, using the droplet diameter as the characteristic length, and an aspect ratio L/D.

To estimate how long it takes for the droplet to cross the channel once motion begins, characteristic time scales must be introduced. These time scales capture the relative influence of gravity, capillarity, and viscosity on the dynamics of the droplet. A characteristic time scale is defined as t_g = √(2H/g), which represents the time required for a free falling object to cross a height H under gravity. The main force opposing motion is capillarity, with characteristic time scale t_γ = √(ρ_l L³/γ). It relates to the Bond number Bo, the aspect ratios L/D and L/H through a ratio between t_γ and t_g,

$$ \frac{t_{\gamma}}{t_{g}} = \sqrt{\frac{\mathrm{Bo} L^3}{2 \mathrm{D}^2 H}} \qquad (9) $$

In addition, the viscous force is related to a timescale t_μ = ρ_l L²/μ_l. The Morton number is

$$ \mathrm{Mo} = \frac{t_{\gamma}^{6}}{t_{\mu}^{4} t_{g}^{2}} \frac{2H}{L} = \frac{g \mu_{l}^{4}}{\rho_{l} \gamma^{3}}. \qquad (10) $$

The Morton number of water at 293 K is Mo = 2.57×10⁻¹¹, meaning that the droplet oscillations and shape changes are dominated by capillarity rather than viscosity. The density ratio is ρ_l/ρ_g = 829 and the dynamic viscosity ratio is μ_l/μ_g = 55.4. The latter three numbers are kept constant at these values throughout.

Several geometric parameters affect the outcome, including the length L, height H, radius r of the solid and the droplet diameter D. However, the focus is limited to L and D, analyzed through the Bond number Bo, the aspect ratio L/D. Moreover, the solid properties are varied with the contact angle θ. To encompass a wide range of
<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Convergence study for a contact angle $\theta = 65^{\circ}$, a Bond number of 1 and an aspect ratio $L/D = 0.5$. Plot of the dimensionless mass center $y_c/H$ as a function of the dimensionless time $t/t_g$ while varying the minimum number of cells inside the width gap $L$. The top and bottom channel surfaces are delineated by the red and black lines, respectively.</figcaption></figure>

scenarios, the Bond number is varied in the range of 0.1 to 10, around the capillary length. The droplet is always larger than the gap width, with $L/D$ ranging from 0.1 to 0.9. Finally, the contact angle $\theta$ is set at both $65^{\circ}$ and $150^{\circ}$, representing hydrophilic and hydrophobic conditions, respectively.

Simulations were also conducted with a drop center positioned above the left wall. The results with the shifted center are similar to those obtained with a central injection point. Since the initial drop position does not alter the conclusions, these findings are omitted from this article.

### 3.3. Convergence study

Given the importance of capillarity in our study and the high curvature of the interface inside the channel, ensuring a high level of refinement near the contact line is crucial. To ensure adequate capture of the contact line, even when the droplet and the channel are large, simulations feature varying numbers of cells. In Fig. 3, the convergence of the dimensionless center of mass of the droplet $y_c/H$ is illustrated. After an initial period of descent, the droplet rebounds when attempting to exit the channel, then oscillates. However, increasing the level of refinement reduced the oscillations caused by the contact line instability. A minimum of 40 cells is necessary to accurately predict the droplet position. Consequently, subsequent simulations utilize a cell count $L/\Delta x$ ranging between 40 to 80 cells in the gap. Simulations are carried out at a refinement level of 10 or 11, depending on the droplet size. This corresponds to a cell count equivalent to a domain size ranging from $2^{10} \times 2^{10}$ to $2^{11} \times 2^{11}$.

## 4. Results and discussions

### 4.1. Description of droplet behavior

Simulations have yielded four penetration scenarios as those depicted in Fig. 1. In hydrophilic cases, clogging occurs at the bottom of the channel, whereas in hydrophobic cases, it occurs at the top of the channel.

The droplet's progression in hydrophilic channels, as illustrated in Figs. 5(a) and 5(b), can be divided into the following three stages:

* **Entry:** Gravity initially pushes the droplet downwards, causing it to contact the channel corners. The capillary forces, directed downward, assist gravity in moving the droplet downwards. Meanwhile, large droplets (characterized by a large Bond number and a small aspect ratio of $L/D$) break. Residual fluid is retained on the upper surface of the channel while the other portion of the droplet descends.

* **Descent:** After the initial descent, the droplet becomes fully enclosed within the channel. Capillary forces acting along the vertical axis at the front and back interfaces are expressed as linear forces with their y-components given by $f_{y,y}^f = -f_{y,y}^b = -\gamma \cos \theta$. Therefore, the projections of these capillary forces along the vertical axis counterbalance each other, leaving gravitational forces as the sole factor driving the droplet's descent.

* **Exit:** As the droplet approaches the bottom of the channel, the curvature of the corner causes the lower interface of the droplet to bend outward, effectively increasing the effective angle $\hat{\theta} = \theta + \alpha$, where $\alpha$ is defined as the angle between the tangent to the wall and the vertical direction. An angle of $\alpha = 0^{\circ}$ indicates that the contact line is located at the beginning of the corner, where the wall is vertical. In contrast, an angle of $\alpha = 90^{\circ}$ indicates that the contact line is located at the end of the corner, where the surface is horizontal (see Fig. 4 for a schematic). Consequently, the vertical component of the capillary forces exerted by the front interface is largest when $\hat{\theta}$ approaches $180^{\circ}$. This means that $\theta$ and $\alpha$ tend to $90^{\circ}$, therefore the contact line tends to approach the bottom flat edge. Depending on the aspect ratio and the Bond number, the resistance of capillarity is sufficient to counteract gravity, leading to a rebound and subsequent clogging. However, in other scenarios, the droplet continues its descent and eventually exits the channel. As the droplet detaches from the solid surface, a thin layer of residual fluid is left behind as a result of the separation process.

For hydrophobic surfaces (Figs. 5(c) and 5(d)), the droplet's evolution follows a similar progression, although it is influenced by opposing capillary effects:

* **Entry:** The droplet initially makes contact with the channel at two points. However, unlike hydrophilic surfaces, where the capillary forces at both interfaces point downwards, hydrophobic surfaces have a front interface that points upward and a back interface pointing downwards. For contact angles $\theta$ greater than $114^{\circ}$ [32], the combined capillary forces point upward, impeding the droplet's spontaneous penetration into the gap. At rounded corners, capillary forces at the front interface peak when $\hat{\theta}$ nears $180^{\circ}$, requiring an external force, such as gravity, for penetration.

* **Descent:** Once the droplet is fully inside the channel, the projections of the capillary forces oppose each other. As in the hydrophilic case, the descent is controlled by gravity. During exit, the capillary forces at the front interface reach their maximum at the rounded corner, whereas the forces at the back interface point downwards, easing the droplet's movement. This creates a brief upward capillary force that can momentarily slow the droplet down, but gravity and inertia continue to drive it downwards.

* **Exit:** Once the front interface exits the capillary, the capillary forces assist in the extraction process [53].

The four scenarios highlight the need to identify the parameter pairs $(L/D, Bo)$ for which clogging occurs. Additional questions arise about the time required for the droplet to pass through the channel.

### 4.2. Crossing time

As noted above, the time it takes for a droplet to pass through the channel is an important parameter. There is a risk that the channel is blocked by subsequent drops. Figs. 6(a) and 7(a) show the normalized
<figure><img src="image_6.png" /><figcaption>Fig. 4. Schematic of the front capillary force $t'_f$ at the front interface as the droplet progress along the bottom corner of the channel. The interface is shown at different positions along the corner arc, corresponding to angles $\alpha = \{0^\circ, 45^\circ, 90^\circ\}$. The illustration assumes a contact angle $\theta = 65^\circ$.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 5. Simulations of a water droplet descending through two-parallel plates for an aspect ratio $L/D = 0.5$. In cases (a) and (b), a contact angle of 65° is employed, while cases (c) and (d) involve a contact angle of 150°. The Bond number is specified as 0.278 for cases (a) and (c), and 2.154 for cases (b) and (d). The colors represent the magnitude of the velocity field.</figcaption></figure>

vertical position of the droplet center of mass over time for hydrophilic ($\theta = 65^\circ$) and hydrophobic ($\theta = 150^\circ$) cases. The crossing time for droplets successfully traversing the channel is depicted in Figs. 6(b) and 7(b). The colors represent various Bond numbers, and the marker shapes indicate the aspect ratio $L/D$.

### 4.2.1. Hydrophilic surface

For $\theta = 65^\circ$, Fig. 6(a) shows that most droplets descend below the upper surface (red line) when entering the channel for $\theta = 65^\circ$, which aligns with theoretical expectations. Droplets with smaller Bond numbers and aspect ratios $L/D$ descend faster. In fact, a small Bond number implies large capillary forces, which, being oriented downwards, result in a rapid initial acceleration of the droplet. The droplets then separate into two groups:

* **Exiting droplets:** These droplets manage to exit the channel (center of mass below the black line), the crossing time averages around $1.3t_g$ and decreases as $L/D$ gets smaller. Droplets with small aspect ratios are flattened, assuming an elongated shape. This elongation increases the crossing time, as the droplet tail enters the channel much later than the droplet front.

* **Trapped droplets:** These droplets, with their center of mass above the black line, cannot exit the channel. When the droplets make contact with the bottom corners, the upward-directed capillary forces cause them to rebound.

For $Bo < 2.154$ and $L/D = 0.1$, the droplets appear to remain above the red line, suggesting that they cannot enter the channel. In fact, these droplets enter the channel. However, their volume is substantial enough to completely fill the channel, resulting in their center of mass being positioned above it. The position of these droplets is poorly represented by their center of mass.

### 4.2.2. Hydrophobic surface

For a contact angle of $\theta = 150^\circ$, the capillary forces exert the opposite effect compared to the hydrophilic scenario. Let us examine how these forces impact the crossing time. In Fig. 7(a), droplets separate into two groups:

* **Exiting droplets:** These droplets pass through the channel with ease. Analyzing the slope of the center of mass reveals a consistent final velocity. However, larger droplets (with a large Bond number and a small aspect ratio $L/D$) show an initial offset, especially for $L/D \le 0.5$. In fact, a contact angle $\theta$ greater than 90° implies that the surface energy of the solid-gas interface is lower than that of the solid-liquid interface, according to Young's equation. This causes the solid to preferentially interact with the gas rather than the liquid. As a result, the droplet center moves downwards in the channel, while its tail spreads, then recedes on the top surface before being pulled through the channel. This process leads to an initial offset. Fig. 7(b) shows that as the droplet is more elongated at smaller $L/D$, the crossing time increases.

* **Repelled droplets:** These droplets struggle to penetrate the channel, evidenced by their center of mass located above the red line, indicating that they remain above the channel. They exhibit oscillations of their center of mass position $y_c$, with amplitudes that increase with aspect ratio $L/D$, due to droplet spreading and receding dynamics.

When droplets exit the channel (Fig. 7(b)), the crossing time is mainly affected by the Bond number. Unlike in hydrophilic cases, there is no initial acceleration burst, making gravity's effect via the Bond number more apparent. Fig. 7(b) shows that the crossing time ranges from $1.1t_g$ to $3.2t_g$. Droplets with small aspect ratios take slightly longer to cross as they require more time to fully enter the channel.
<figure><img src="image_8.png" /><figcaption>**Fig. 6.** For a contact angle of 65°: (a) shows the droplet's center of mass along the y axis, relative to the channel height *H*, plotted against the dimensionless time *t*/*t*<sub>g</sub> for various Bond numbers Bo and L/D ratios. The superior and inferior channel surfaces are delineated by the red and black lines, respectively. (b) depicts the dimensionless crossing time *t*<sub>c</sub>/(*t*<sub>g</sub>) for varying L/D and Bo values. Colors and shapes represent, respectively, the Bond number and the aspect ratio.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. 7.** For a contact angle of 150°: (a) illustrates of the droplet's center of mass along the y axis, relative to the channel height *H*, plotted against the dimensionless time *t*/*t*<sub>g</sub> for various Bond numbers Bo and L/D ratios. The superior and inferior channel surfaces are delineated by the red and black lines, respectively. (b) depicts the dimensionless crossing time *t*<sub>c</sub>/(*t*<sub>g</sub>) for varying L/D and Bo values. Colors and shapes represent, respectively, the Bond number and the aspect ratio.</figcaption></figure>

### 4.3. Droplet breakup on entry

The mechanisms that lead to clogging were highlighted, emphasizing the need to quantify clogged states, while the role of droplet breakup still requires further investigation. As breakup reduces droplet volume, which reduces the hydrostatic pressure, potentially impacting the clogging scenarios. Figs. 8(a) and 8(b) show the ratio of the final volume to the initial volume $V_f/V_i$ for hydrophilic and hydrophobic surfaces. The initial volume (per unit length in 2D) is the unit disk area $V_i = \pi/4D^2$, the final volume is the volume of the main droplet in the final state and the main droplet is selected as the largest and most centered droplet. In most cases, less than 50% of the mass is lost, but for large droplets (large Bond number Bo ≥ 5.995 and small aspect ratio L/D ≤ 0.2), the ejected mass can reach 80%. Fragmentation occurs mainly during initial penetration, in the spreading phase. High contact angles reduce fragmentation by causing the droplet tail to recede. Consequently, 78% of the hydrophilic cases show fragmentation, compared to only 5% of the hydrophobic cases. As the Bond number increases,

the fragmentation also increases. The smallest volume is approximated by $V_i/V_f = 0.59L/D + 0.26$ for Bo = 10 and $\theta = 65^\circ$. This sets a limit for Bond numbers that ranges from 0.1 to 10. Fig. 9 shows the dimensionless shear stress $\tau_{xy}/(\rho_l g D)$ and the velocity field of a droplet undergoing fragmentation. Here, $\tau_{xy} = \mu_l(\partial_j u_x + \partial_x u_y)$, and $\rho_l g D$ represents the hydrostatic pressure of a water column of height $D$. As the droplet moves downwards, its sides extend outward, while its center continues to descend. The droplet separates at two stagnation points located at approximately $x = \pm(L+r)$, the junction of the horizontal surface and the rounded corner, where the streamlines separate and the shear stress is zero. We note that in hydrophilic cases, fragmentation significantly impacts the exit process.

### 4.4. Critical conditions for clogging

The clogged states, depicted as a function of L/D and Bo, are visualized in Figs. 10(a) and 10(b), with marker shapes and the black line delineating the free states from the clogged. Among the total
<figure><img src="image_10.png" /><figcaption>**Fig. 8.** Droplet volume ratio between the final $V_f$ and initial $V_i$ state for hydrophilic $\theta = 65^\circ$ (a) and hydrophobic $\theta = 150^\circ$ (b) cases, considering various values of Bo and $L/D$. The red dotted line indicates a linear regression of the curve for Bo = 10.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 9.** Dimensionless shear stress $\tau_{xy}/(\rho_l g D)$ and velocity field at $t = 0.28 t_g$ for parameters $L/D = 0.4$, Bo = 6 and $\theta = 65^\circ$. The droplet is in the process of fragmenting as it descends into the channel.</figcaption></figure>

90 simulations conducted for hydrophilic cases, 54 drops successfully crossed the channel, while 36 failed to do so. Similarly, for hydrophobic cases, 43 drops crossed, while 47 failed out of the 90 simulations. Droplets are more likely to pass through the channel when the Bond number and aspect ratio are large, as capillarity hinders this process. In fact, the Laplace pressure is given by $\gamma\kappa(1/L)$, linearly depending on $1/L$. Therefore, as $L/D$ increases, the capillary effects, which contribute to clogging, decrease. This trend is consistent in both hydrophilic and hydrophobic scenarios. In the hydrophilic scenario, there is a moderate increase in the occurrence of clogging when $L/D$ is small and Bo is large, while in the hydrophobic scenario, clogging occurs less frequently when $L/D$ is large and Bo is small. In particular, the channel is always clogged when $Bo < 0.6$ for hydrophilic surfaces.

The colors in Figs. 10(a) and 10(b) represent the percentage of ejected mass, seen primarily in larger droplets with large Bond numbers and aspect ratios. On the hydrophilic surface, when the Bond number is small ($Bo < 2.154$) and the aspect ratio is $L/D = 0.1$, droplet

fragmentation does not occur because the droplet fills the channel before its top interface reaches the stagnation points. However, for $L/D = 0.2$, the droplet is not large enough to fill the channel, leading to breakup even at small Bond numbers. As the aspect ratio increases, a greater portion of the droplet mass is initially located above the channel gap, resulting in a decrease in the ejected mass.

To further incorporate the impact of the ejected mass, an effective Bond number $\tilde{B}_o$ and an aspect ratio $L/\tilde{D}$ where $\tilde{D}$ is the droplet diameter derived from the remaining volume $\tilde{V}$ after fragmentation. The small graphs in Figs. 10(a) and 10(b) show that the boundary between the clogged and free states has a less steep slope when considering fragmentation, especially for hydrophilic surfaces. Since ejected mass occurs mainly for small $L/D$ and large Bo, clogging situations remain unchanged for large $L/D$ and small Bo. Thus, droplet fragmentation tends to favor clogged states by reducing the effective volume of the droplet within the channel, which diminishes gravitational forces. Therefore, models predicting droplet exit from hydrophilic channels should consider fragmentation for aspect ratios $L/D$ below 0.5.

Finally, our findings are consistent with those of Hagg [59]. His analysis and simulations show that the transition from free-flowing to clogged states occurs for an aspect ratio $L/D$ between 0.4 to 0.5 at a contact angle of 150° and $Bo = 1$.

## 5. Conclusion

This study was initiated to improve the efficiency of cooling heat exchangers using sprays. Although spray cooling improves heat transfer, the potential for channel clogging emphasizes the need for a thorough understanding of its dynamics. The investigation focused on the behavior of liquid droplets in capillary systems, particularly on hydrophilic and hydrophobic surfaces.

Numerical simulations were performed using an adaptive PLIC-VOF solver, with the channel modeled as an embedded boundary. Droplet dynamics were analyzed in channels formed by two parallel plates. The analysis considered key parameters, including the Bond number $Bo$, the ratio of channel length to droplet diameter $L/D$, and contact angle $\theta$, each assessed for its impact on droplet clogging and crossing time. The key findings are summarized as follows:

1. **Surface effect:** Distinct behaviors were observed on hydrophilic and hydrophobic surfaces. On hydrophilic surfaces, small droplets quickly entered the channel because of capillary forces,
<figure><img src="image_12.png" /><figcaption>**Fig. 10.** Clogging state as a function of L/D and Bo for both hydrophilic $\theta = 65^{\circ}$ (a) and hydrophobic $\theta = 150^{\circ}$ (b). Squares and triangles represent respectively clogged channels and free channels. The black line outlines the clogging boundary, and colors indicate the corresponding ejected mass due to fragmentation. The small graphs in the corners reflect the ejected mass values in $\bar{B}_o$ and $L/\bar{D}$.</figcaption></figure>

but were at risk of clogging when exiting. In contrast, hydrophobic surfaces made it harder for droplets to enter the channel while facilitating their exit. This means that capillarity acts as an energy barrier that affects how droplets move in and out of the channel. Only kinetic energy or external forces such as gravity can overcome this barrier.

2. **Entry and exit dynamics:** During the entry process, capillary forces help droplets accelerate into hydrophilic channels. However, these same forces make it difficult for the droplets to exit. Once a droplet is inside the channel, the capillary forces at both ends balance each other out. As a result, the downward movement of the droplets depends only on gravity. The time it takes to cross the channel is influenced by gravity, the contact angle, the ratios L/D and L/H, ranging from 1.1$t_g$ to 3.2$t_g$. Surface type also plays a role: hydrophobicity introduces an initial delay as the droplet tail spreads and then recedes.

3. **Fragmentation impact:** Larger droplets, characterized by large Bond numbers and aspect ratios, tend to fragment when entering the channel, particularly on hydrophilic surfaces. This fragmentation can lead to an ejected mass of up to 80%, significantly reducing the amount of liquid that enters the channel. This, in turn, greatly affects the exit process of droplets on hydrophilic surfaces. Therefore, models predicting exit dynamics on hydrophilic surfaces should take fragmentation into account.

Our study used simulations to explore droplets behavior in capillary channels, providing insights into the dynamic of clogging, penetration, and exit. The 2D approach was computationally efficient and captured essential trends, especially the effects of surface properties. However, 2D simulations have limitations compared to 3D. In 3D, phenomena such as droplet shape asymmetry, surface curvature, and more accurate fragmentation would be better captured. Although 2D results are valuable for delineating the phase diagram, future 3D simulations would offer a more complete and realistic picture of droplet dynamics in heat exchangers, especially for more complex cases.

### CRediT authorship contribution statement

**Damien Thomas:** Writing – review & editing, Writing – original draft, Visualization, Investigation, Formal analysis, Conceptualization. **Stéphanie Lacour:** Validation, Supervision, Investigation, Conceptualization. **Stéphane Zaleski:** Validation, Supervision, Formal analysis.

### Declaration of competing interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Damien THOMAS reports financial support was provided by European Research Council. If there are other authors, they declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

### Acknowledgments

This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No. 883849). We thank the European PRACE group, the Swiss supercomputing agency CSCS, the French national GENCI supercomputing agency and the relevant supercomputer centers for their grants of CPU time on massively parallel machines, and their teams for assistance and the use of Irene-Rome at TGCC and ADASTRA at CINES.

### Data availability

The source code [77] used to run the simulations can be found on Basilisk website.

### References

[1] Riofrío MC, Caney N, Gruss J-A. State of the art of efficient pumped two-phase flow cooling technologies. Appl Therm Eng 2016;104:333–43. http://dx.doi.org/10.1016/j.applthermaleng.2016.05.061.

[2] Hamlin S, Hunt R, Tassou SA. Enhancing the performance of evaporative spray cooling in air cycle refrigeration and air conditioning technology. Appl Therm Eng 1998;18(11):1139–48. http://dx.doi.org/10.1016/S1359-4311(98)00028-3.

[3] Vrachopoulos MG, Filios AE, Kotsiovelos GT, Kravvaritis ED. Incorporated evaporative condenser. Appl Therm Eng 2007;27(5):823–8. http://dx.doi.org/10.1016/j.applthermaleng.2006.09.021.

[4] Youbi-Idrissi M, Macchi-Tejeda H, Fournaison L, Guilpart J. Numerical model of sprayed air cooled condenser coupled to refrigerating system. Energy Convers Manage 2007;48(7):1943–51. http://dx.doi.org/10.1016/j.enconman.2007.01.025.

[5] Hajidavalloo E, Eghtedari H. Performance improvement of air-cooled refrigeration system by using evaporatively cooled air condenser. Int J Refrig 2010;33(5):982–8. http://dx.doi.org/10.1016/j.ijrefrig.2010.02.001.
[6] Yu FW, Chan KT. Improved energy performance of air-cooled chiller system with mist pre-cooling. Appl Therm Eng 2011;31(4):537-44. http://dx.doi.org/10.1016/j.applthermaleng.2010.10.012.

[7] Yang J, Chan KT, Wu X, Yang X, Zhang H. Performance enhancement of air-cooled chillers with water mist: Experimental and analytical investigation. Appl Therm Eng 2012;40:114-20. http://dx.doi.org/10.1016/j.applthermaleng.2012.02.001.

[8] Tissot J, Boulet P, Trinquet F, Fournaison L, Lejeune M, Liaudet F. Improved energy performance of a refrigerating machine using water spray upstream of the condenser. Int J Refrig 2014;38:93-105. http://dx.doi.org/10.1016/j.ijrefrig.2013.08.025.

[9] Zhang F, Bock J, Jacobi AM, Wu H. Simultaneous heat and mass transfer to air from a compact heat exchanger with water spray precoaling and surface deluge cooling. Appl Therm Eng 2014;63(2):528-40. http://dx.doi.org/10.1016/j.applthermaleng.2013.11.046.

[10] Wang T, Sheng C, Nnanna AGA. Experimental investigation of air conditioning system using evaporative cooling condenser. Energy Build 2014;81:435-43. http://dx.doi.org/10.1016/j.enbuild.2014.06.047.

[11] Huang C-N, Ye Y-H. Development of a water-mist cooling system: A 12,500 Kcal/h air-cooled chiller. Energy Rep 2015;1:123-8. http://dx.doi.org/10.1016/j.egyr.2015.04.002.

[12] Liu H, Zhou Q, Liu Y, Wang P, Wang D. Experimental study on cooling performance of air conditioning system with dual independent evaporative condenser. Int J Refrig 2015;55:85-92. http://dx.doi.org/10.1016/j.ijrefrig.2015.03.012.

[13] Martínez P, Ruiz J, Cutillas CG, Martínez PJ, Kaiser AS, et al. Experimental study on energy performance of a split air-conditioner by using variable thickness evaporative cooling pads coupled to the condenser. Appl Therm Eng 2016;105:1041-50. http://dx.doi.org/10.1016/j.applthermaleng.2016.01.067.

[14] Dreyer AA, Kriel DE, Erens P. Analysis of spray-cooled finned-tube heat exchangers. Heat Transf Eng 1992;13(4):53-71. http://dx.doi.org/10.1080/01457639208939788.

[15] Chen C-W, Yang C-Y, Hu Y-T. Heat transfer enhancement of spray cooling on flat aluminum tube heat exchanger. Heat Transf Eng 2013;34(1):29-36. http://dx.doi.org/10.1080/014576372.2013.694742.

[16] Wataru N, Heikichi K, Shigeki H. Heat transfer from tube banks to air/water mist flow. Int J Heat Mass Transfer 1988;31(2):449-60. http://dx.doi.org/10.1016/0017-9310(88)90027-0.

[17] Kim M-H, Song S, Bullard CW. Effect of inlet humidity condition on the air-side performance of an inclined brazed aluminum evaporator. Int J Refrig 2002;25(5):611-20. http://dx.doi.org/10.1016/S0140-7007(01)00061-5.

[18] Hasan A, Sirén K. Performance investigation of plain and finned tube evaporatively cooled heat exchangers. Appl Therm Eng 2003;23(3):325-40. http://dx.doi.org/10.1016/S1359-4311(02)00194-1.

[19] Elsherbini Al, Jacobi AM. A model for condensate retention on plain-fin heat exchangers. J Heat Transf 2005;128(5):427-33. http://dx.doi.org/10.1115/1.2175091.

[20] Popli S, Hwang Y, Radermacher R. Performance enhancement of enhanced herringbone wavy-fin round tube inclined heat exchangers with and without hydrophilic coating using evaporative spray and deluge cooling. In: ASME international mechanical engineering congress and exposition. 56352, American Society of Mechanical Engineers; 2013. V08BT09A035.

[21] Shin J, Ha S. The effect of hydrophilicity on condensation over various types of fin-and-tube heat exchangers. Int J Refrig 2002;25(6):688-94. http://dx.doi.org/10.1016/S0140-7007(01)00065-1.

[22] Junior ICA, Smith-Schneider P. Consolidated experimental heat and mass transfer database for a reduced scale evaporative condenser. Int J Refrig 2016;66:21-31. http://dx.doi.org/10.1016/j.ijrefrig.2015.12.006.

[23] Vendé P-E. Incidence de l'implémentation d'un Système de Pulvérisation de gouttelettes d'eau Sur l'efficacité d'un groupe de production de Froid Embarqué (Ph.D. thesis). Sorbonne université; 2018.

[24] de Gennes PG. Wetting: Statics and dynamics. Rev Modern Phys 1985;57(3):827-63. http://dx.doi.org/10.1103/RevModPhys.57.827.

[25] Yarin AL. Drop impact dynamics: Splashing, spreading, receding, bouncing.... Annu Rev Fluid Mech 2006;38(Volume 38, 2006):159-92. http://dx.doi.org/10.1146/annurev.fluid.38.050304.092144.

[26] Josserand C, Thoroddsen ST. Drop impact on a solid surface. Annu Rev Fluid Mech 2016;48(Volume 48, 2016):365-91. http://dx.doi.org/10.1146/annurev-fluid-122414-034401.

[27] Zhao S-C, de Jong R, van der Meer D. Raindrop impact on sand: A dynamic explanation of crater morphologies. Soft Matter 2015;11(33):6562-8. http://dx.doi.org/10.1039/C5SM00957J.

[28] Vaikuntanathan V, Sivakumar D. Maximum spreading of liquid drops impacting on groove-textured surfaces: Effect of surface texture. Langmuir 2016;32(10):2399-409. http://dx.doi.org/10.1021/acs.langmuir.5b04639.

[29] Yamamoto K, Takezawa H, Ogata S. Droplet impact on textured surfaces composed of commercial stainless razor blades. Colloids Surfaces A: Physicochem Eng Asp 2016;506:363-70. http://dx.doi.org/10.1016/j.colsurfa.2016.06.024.

[30] Lorenceau É, Quéré D. Drops impacting a sieve. J Colloid Interface Sci 2003;263(1):244-9. http://dx.doi.org/10.1016/S0021-9797(03)00126-7.

[31] Bouchard DJ, Chandra S. Droplet impact and flow into a gap between parallel plates. Phys Fluids 2019;31(6):062104. http://dx.doi.org/10.1063/1.5097944.

[32] Marmur A. Penetration of a small drop into a capillary. J Colloid Interface Sci 1988;122(1):209-19. http://dx.doi.org/10.1016/0021-9797(88)90304-9.

[33] Choi H, Ma L, Liang H. A kinetic study of the spontaneous penetration of a water drop into a hydrophobic pore. Surf Topogr: Metrol Prop 2017;5(1):014003. http://dx.doi.org/10.1088/2051-672X/aa6126.

[34] Delbos A, Lorenceau E, Pitois O. Forced impregnation of a capillary tube with drop impact. J Colloid Interface Sci 2010;341(1):171-7. http://dx.doi.org/10.1016/j.cis.2009.09.016.

[35] Hosseini S, Dalili A, Ashgriz N, Chandra S. Droplet impact and penetration on series of parallel tubes. Engineering Conferences International — ECI Digital Archives; 2014, p. "".

[36] Willmott GR, Neto C, Hendy SC. An experimental study of interactions between droplets and a nonwetting microfluidic capillary. Faraday Discuss 2010;146:233. http://dx.doi.org/10.1039/b925588e.

[37] Willmott GR, Briole A, Szczepaniak F. Inertial capillary uptake of drops. Phys Rev E 2020;101(4):043109. http://dx.doi.org/10.1103/PhysRevE.101.043109.

[38] Dhondi S, Pereira GG, Hendy SC. Effect of molecular weight on the capillary absorption of polymer droplets. Langmuir 2012;28(27):10256-65. http://dx.doi.org/10.1021/la300903w.

[39] Wang Y, Chen S, Liu Y. Spontaneous uptake of droplets into non-wetting capillaries. Comput & Fluids 2016;134-135:190-5. http://dx.doi.org/10.1016/j.compfluid.2016.05.017.

[40] Wang Y, Chen S, Wu B. Self-driven penetration of droplets into non-wetting capillaries. In: ICCFD8, Comput & Fluids In: ICCFD8, 2017;vol. 154:211-5. http://dx.doi.org/10.1016/j.compfluid.2017.06.006.

[41] Zhang D, Shangguan Q, Wang Y. Apply surface wettability gradient to non-wetting capillary: A simulation study on spontaneous droplet flow. AIP Adv 2018;8(9):095223. http://dx.doi.org/10.1063/1.5047450.

[42] Wu Y, Wang F, Huang W, Selzer M, Nestler B. Capillary adsorption of droplets into a funnel-like structure. Phys Rev Fluids 2022;7(5):054004. http://dx.doi.org/10.1103/PhysRevFluids.7.054004.

[43] Yue P, Renardy Y. Spontaneous penetration of a non-wetting drop into an exposed pore. Phys Fluids 2013;25(5):052104. http://dx.doi.org/10.1063/1.4804957.

[44] Ferdowsi P, Bussmann M. Axisymmetric drop penetration into an open straight capillary. In: 20th annual conference of the CFD society of Canada. 2012, p. "".

[45] Bachmann J, Ellies A, Hartge KH. Development and application of a new sessile drop contact angle method to assess soil water repellency. J Hydrol 2000;231-232:66-75. http://dx.doi.org/10.1016/S0022-1694(00)00184-0.

[46] Choi H, Liang H. Wettability and spontaneous penetration of a water drop into hydrophobic pores. J Colloid Interface Sci 2016;477:176-80. http://dx.doi.org/10.1016/j.cis.2016.05.029.

[47] Denesuk M, Smith GL, Zelinski BJJ, Kreidl NJ, Uhlmann DR. Capillary penetration of liquid droplets into porous materials. J Colloid Interface Sci 1993;158(1):114-20. http://dx.doi.org/10.1006/jcis.1993.1235.

[48] Miyanaji H, Zhang S, Yang L. A new physics-based model for equilibrium saturation determination in binder jetting additive manufacturing process. Int J Mach Tools Manuf 2018;124:1-11. http://dx.doi.org/10.1016/j.ijmachtools.2017.09.001.

[49] Seveno D, LeDauphin V, Martic G, Voué M, De Coninck J. Spreading drop dynamics on porous surfaces. Langmuir 2002;18(20):7496-502. http://dx.doi.org/10.1021/la025520h.

[50] Ding H, Theofanous TG. The inertial regime of drop impact on an anisotropic porous substrate. J Fluid Mech 2012;691:546-67. http://dx.doi.org/10.1017/jfm.2011.491.

[51] Bordoloi AD, Longmire EK. Drop motion through a confining orifice. J Fluid Mech 2014;759:520-45. http://dx.doi.org/10.1017/jfm.2014.572.

[52] Marmur A. The radial capillary. J Colloid Interface Sci 1988;124(1):301-8. http://dx.doi.org/10.1016/0021-9797(88)90351-7.

[53] Piroird K, Clanet C, Quéré D. Capillary extraction. Langmuir 2011;27(15):9396-402. http://dx.doi.org/10.1021/la201490m.

[54] Schebarchov D, Hendy SC. Capillary absorption of metal nanodroplets by single-wall carbon nanotubes. Nano Lett 2008;8(8):2253-7. http://dx.doi.org/10.1021/nl080875s.

[55] Schebarchov D, Hendy SC. Uptake and withdrawal of droplets from carbon nanotubes. Nanoscale 2011;3(1):134-41. http://dx.doi.org/10.1039/CONR00477D.

[56] Bormashenko E. Wetting of flat and rough curved surfaces. J Phys Chem C 2009;113(40):17275-7. http://dx.doi.org/10.1021/jp905237v.

[57] Bormashenko E. Wetting transitions on biomimetic surfaces. Philos Trans R Soc A: Math Phys Eng Sci 2010;368(1929):4695-711. http://dx.doi.org/10.1098/rsta.2010.0121.

[58] Bormashenko E. Progress in understanding wetting transitions on rough surfaces. In: Reinhard miller, honorary issue, Adv Colloid Interface Sci In: Reinhard miller, honorary issue, 2015;222:92-103. http://dx.doi.org/10.1016/j.cis.2014.02.009.

[59] Hagg D, Eifert A, Dörr A, Bodziony F, Marschall H. Counter-intuitive penetration of droplets into hydrophobic gaps in theory and experiment. Sci Rep 2023;13(1):16518. http://dx.doi.org/10.1038/s41598-023-43138-2.
[60] Popinet S. Gerris: A tree-based adaptive solver for the incompressible Euler equations in complex geometries. J Comput Phys 2003;190(2):572-600. http://dx.doi.org/10.1016/S0021-9991(03)00298-5.

[61] Popinet S. An accurate adaptive solver for surface-tension-driven interfacial flows. J Comput Phys 2009;228(16):5838-66. http://dx.doi.org/10.1016/j.jcp.2009.04.042.

[62] Lagrée P-Y, Staron L, Popinet S. The granular column collapse as a continuum: Validity of a two-dimensional Navier-Stokes model with a $\mu(I)$-Rheology. J Fluid Mech 2011;686:378-408. http://dx.doi.org/10.1017/jfm.2011.335.

[63] Popinet S. conserving.h: Momentum-conserving advection of velocity. Tech. rep., Basilisk documented source code; 2017-2024, URL http://basilisk.dalembert.upmc.fr/src/navier-stokes/conserving.h.

[64] Arrufat T, Crialesi-Esposito M, Fuster D, Ling Y, Malan L, Pal S, et al. A mass-momentum consistent, volume-of-fluid method for incompressible flow on staggered grids. Comput & Fluids 2021;215:104785.

[65] Tryggvason G, Scardovelli R, Zaleski S. Direct numerical simulations of gas-liquid multiphase flows. Cambridge: Cambridge University Press; 2011, http://dx.doi.org/10.1017/CBO9780511975264.

[66] Schwartz P, Barad M, Colella P, Ligocki T. A cartesian grid embedded boundary method for the heat equation and Poisson's equation in three dimensions. J Comput Phys 2006;211(2):531-50. http://dx.doi.org/10.1016/j.jcp.2005.06.010.

[67] Johansen H, Colella P. A cartesian grid embedded boundary method for Poisson's equation on irregular domains. J Comput Phys 1998;147(1):60-85. http://dx.doi.org/10.1006/jcph.1998.5965.

[68] Tavares M, Josserand C, Limare A, Lopez-Herrera JM, Popinet S. A coupled VOF/Embedded boundary method to model two-phase flows on arbitrary solid surfaces. Comput & Fluids 2024;278:106317. http://dx.doi.org/10.1016/j.compfluid.2024.106317.

[69] Popinet S. Sessile drop on an embedded boundary. 2021, http://basilisk.fr/sandbox/popinet/contact/sessile.c.

[70] Popinet S. Sessile drop on an embedded boundary inclined at 45 degrees. 2021, http://basilisk.fr/sandbox/popinet/contact/sessile-inclined.c.

[71] O'Brien A, Afkhami S, Bussmann M. Pore-scale direct numerical simulation of haines jumps in a porous media model. Eur Phys J Spec Top 2020;229(10):1785-98. http://dx.doi.org/10.1140/epjst/e2020-000008-0.

[72] Hirt CW, Nichols BD. Volume of fluid (VOF) method for the dynamics of free boundaries. J Comput Phys 1981;39(1):201-25. http://dx.doi.org/10.1016/0021-9991(81)90145-5.

[73] Rudman M. A volume-tracking method for incompressible multifluid flows with large density variations. Internat J Numer Methods Fluids 1998;28(2):357-78. http://dx.doi.org/10.1002/(SICI)1097-0363(19980815)28:2<357::AID-FLD750>3.0.CO;2-D.

[74] Huh C, Scriven L. Hydrodynamic model of steady movement of a solid/liq-uid/fluid contact line. J Colloid Interface Sci 1971;35(1):85-101. http://dx.doi.org/10.1016/0021-9797(71)90188-3.

[75] Afkhami S, Zaleski S, Bussmann M. A mesh-dependent model for applying dynamic contact angles to VOF simulations. J Comput Phys 2009;228(15):5370-89. http://dx.doi.org/10.1016/j.jcp.2009.04.027.

[76] Chen X, Lu J, Tryggvason G. Numerical simulation of self-propelled non-equal sized droplets. Phys Fluids 2019;31(5):052107. http://dx.doi.org/10.1063/1.5094757.

[77] Thomas D, Lacour S, Zaleski S. Clogging of a dual-plate channel by a droplet. 2024, http://basilisk.fr/sandbox/thomas/run/channel_clogging.c.