Delamination and out-of-plane deformation in drying colloidal suspensions

Abstract

A drop of a colloidal suspension placed on a substrate forms a solid particle deposit as it dries. As water evaporates, large gradients in pore pressure inside the porous deposit cause shrinkage and stresses. The deposit cracks, then delaminates from the substrate, and bends out of plane, creating a striking three-dimensional structure. Previous models have attributed the out-of-plane deformation to pore pressure gradients through the deposit's thickness, a hypothesis our findings contradict. Through a combination of interference and confocal microscopy, we show that the final curvature strongly depends on the deposit thickness, with thinner deposits curving more. We propose a mechanism where the curvature is driven not by vertical pressure gradients, but by much larger radial pressure gradients across the length of the deposit. The resulting in-plane differential shrinkage creates geometric frustration that is resolved through out-of-plane buckling. We validate this mechanism using non-Euclidean plate simulations, which successfully reproduce the buckling behavior and the observed dependence of curvature on thickness.

---

Introduction

Lilies curve as they bloom, the scales of pine cones bend to open when the weather is dry, and mud in Death Valley curls up as it dries; thin objects that dynamically curve are ubiquitous in the natural world.^1–3 This complex bending behavior can be probed in evaporating drops of aqueous colloidal suspensions. As a drop evaporates on a solid substrate, the colloidal particles assemble into a solid deposit that cracks, delaminates from the substrate, and curves out of plane.^4–9 While the formation of the particle deposit and the conditions for crack formation are well documented,^10,11 an understanding of the out-of-plane deformation of the deposit is lacking.

The contact line of a drop of a colloidal suspension deposited on a hydrophilic substrate remains pinned to the substrate.¹² A capillary-driven flow brings particles to the edge of the drop, leading to the formation of a deposit of close-packed particles that grows radially inward and eventually covers almost the entire initially wetted area.^9,13–15 As the deposit grows, it remains saturated with water, and water lost to evaporation at the top of the deposit is replaced by water from the liquid region.^16,17 The flow of water is sustained by a lower pressure in the pores of the deposit compared to that in the liquid region.^18–20 This negative pressure in the pores causes the deposit to shrink, but shrinkage is initially prevented by the adhesion of the deposit to the substrate, which leads to in-plane tensile drying stresses.^{10,11,21–25} The in-plane stresses increase until a critical stress is reached and cracks form.^20,26,27 As evaporation continues, the drying stresses continue to increase. This prompts two different modes of stress relaxation: if the deposit remains adhered to the substrate, the stresses relax by the formation of additional cracks.^28–30 Alternatively, if the deposit delaminates from the substrate, it becomes free to shrink to relax the drying stresses.^31,32 The delaminated deposit curves out of plane, creating three-dimensional structures that resemble blooming flowers.^4–9 The curvature of the deposit indicates that the shrinkage is not uniform, since a deposit that shrunk uniformly would remain flat.³³ The mechanism by which the delaminated deposit shrinks, whether curvature increases locally or uniformly across the deposit, and the parameters that govern the curvature are largely unknown. Here, we measure the curved shape of the deposit during and after drying, and demonstrate that the curvature is set locally at the advancing solidification front. We show that there is a strong dependence of the curvature on the thickness of the deposit, in contradiction with recent models that account for the bending in terms of a differential shrinkage across the thickness of the deposit.^5,6,8 Instead, we propose that the curvature is set by a differential shrinkage across the length of the deposit, which causes geometric frustration that is resolved by a large out-of-plane deformation.

Materials and methods

Drop drying experiments are performed using charge-stabilized suspensions of colloidal silica particles (Ludox AS-40, Sigma-Aldrich) of diameter 2a = 22 nm. The Ludox AS-40 particles are negatively charged, with a zeta potential of −47 mV reported in ref. 34. The suspensions are diluted with deionized water to particle volume fractions ranging from ϕ₀ = 0.08–0.23. The diluted suspensions are stable. Drops of volume Ω₀ = 0.3–1 µL are deposited on clean microscope glass slides using a micropipette (Gilson 0.2–2 µL), where they adopt a contact angle 20 ± 2°. Prior to use, the glass slides are cleaned by sonication in acetone for five minutes followed by rinsing with isopropanol, and the colloidal suspensions are sonicated for twenty minutes to ensure complete dispersion. The drops dry at ambient temperature T = 22 ± 2 °C and constant relative humidity ranging from 6–60%.

In side-view and bottom-view experiments, we image the drops from below using an inverted microscope (Eclipse TE2000-U, Nikon, with Lumix GH5 camera) with magnification 6×, filming at 30 fps. We obtain interference images by shining a single wavelength λ = 532 nm light source from below the sample.⁹ We image the drops from the side using a Ximea xiQ camera with a f = 50 mm objective and a 50 mm extension tube.

In confocal microscopy experiments, we use a laser scanning confocal microscope (VK-X 3D, Keyence) to measure the upper interface of the drying drop. By scanning in the vertical direction and recording the maximum reflected laser intensity, the location of the interface is identified. This process is repeated along a radius of the drop to measure the top interface z(r).

To identify the mechanisms governing the out-of-plane deformation of the deposit, we perform non-Euclidean plate simulations.³⁵ The deposit is modeled as a thin elastic sheet of thickness h and radius R₀, subjected to a target strain ε_t(r, z) = P(r, z)/Y. This strain represents the local shrinkage the material would undergo to reach a stress-free state, where P(r, z) is the pore pressure determined from Darcy's law and Y is the Young's modulus of the deposit. The sheet's configuration is determined by minimizing the total elastic energy E, described by the Koiter shell model as the sum of stretching energy E_s and bending energy E_b:

In this formulation, a_c and b_c are the first and second fundamental forms of the midsurface in the current configuration, while a_r and b_r represent the rest (target) configuration. The elastic energy norm is defined as ‖A‖²_e = α Tr²(A) + 2β Tr(A²), with α and β the plane-stress Lamé parameters,^36,37 and U denotes the domain of the mid-surface parametrization. To account for vertical gradients in pore pressure, we use a bilayer model where each layer of thickness h/2 is assigned a specific rest metric. The numerical simulations consider a wedge-shaped segment of the deposit – representing the region between two radial cracks – discretized using a triangulated mesh. The equilibrium deformed configuration is then computed by minimizing the discretized elastic energy over all possible embeddings of this mesh, given the prescribed rest configuration of each triangle. We systematically vary the dimensionless groups P*/Y (target strain amplitude) and h/R₀ (thickness ratio) to observe how the sheet resolves the geometric frustration resulting from the non-uniform shrinkage.

Results

Deposit delamination and out-of-plane deformation

We deposit a drop of an aqueous suspension of colloidal silica particles on a clean glass slide. The drop adopts a spherical cap shape of radius R₀, as shown in Fig. 1A. As water evaporates, the initial contact line remains pinned. A deposit of close-packed particles with width w forms from the edge of the drop, enclosing a liquid spherical cap of radius r_l = R₀ − w that shrinks as the solidification front moves inward.^19,38 At a deposit width w_crack, a first radial crack appears in the deposit. The deposit delaminates from the glass substrate around the crack, as evidenced by the interference fringes in Fig. 1B and in Movie S1 (SI).²⁰ The first crack is rapidly followed by a succession of regularly spaced radial cracks that separate the deposit into segments whose shapes resemble the petals of a flower. As the liquid region recedes, the radial cracks propagate inward and the deposit delaminates at a propagating delamination front of radius r*. The drying dynamics are thus characterized by two circular fronts that propagate toward the center of the drop: a solidification front that separates the particle deposit and the liquid region at radius r_l, and a delamination front at radius r*. Both fronts eventually reach the center of the drop as the liquid region disappears and almost the entire deposit is delaminated, with only a small area of contact with the substrate in the very center of the deposit. The petals curve as they delaminate, leading to the out-of-plane deformation visible in Fig. 1A.^9,39

Fig. 1 Drying and out-of-plane deformation of a colloidal suspension drop. (A) Side-view images of a drop of a silica particle suspension with initial volume Ω₀ = 0.3 μL and initial particle volume fraction ϕ₀ = 0.18 drying on a glass slide at a relative humidity of RH = 49% at different times t since drop deposition. A solid deposit forms at the edge of the drop and grows inward as the liquid region of radius r_l in the center recedes. When the deposit width w reaches a critical value w_crack, a first radial crack forms in the deposit, followed by regularly spaced radial cracks. The deposit delaminates from the substrate and bends upward, reaching a vertical deflection δ at the edge of the deposit. (B) Simultaneous bottom-view images of the drop using interference microscopy with a monochromatic light of wavelength λ = 532 nm. The black and white interference fringes indicate that the deposit delaminates from the substrate around the cracks. The delamination front propagates at a radius r* > r_l from the center of the drop. r* is measured at the centerline of the deposit segment. The scale bars represent 0.5 mm.

Dynamics of out-of-plane deformation

We measure the curvature of the petals close to the delamination front using interference microscopy. A monochromatic light source of wavelength λ illuminates the drop from below and induces interference of light reflected at the top of the substrate and light reflected at the bottom of the deposit. The successive white and black interference fringes visible in Fig. 2A denote an increase in the thickness of the air layer between the substrate and the deposit by λ/4. We reconstruct the deflection z of the centerline of a segment of deposit between two cracks close to the delamination front by measuring the location of the interference fringes, knowing that the deposit is adhered to the substrate with z = 0 at the delamination front r*, as shown in Fig. 2B. Interference fringes are visible for deflections z < 6 µm, which corresponds to the region close to the delamination front. This limit is due to the coherence length of the light source and to the angle that the delaminated deposit adopts with the substrate.⁴⁰ The deflection profile close to the delamination front, for z < 6 µm, can be fit to a parabola defined aswhere κ* is the local curvature at the delamination front. The deflection in this region is smaller than the thickness of the deposit. The evolution of κ* with deposit width w = R₀ − r_l is shown in Fig. 2C. Before cracks form at a deposit width w_crack, the deposit is adhered to the substrate and the curvature is everywhere equal to zero. At the onset of delamination, the curvature of the delaminated deposit is small. The curvature at the delamination front progressively increases as the deposit width increases until it reaches a maximum . κ* then remains close to . The increase in curvature with deposit width collapses onto a master curve when we report versus (w − w_crack)/R₀, for experiments at different initial particle volume fractions ϕ₀ and relative humidities RH, as shown in Fig. 2D. The delamination dynamics thus only depend on the deposit width w_crack at which the first crack forms and delamination starts, and on the final curvature .

Fig. 2 Dynamics of out-of-plane bending. (A) Bottom-view interference microscopy image of a drying drop with initial volume Ω₀ = 0.3 µL and initial particle volume fraction ϕ₀ = 0.09 at a relative humidity of RH = 18%. We identify the radial location of the interference fringes along the dashed green line to measure the deflection profile of the delaminated deposit at r > r*. The scale bar represents 0.5 mm. (B) Vertical deposit deflection z versus radius r along a radial line at different times Δt from the onset of delamination. Each data point corresponds to an interference fringe and a vertical increment in z of λ/4 = 133 nm. The curvature at the delamination front κ* is obtained by fitting the deflection to a parabolic profile z = (κ*/2)(r − r*)². (C) Curvature at the delamination front κ* versus normalized deposit width w/R₀ for drops with Ω₀ = 0.3 µL and varying initial particle volume fraction ϕ₀ drying at RH = 52 ± 3%. For increasing ϕ₀, crack formation and delamination occur for larger deposit width w_crack and the maximum curvature is lower. (D) Normalizing the curvature at the delamination front κ* with the maximum curvature and subtracting the deposit width at crack formation w_crack from the deposit width w leads to a master curve of all the data. In addition to the data shown in panel C, experiments at relative humidity RH = 15 ± 3% are also shown. The volume fractions are identified in the legend in panel C.

Curvature set locally

While the curvature κ* is defined locally at the delamination front, the entire delaminated deposit deflects as the delamination front propagates, leading to a large vertical deflection δ at the edge of the deposit. To connect the local curvature at the delamination front with the global shape of the deposit, we determine the shape of the deposit using confocal microscopy. We measure the location of the top interface of the deposit, z_top ≈ z + h, where h is the thickness of the deposit, as the deposit grows and bends upward, as shown in Fig. 3A. The separation between the liquid region and the solid deposit at r_l is clearly visible, as well as the separation between the adhered and bent parts of the deposit. The curves for the top interface of the deposit at different times are superimposed in the region r_l < r < r* where the deposit is formed but not delaminated, indicating that the deposit thickness remains constant before delamination.

Fig. 3 Deposit curvature is set locally. (A) Location of the upper interface z_top of the liquid region and the solid deposit along a fixed radial line at different times t since drop deposition, measured using a laser confocal microscope. The initial drop volume is Ω₀ = 0.3 µL, the initial particle volume fraction is ϕ₀ = 0.14 and the relative humidity is RH = 14%. (B) Superposition of the upper interface z_top in the curved region of the deposit obtained by a rotation and a translation. The superposition indicates that the deposit only deforms in a small region close to the delamination front and then maintains its curvature.

Remarkably, the curvature of the delaminated deposit remains fixed as the deposit continues to grow. The delaminated section of the profiles at each time can be superimposed onto the final profile by a combination of a translation and a rotation that leaves the curvature unchanged, as shown in Fig. 3B and in Movie S2 (SI). This superposition indicates that the deposit only acquires curvature in a small region close to the delamination front. Consequently, the increase in the vertical deflection δ at the tip of the deposit is only due to the propagation of the delamination front toward the center of the drop, and not to a change in curvature in the part of the deposit that is already delaminated. When the delamination front reaches the center of the drop, δ reaches a maximum value and then remains close to this maximum value, as shown in Fig. S1 (SI). Once the liquid region disappears, air replaces water in the pores of the deposit and the drying stresses caused by the water pore pressure vanish.^11,20 The deposit remains curved, indicating that the deformation is plastic and originates from particle rearrangements.⁵ Fig. 3B shows that this plastic deformation is localized in a small region close to the delamination front, and that the curvature obtained when the deposit delaminates is maintained as drying continues.

Effect of deposit thickness on curvature

To connect the local curvature κ* with the global shape of the deformed deposit, we measure the maximum deflection of the tip of the deposit δ_max and define the maximum average curvature [small kappa, Greek, macron]_max = 2δ_max/R₀². The maximum average curvature [small kappa, Greek, macron]_max measured in side-view images is close to the maximum local curvature measured in interference microscopy, as shown in Fig. S2 (SI).

The curvature depends strongly on the initial particle volume fraction ϕ₀, as shown in Fig. 4A and in Movie S3 (SI). The initial particle volume fraction sets the average thickness h̄ of the particle deposit. Particle volume conservation yields h̄ = ϕ₀Ω₀/(ϕ_depositπR₀²), where Ω₀ is the initial drop volume and ϕ_deposit ≈ 0.6 is the particle volume fraction in the deposit, which is close to the random close-packing volume fraction.^9,19 We vary h̄ by changing the initial drop volume between 0.3 µL and 1 µL, and the initial particle volume fraction between 0.08 and 0.23. In addition, we vary the relative humidity and calculate the resulting average evaporative flux j₀ as j₀ = Ω₀/(πR₀²t_f)(1 − ϕ₀/ϕ_deposit), where the drying time t_f is the time when r_l = 0 and the liquid region disappears.²⁰ The deposit thickness h̄ governs the curvature [small kappa, Greek, macron]_max, as shown in Fig. 4B. For thin deposits, the curvature is considerable, of the order of 1 mm⁻¹ ≈ 1/R₀. The average evaporative flux j₀ does not have a noticeable effect on the curvature. The decrease in ionic strength due to diluting the stock suspension with deionized water also does not affect the curvature, as shown in Fig. S3 and Supplementary Text (SI).

Fig. 4 Dependence of deposit curvature on deposit thickness. (A) Side-view image of dried drops with initial volume Ω₀ = 0.3 µL and initial particle volume fractions (from left to right) ϕ₀ = 0.08, ϕ₀ = 0.18 and ϕ₀ = 0.23, corresponding to average deposit thicknesses h̄ = 11 µm, h̄ = 25 µm and h̄ = 32 µm. (B) Maximum average curvature [small kappa, Greek, macron]_max = 2δ_max/R₀² measured from side-view images versus average deposit thickness h̄ for drops with initial volume Ω₀ = 0.3–1 µL, ϕ₀ = 0.08–0.23, drying at relative humidity RH = 6–60%. The color indicates the calculated evaporative flux j₀, which is larger for a lower relative humidity. The green horizontal line at the bottom indicates the target bending curvature κ_t = 0.003–0.01 mm⁻¹ as predicted using the equation in ref. 6. δ_max and R₀ are measured separately for the furthest left and right petals to account for potential asymmetry in the drying. (C) Upper and lower interfaces of the dried deposit for a drop with Ω₀ = 0.3 µL and ϕ₀ = 0.14, measured using a laser confocal microscope in film mode.⁴⁰ The color indicates the local curvature κ of the lower deposit interface. Regions of the deposit with a larger local thickness h have a smaller curvature, while thinner regions of the deposit have a larger curvature.

So far, we have seen that the final curvature of the deposit is set by the deposit thickness, and that the deposit acquires curvature as it deforms locally next to the delamination front. Remarkably, the curvature of the deposit even locally depends on the deposit thickness, as shown in Fig. 4C, where we report both the upper and lower interfaces of the deposit measured using confocal microscopy. Thicker regions of the deposit exhibit a smaller curvature, while thinner regions exhibit a larger curvature.

Discussion

From flow in the deposit to curvature

What sets the drying-induced curvature? The solid particle deposit remains saturated with water during drying, and water lost to evaporation with evaporative flux j₀ is replaced by water from the liquid region in the center of the drop.^16,19 Water flows inside the pores of the deposit in the vertical direction with a vertical flow velocity u_z ∼ j₀ and an in-plane radial flow velocity set by volume conservation u_r ∼ j₀w/h̄, as shown in Fig. 5A.⁸ To overcome viscous dissipation in the porous deposit, this flow velocity requires a pressure gradient given by Darcy's law as ∇P(r, z) = −(μ/k)u, where μ is the dynamic viscosity of water and k is the permeability of the deposit estimated using the Kozeny–Carman equation, as shown in Fig. 5B.^17,41 Integrating Darcy's law, the vertical pressure difference between the top and the bottom of the deposit scales as ΔP_z ∼ −P*h̄²/R₀², where P* = μj₀R₀²/(kh̄) is the porous pressure scale, while the radial pressure difference between the edge of the deposit and the solidification front scales as ΔP_r ∼ −P*w²/R₀². The pressure in the pores of the deposit is balanced against atmospheric pressure by the curvature of the microscopic menisci between particles at the top of the deposit, as shown in the inset of Fig. 5B.^14,17 If the pressure difference across the pores of the menisci surpasses the maximum capillary pressure ∼γ/a with γ the surface tension and a the particle radius, air displaces the microscopic menisci and invades the deposit.^14,17 Due to the nanometer size of a that sets both the maximum negative pressure in the deposit and the gradients in pore pressure via the deposit permeability k ∝ a², the pore pressure is negative and reaches values of −1 to −10 MPa, as detailed in supplementary information (SI).²⁰ This negative pressure induces shrinkage in the deposit: P sets an isotropic target strain ε_t(r, z) = P(r, z)/Y < 0, where Y is the Young's modulus of the deposit. The target strain is the strain that would release all stresses.^30,31,42 When the deposit is adhered to the substrate, that target strain cannot be reached because of the no-strain boundary condition ε = 0 imposed at the bottom of the deposit. This leads to tensile drying stresses σ = Y(ε − ε_t) and to the formation of cracks. When the deposit delaminates, however, it becomes free to deform. Consequently, no new cracks form in the delaminated deposit.

Fig. 5 Flow in the porous deposit causes bending and buckling. (A) Flow in the porous deposit. Water evaporates from the surface of the deposit with evaporative flux j₀, inducing a flow velocity u in the porous deposit. The vertical component of the velocity is exaggerated for visualization purposes. (B) Pressure in the porous deposit. The flow velocity in the porous deposit requires a gradient in pore pressure P(r, z), with a lower pore pressure at the top and at the outer edge of the deposit. The pressure difference between z = 0 and z = h is exaggerated in the sketch. The difference between the pore pressure in the deposit and atmospheric pressure is sustained by the capillary pressure at the microscopic menisci between particles at the top surface of the deposit, as shown in the inset. (C) Simulation result for a sheet with thickness h/R₀ = 15 × 10⁻³ and target strain amplitude P*/Y = 0.3 with r_l/R₀ = 0.4. The vertical deflections are magnified by a factor 100. The dashed line denotes the midline of the sheet, which has a constant curvature κ. (D) Simulation result for a sheet with thickness h/R₀ = 3 × 10⁻³ and target strain amplitude P*/Y = 0.3 with r_l/R₀ = 0.4. The deformed sheet has a non-zero Gaussian curvature and reaches vertical deflections 100 times larger than the sheet in (C). The curvature is measured in the region where r > r_l and z < 6 × 10⁻³R₀. (E) Normalized curvature κR₀ of the midline of the sheet close to r_l versus the target strain amplitude P*/Y. Simulations where the sheet is in the bending regime are represented by circles and simulations where the sheet buckles are represented by triangles. The dashed blue line represents a slope of 1, the dashed purple line represents a slope of 2/3. (F) Regime map of bending (blue circles) and buckling (purple triangles) as a function of simulation parameters h/R₀ and P*/Y. The black line indicates P*/Y = 3 × 10³(h/R₀)².

The target strain in the delaminated deposit is not uniform, which leads to the out-of-plane deformation of the deposit. Previous works predict that the curvature depends linearly on the average evaporative flux j₀ and is independent of the average deposit thickness h̄.^6,8 Both of these predictions are in disagreement with our experimental observations; the curvature we measure is independent of j₀ and decreases with h̄. In addition, the predicted curvature κ_t ∼ μj₀/(kY) is much smaller than the curvature κ measured experimentally: with our experimental values of j₀ = 0.2–0.8 µm s⁻¹ and the value of Y = 150 MPa measured in ref. 6, κ_t = 0.003–0.01 mm⁻¹, a factor of 10 to 300 times less than the curvature in our experiments, as shown in Fig. 4B. The κ_t model only considers vertical flows in the deposit, which cause a vertical pressure difference ΔP_z that leads to a vertical gradient in target strain, with the top of the deposit shrinking more than the bottom. Since κ ≫ κ_t, vertical flows alone do not explain the curvature of the drying drops. In addition to the vertical pressure difference, there is a radial pressure difference ΔP_r due to in-plane radial flows in the deposit. Regions of the deposit further away from the solidification front experience a more negative pore pressure and shrink more. This radial pressure difference is on the order of (w/h̄)² ≈ 400 to 10 000 times larger than the vertical pressure difference; the pressure and the target strain thus vary mainly in the radial direction.⁸ However, an in-plane gradient in shrinkage does not directly correspond to an out-of-plane deformation. To elucidate the discrepancy between the predicted target curvature κ_t and the curvature κ measured in experiments, we conduct non-Euclidian plate simulations that consider the variations in target strain across both the vertical and in-plane radial directions. The simulations model the deposit as a thin elastic sheet of uniform thickness h subject to an isotropic in-plane target strain ε_t(r, z) = P(r, z)/Y that varies in both the vertical and radial directions.³⁷ We model the wedge-shaped region of the deposit between two cracks. The elastic energy per unit area of the sheet is represented using the Koiter shell model containing both stretching energy E_s and bending energy E_b, allowing for a non-Euclidean rest configuration associated with the strain due to shrinkage, as discussed in Fig. S4 and Supplementary Text (SI).^35,37 We minimize the total elastic energy to obtain the shape of the deposit.³⁷

Bending and buckling regimes in simulations

The simulations highlight two different modes of deformation. For low target strain amplitudes P*/Y or large thickness ratios h/R₀, the sheet bends with a constant curvature in the r direction and remains flat with no curvature in the perpendicular direction, as illustrated in Fig. 5C. In contrast, for larger target strains and thinner sheets, the sheet is doubly curved with a non-zero Gaussian curvature, as illustrated in Fig. 5D. The curvature in the radial direction is largest close to r_l and the maximum deflection of the sheet is two orders of magnitude larger than in the first mode of deformation. The curvature of sheets that maintain zero Gaussian curvature follows κR₀ = κ_tR₀ ∼ μj₀R₀/(kY) ∼ (P*/Y)(h/R₀), as shown by the circles in Fig. 5E and Fig. S5 (SI), and in disagreement with our experimental results. The deformation in this first regime comes from vertical gradients in target strain and minimizes the bending energy E_b ∝ (κ − κ_t)². Because the curvature stems from minimizing the bending energy, we refer to this regime as the bending regime. In contrast, sheets that exhibit non-zero Gaussian curvature have a much larger normalized curvature κR₀ along the radial direction, as shown in Fig. 5E. In this second regime, the deformation stems from geometric frustration due to in-plane radial variations in target strain: the sheet cannot shrink more for larger r than for smaller r while staying flat or cylindrically curved, and instead it must acquire a Gaussian curvature as a result of Gauss' Theorema Egregium.^33,37,42 Since geometry prevents a flat or cylindrically curved sheet from satisfying the large variations in target strain along the radial direction, the sheet in its initial configuration has a large stretching energy E_s. Minimizing this stretching energy results in a large out-of-plane deformation with κ ≫ κ_t and non-zero Gaussian curvature. The same mechanism is at play in the blooming of a lily, where the growth of plant tissues at different rates leads to Gaussian curvature.² Because the out-of-plane deformation stems from in-plane effects, we denote this regime the buckling regime. Note that the bending energy E_b plays opposite roles in the bending regime and in the buckling regime: in the bending regime, minimizing E_b causes the sheet to curve to κ_t; in the buckling regime, minimizing E_b penalizes the large curvature κ ≫ κ_t and reduces the curvature.

The simulations recover the experimentally observed dependence of the curvature on the thickness of the deposit in two ways. First, the local curvature in the buckling regime increases with the target strain amplitude P*/Y = μj₀R₀²/(Ykh), as shown in Fig. 5E, and thus decreases with h. Second, the crossover from the bending to the buckling regime occurs for a critical (P*/Y)_crit ∼ (h/R₀)², as shown in Fig. 5F, which can be written in terms of a critical thickness h_crit/R₀ ∼ (μj₀R₀/k)^1/3. For h > h_crit, the sheet bends with a small curvature κ_t; for h < h_crit, the sheet buckles with a large curvature. A similar boundary between bending and buckling for a different variation in target strain was obtained in ref. 43 by setting E_b = E_s. In addition, the simulations indicate that only sheets in the buckling regime can reach curvatures comparable to the normalized curvatures κR₀ = 0.1–1 observed in experiments. Because the vertical variation in target strain ΔP_z/Y that causes bending and the in-plane variation in target strain ΔP_r/Y that causes buckling are proportional, a sheet in the bending regime subjected to a progressively increasing target strain buckles before reaching a curvature comparable to the one observed in experiments. This confirms that bending purely due to vertical flows in the porous deposit cannot explain the curvatures observed in the experiments.

The simulations model an elastic sheet and do not capture the local and plastic nature of the deposit deformation highlighted in Fig. 3. In addition, the shrinkage in simulations for r close to R₀ is large, up to ε = −0.2, while the shrinkage we observe in experiments is much smaller. Since the buckling curvature in the simulations is localized in a small region where the strain is small, we perform additional simulations setting a maximum negative target strain of −0.05. These simulations recover the same localized curvature, indicating that the buckling curvature is set locally, as shown in Fig. S6 (SI). This suggests the scenario for buckling-induced curvature illustrated in Fig. 6: in the region between r_l and r*, the deposit has formed but remains adhered to the substrate, preventing any in-plane strain and curvature. In a small intermediate region between the delamination front at r* and a hypothetical radius r_p, the target strain ε_t varies rapidly with r causing the deposit to buckle and acquire a large curvature. This deformation is plastic and the deposit retains the curvature it acquires in the intermediate region. Further away from the delamination front for r > r_p, the deposit does not deform, likely because a maximum shrinkage −ε_max has been reached.

Fig. 6 Buckling occurs locally. (A) Sketch of the local buckling mechanism with three distinct regions. The deposit in the region between the solidification front r_l and the delamination front r* is adhered to the substrate and does not shrink in plane, ε = 0. The deposit in the intermediate delaminated region just past the delamination front shrinks in plane with large radial variation in target strain, resulting in buckling and curvature. This deformation is plastic and the deposit maintains its curvature for r > r_p as defined in the sketch. (B) Representation of the three regions on a side-view image of a drying drop. The width of the intermediate region where the strain varies and the deposit curvature increases is indicated approximately.

Buckling explains the large curvatures observed in experiments and the dependence on the deposit thickness. However, because buckling depends on the pore pressure, which is set by flows in the deposit, the buckling curvature is expected to depend on the evaporative flux j₀ – a dependence that is not observed experimentally and points to limitations in our simulations. Two important aspects of drop drying are not captured in the simulations. The first aspect is the adhesion of the deposit to the substrate, which penalizes delamination. The adhesion energy per unit area E_adh is minimized in addition to the elastic energies E_b and E_s. The ratio of the adhesion energy and the elastic energy sets the location of the delamination front r*, which is distinct from the solidification front r_l, while in simulations without adhesion energy r* = r_l. The master curve of versus (w − w_crack)/R₀ shown in Fig. 2 can be interpreted as a transition from a regime dominated by adhesion energy where κ* = 0 to a regime dominated by the in-plane stretching energy where . When the deposit first cracks, it delaminates but remains almost flat. As (w − w_crack)/R₀ increases, the pore pressure becomes more negative away from the solidification front and the in-plane gradients in target strain increase, causing the deposit to buckle out of plane with a curvature that approaches as E_adh/E_s decreases towards zero. The second aspect of drop drying that could explain the discrepancy between our model and the experiments is the microstructure of the deposit, which we begin to consider by imposing a maximum negative target strain in simulations. The deposit is composed of nanoparticles pushed together by the radial flow from the liquid region to the deposit. The drag force of the flow through the deposit causes anisotropy, in-plane plastic deformation and shear banding in the deposit.^44–46 Modeling the deposit directly as an assembly of spheres in Hertzian contact leads to a non-linear stress–strain constitutive relation for the deposit, which could affect the scaling of the curvature obtained in simulations.^47,48 The mechanical properties of the deposit also vary from the solidification front to the edge of the deposit.¹⁹ In centimeter-thick samples of drying clay, the drying-induced curvature has been associated with vertical gradients in particle size, and for drops of a nanoparticle suspension on hydrophobic surfaces vertical gradients in particle concentration have been observed.^49–51 A vertical gradient in the microstructure could explain part of the large curvatures we observe. For our monodisperse suspensions and thin deposits, however, we expect particles to smooth out vertical gradients in particle concentration by diffusion.⁵² Experiments that vary the particle modulus, particle size and adhesion between the deposit and the substrate could provide further insights.

Conclusions

Out-of-plane deformation of thin films is seen from peeling paint to the morphogenesis of biological tissues. Our findings fundamentally change the understanding of this process in drying deposits of colloidal suspensions. We find that the deposit's curvature is set locally at the delamination front, with a final shape that is strongly dependent on thickness—a key parameter previously unaccounted for. We demonstrate that the deformation is a buckling phenomenon, driven by dominant in-plane differential shrinkage from radial pressure gradients. The stunning array of colloidal flowers created by deposits of various thicknesses serves as transient testimony of the intricate nature of drying phenomena.

Conflicts of interest

There are no conflicts to declare.

Data availability

All data needed to evaluate the conclusions in the paper are present in the paper and/or the supplementary information (SI). Supplementary information: details of the simulation and supplementary movies. See DOI: https://doi.org/10.1039/d5sm01099c.

Acknowledgements

We thank Rohan Abeyaratne, Philippe Bourrianne, Efi Efrati and Matthew G. Hennessy for interesting discussions. I. B. acknowledges support from the American Chemical Society Petroleum Research Fund, grant 66034-ND9. P. L. acknowledges funding from the MathWorks Foundation and the Ben Gold Fellowship Fund. M. I. acknowledges funding from the Meryl and Stewart Robertson Fund.

Notes and references

1. C. Dawson, J. F. V. Vincent and A.-M. Rocca, Nature, 1997, 390, 668.
2. H. Liang and L. Mahadevan, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 5516–5521.
3. R. W. Style, S. S. L. Peppin and A. C. F. Cocks, J. Geophys. Res.: Earth Surf., 2011, 116, F01025.
4. F. Parisse and C. Allain, Phys. Fluids, 1996, 8, S6.
5. F. Giorgiutti-Dauphiné and L. Pauchard, Colloids Surf., A, 2015, 466, 203–209.
6. H. Lama, L. Pauchard, F. Giorgiutti-Dauphine and S. Khawas, Phys. Rev. Mater., 2023, 7, 025604.
7. P. Lilin, P. Bourrianne, G. Sintès and I. Bischofberger, Phys. Rev. Fluids, 2020, 5, 110511.
8. A. Osman, L. Goehring, H. Stitt and N. Shokri, Soft Matter, 2020, 16, 8345–8351.
9. P. Bourrianne, P. Lilin, G. Sintès, T. Nîrca, G. H. McKinley and I. Bischofberger, Soft Matter, 2021, 17, 8832–8837.
10. A. F. Routh, Rep. Prog. Phys., 2013, 76, 046603.
11. L. Goehring, A. Nakahara, T. Dutta, S. Kitsunezaki and S. Tarafdar, Desiccation Cracks and Their Patterns: Formation and Modelling in Science and Nature, John Wiley & Sons, 2015.
12. R. D. Deegan, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2000, 61, 475–485.
13. R. D. Deegan, O. Bakajin, T. F. Dupont, G. Huber, S. R. Nagel and T. A. Witten, Nature, 1997, 389, 827–829.
14. H. M. van der Kooij, G. T. van de Kerkhof and J. Sprakel, Soft Matter, 2016, 12, 2858–2867.
15. G. Berteloot, A. Hoang, A. Daerr, H. P. Kavehpour, F. Lequeux and L. Limat, J. Colloid Interface Sci., 2012, 370, 155–161.
16. J. M. Salamanca, E. Ciampi, D. A. Faux, P. M. Glover, P. J. McDonald, A. F. Routh, A. C. I. A. Peters, R. Satguru and J. L. Keddie, Langmuir, 2001, 17, 3202–3207.
17. E. R. Dufresne, E. I. Corwin, N. A. Greenblatt, J. Ashmore, D. Y. Wang, A. D. Dinsmore, J. X. Cheng, X. S. Xie, J. W. Hutchinson and D. A. Weitz, Phys. Rev. Lett., 2003, 91, 224501.
18. A. F. Routh and W. B. Russel, AIChE J., 1998, 44, 2088–2098.
19. L. Goehring, W. J. Clegg and A. F. Routh, Langmuir, 2010, 26, 9269–9275.
20. P. Lilin and I. Bischofberger, Langmuir, 2022, 38, 7442–7447.
21. L. Goehring, W. J. Clegg and A. F. Routh, Phys. Rev. Lett., 2013, 110, 024301.
22. R. C. Chiu, T. J. Garino and M. J. Cima, J. Am. Ceram. Soc., 1993, 76, 2257–2264.
23. H. Lei, J. A. Payne, A. V. McCormick, L. F. Francis, W. W. Gerberich and L. E. Scriven, J. Appl. Polym. Sci., 2001, 81, 1000–1013.
24. D. M. Holmes, R. Vasant Kumar and W. J. Clegg, J. Am. Ceram. Soc., 2006, 89, 1908–1913.
25. K. Zhang, J. Yan, Q. He, C. Xu, J. van de Koppel, B. Wang, B. Cui and Q.-X. Liu, Sci. Adv., 2023, 9, eabq3520.
26. W. Man and W. B. Russel, Phys. Rev. Lett., 2008, 100, 198302.
27. H. N. Yow, M. Goikoetxea, L. Goehring and A. F. Routh, J. Colloid Interface Sci., 2010, 352, 542–548.
28. A. Kumar and G. U. Kulkarni, Langmuir, 2021, 37, 13141–13147.
29. P. Lilin and I. Bischofberger, Phys. Rev. Fluids, 2022, 7, 110505.
30. P. Lilin, M. Ibrahim and I. Bischofberger, Sci. Adv., 2024, 10, eadp3746.
31. V. Lazarus, EPL, 2017, 117, 24002.
32. M. Liu, S. Yu, L. He and Y. Ni, Soft Matter, 2022, 18, 5906–5927.
33. P. G. Ciarlet, J. Elasticity, 2005, 78, 1–215.
34. K. Manne, S. Magkiriadou, P. Yazhgur, A. F. Demirörs and F. Scheffold, Front. Soft Matter, 2023, 3, 1224158.
35. E. Efrati, E. Sharon and R. Kupferman, J. Mech. Phys. Solids, 2009, 57, 762–775.
36. C. Weischedel, A. Tuganov, T. Hermansson, J. Linn and M. Wardetzky, Construction of Discrete Shell Models by Geometric Finite Differences, Fraunhofer ITWM Technical Report 220, 2012.
37. W. M. van Rees, E. Vouga and L. Mahadevan, Proc. Natl. Acad. Sci. U. S. A., 2017, 114, 11597–11602.
38. C. N. Kaplan and L. Mahadevan, J. Fluid Mech., 2015, 781, R2.
39. F. Giorgiutti-Dauphiné and L. Pauchard, Eur. Phys. J. E: Soft Matter Biol. Phys., 2014, 37, 1–7.
40. J. Schmit, K. Creath and J. C. Wyant, Optical Shop Testing, John Wiley & Sons, Ltd, 2007, ch. 15, pp. 667–755.
41. C. J. Brinker and G. W. Scherer, Sol-Gel Science: The Physics and Chemistry of Sol-Gel Processing, Academic Press, New York, 1990.
42. Y. Klein, E. Efrati and E. Sharon, Science, 2007, 315, 1116–1120.
43. Y. Klein and E. Sharon, Phys. Rev. Lett., 2021, 127, 105501.
44. F. Boulogne, L. Pauchard, F. Giorgiutti-Dauphiné, R. Botet, R. Schweins, M. Sztucki, J. Li, B. Cabane and L. Goehring, EPL, 2014, 105, 38005.
45. B. Yang, J. S. Sharp and M. I. Smith, ACS Nano, 2015, 9, 4077–4084.
46. P.-C. Kiatkirakajorn and L. Goehring, Phys. Rev. Lett., 2015, 115, 088302.
47. W. B. Russel, N. Wu and W. Man, Langmuir, 2008, 24, 1721–1730.
48. O. Prakash Bamboriya and M. S. Tirumkudulu, Soft Matter, 2025, 21, 2986–2993.
49. M. Zielinski, M. Sánchez, E. Romero and A. Atique, Geoderma, 2014, 226–227, 85–93.
50. M. Maedo, M. Sánchez, D. Aljeznawi, O. Manzoli, L. J. N. Guimarães and P. R. Cleto, Acta Geotechnol., 2020, 15, 2619–2635.
51. S. Basu, L. Bansal and A. Miglani, Soft Matter, 2016, 12, 4896–4902.
52. F. Giorgiutti-Dauphiné and L. Pauchard, Eur. Phys. J. E: Soft Matter Biol. Phys., 2018, 41, 32.

---