Time-dependent modelling of thin poroelastic films drying on deformable plates

Abstract

Understanding the generation of mechanical stress in drying, particle-laden films is important for a wide range of industrial processes. One way to study these stresses is through the cantilever experiment, whereby a thin film is deposited onto the surface of a thin plate that is clamped at one end to a wall. The stresses that are generated in the film during drying are transmitted to the plate and drive bending. Mathematical modelling enables the film stress to be inferred from measurements of the plate deflection. The aim of this paper is to present simplified models of the cantilever experiment that have been derived from the time-dependent equations of continuum mechanics using asymptotic methods. The film is described using nonlinear poroelasticity and the plate using nonlinear elasticity. In contrast to Stoney-like formulae, the simplified models account for films with non-uniform thickness and stress. The film model reduces to a single differential equation that can be solved independently of the plate equations. The plate model reduces to an extended form of the Föppl-von Kármán (FvK) equations that accounts for gradients in the longitudinal traction acting on the plate surface. Consistent boundary conditions for the FvK equations are derived by resolving the Saint-Venant boundary layers at the free edges of the plate. The asymptotically reduced models are in excellent agreement with finite element solutions of the full governing equations. As the Péclet number increases, the time evolution of the plate deflection changes from $t$ to $t^{1/2}$, in agreement with experiments.

1. Introduction

The drying of thin films that consist of a volatile solvent and a particulate phase is relevant to a wide range of industrial applications [24], such as the fabrication of electrodes for lithium-ion batteries [41] and flexible conductive coatings [35]. As solvent is removed from the mixture by evaporation, the particulate phase aggregates to form a porous solid. The evaporation-driven flow of solvent through the pores will generate a pressure gradient that, in turn, can deform the solid. The removal of solvent from the pore space will also cause the solid to contract. If the film has been placed on a rigid substrate, then the adhesion of the solid to the substrate will restrict contractions in the horizontal plane. As a result, mechanical stresses will develop within the film [9]. The film stress can be relieved through a myriad of mechanical instabilities [16], including fracture [5, 10], buckling [31] and delamination [15, 30]. Although drying-induced instabilities have traditionally been viewed as detrimental, there is increasing interest in understanding how these can be harnessed in applications such as medical diagnostics [34] and lithography [23].

Avoiding or triggering instabilities during film drying requires quantitative knowledge about the evolution of the film stress. However, directly measuring the mechanical stress during drying is difficult. A number of innovative experiments have therefore been designed around the aim of indirectly measuring the film stress [4, 40]. The cantilever technique is one such experiment, whereby a thin film (or drop) is deposited onto the upper surface of a thin plate that is supported at one end by a wall [12]. Due to the adhesion between the film and the plate, the mechanical stresses that are generated in the film during drying are transmitted to the plate. In response, the plate will bend, causing a deflection that can be measured. The film stress can then be inferred from the deflection using mathematical modelling.

In the seminal work of Stoney [36], a simple model is presented for the cantilever experiment. Using classical beam theory, Stoney was able to determine a relationship between the film stress and the beam deflection. Stoney’s derivation relies on a number of key assumptions, namely the film is thin compared to the beam, the film has a uniform and constant thickness, and bending is driven by a uniform longitudinal film stress. In the context of drying, the latter point implies that the film has a homogeneous composition. Moreover, Stoney’s derivation requires the explicit definition of a neutral axis, which is a longitudinal axis in the beam where no strain occurs. A number of authors have extended Stoney’s model by utilising more sophisticated beam theories [14]. For example, Petersen et al. [32] accounted the finite thickness of the film. However, Chiu [6] has criticised several of these models and argues that Stoney’s choice of neutral axis is only correct when the film and beam have the same Young’s modulus. As discussed by Francis et al. [14], the Young’s modulus of the film is typically one or two orders of magnitude smaller than that of the plate. At the other end of the modelling spectrum, Lei et al. [26] used the finite element method to simulate the drying of an elasto-viscoplastic film on a cantilever beam. To the best of our knowledge, no comparison has been made between simple models based on beam theory and finite element solutions of the full equations of continuum mechanics.

By analysing a cantilever experiment with a Stoney-like model, Croll [8] observed that the stress in the film is independent of the film height. Similar observations were recently made for drying polymer films by Tomar et al. [39], who also noted that the kinetics of stress generation are sensitive to the Péclet number. The film stress was found to increase linearly with time for small Péclet numbers and nonlinearly for larger Péclet numbers. The nonlinear evolution of the stress was attributed to the formation of a polymer-rich skin at the film surface caused by a large rate of evaporation. When a skin forms, the longitudinal stress across the film will be highly non-uniform and the applicability of Stoney’s model to this scenario comes into question. Skin formation is particularly problematic when manufacturing electrodes for lithium-ion batteries, as stress and composition gradients can lead to mechanical failure and poor electrical conductivity [13, 20, 21]. Understanding the mechanics of skin formation is therefore crucial to improving battery performance.

The aim of this paper is to revisit the derivation of Stoney-like formulae using asymptotic methods. Starting from the time-dependent equations of three-dimensional continuum mechanics, a hierarchy of simplified models for the cantilever experiment will be derived. More specifically, the film will be modelled using nonlinear poroelasticity [3, 29] and the plate using nonlinear elasticity [19]. The asymptotic analysis is based on the small aspect ratios of the film and the plate. By using asymptotic methods to simplify the equations of elasticity rather than directly appealing to beam theory, the explicit definition of a neutral axis can be avoided. The reduced models that are presented here extend Stoney-like models by capturing the evolution of the film thickness, composition and stress, all of which can be highly non-uniform. Moreover, the analysis reveals parameter regimes where the plate deflection is driven by both the longitudinal and transverse stress in the film. The reduced models are validated by comparing them against finite element solutions of the full system of equations and found to be very accurate.

A key feature of the asymptotic reduction is that it leads to a one-way decoupling of the models, thus allowing the film and plate equations to be solved sequentially. Moreover, the dimensionality of the film and plate models are reduced. The combined result is a computationally efficient model that is able to resolve the spatio-temporal dynamics of drying and bending. In asymptotically reducing the equations for the film, we invoke a poroelastic version of lubrication theory [11, 18, 22]. Furthermore, we extend our previous work on drying-induced stresses in poroelastic drops [17] by accounting for a wide range of drying regimes that are determined by the magnitude of the Péclet number. In all cases, the film model reduces to either a nonlinear ordinary differential equation or a nonlinear diffusion equation for the solvent concentration. The deformation of the plate can be described by an extended form of the Föppl-von Kármán (FvK) equations [25] that accounts for gradients in the longitudinal traction exerted on the upper plate surface by the film. In order to derive a consistent set of boundary conditions at the free edges of the plate, a careful analysis of the Saint-Venant boundary layers is required [19]. If considering a beam rather than a plate, we show that the FvK equations reduce to a linear second-order differential equation.

The governing equations for the film and the plate are presented in Section 2. A scaling analysis and preliminary reduction of the plate equations are carried out in Section 3. The equations are non-dimensionalised in Section 4. A detailed asymptotic reduction is then carried out in Section 5. In Section 6, the solutions of the asymptotically reduced models are compared against finite element solutions of the full model. A parametric study is also conducted. The paper concludes in Section 7.

2. Modelling

The formation of a coating is a complex process that involves the transformation of a liquid-like mixture into a solid. The generation of mechanical stress occurs after the gel point of the mixture has been crossed, and the mechanical response of the film becomes more solid-like than liquid-like. In cantilever experiments, the crossing of the gel point can be identified as the time at which the plate begins to bend. Before this time, film drying takes place but the deflection of the plate is negligible.

Given that we are interested in modelling the deflection of the plate, we will assume that the gel point has just been crossed so that a poroelastic matrix has formed throughout the film. One-dimensional models that couple solidification and poromechanics have been considered by Style and Peppin [37] and Punati and Tirumkudulu [33]. We therefore consider a poroelastic film that is drying on a thin, flexible plate, as illustrated in Figure 1. The plate is assumed to be rectangular with length $L$ , width $W$ and thickness $H_p$ . The length and width are assumed to be comparable in size. However, the thickness of the plate is assumed to be much smaller than the length and width. The small aspect ratio of the plate is denoted by $\delta = H_p/ L \ll 1$ .

Figure 1. Bending of a plate during film drying. The plate has length $L$ , width $W$ and height $H_p$ . The initial film thickness is given by $H_f(X_1, X_2)$ , where $X_1$ and $X_2$ are Lagrangian coordinates that lie in the plane spanned by the plate centreline (dashed line). The plate is clamped to a wall at $X_1 = 0$ . The origin of the vertical Lagrangian coordinate $Z$ coincides with the plate centreline. Panels (a) and (b): Cross-sections of the Lagrangian (undeformed) configuration. Panel (c): Cross-section of the Eulerian (deformed) configuration.

We let $X_1$ , $X_2$ , and $X_3 = Z$ denote Cartesian coordinates associated with the initial (Lagrangian) configuration of the system. The coordinates $X_1$ and $X_2$ are chosen to lie in the plane that is parallel to the upper and lower surfaces of the plate. The vertical (or transverse) Lagrangian coordinate, $Z$ , is chosen such that $Z = 0$ corresponds to the centre surface of the plate. Consequently, $Z = \pm H_p/2$ denotes the upper and lower surfaces. The plate is envisioned as being clamped to a wall at $X_1 = 0$ and having free edges at $X_1 = L$ and $X_2 = \pm W/ 2$ .

The film has a non-uniform thickness that evolves in time. The initial film thickness is denoted by $H_f(X_1, X_2)$ and has a maximum value given by $H_f^0 = \max (H_f)$ . The film is also assumed to be thin relative to the length of the plate so that $\epsilon = H_f^0/ L \ll 1$ . We assume that the film thickness tends to zero as the edges of the plate are approached. That is, the contact line of the film is assumed to be static, in agreement with experiments [5, 30], and pinned to the edges of the plate. Recession of the contact line during drying is typically prohibited by the strong adhesion of the solid constituents to the substrate [16, 38]. Our modelling setup is similar to the experiments of Tomar et al. [39], where a drop of polymer solution was deposited into the middle of a cantilever plate and allowed to dry.

The film is assumed to remain bonded to the plate during drying. The upper surface of the plate, originally located at $Z = H_p/2$ , becomes vertically displaced to the Eulerian position $z = H_p/2 + w(x_1, x_2, t)$ , where $x_1$ and $x_2$ are in-plane coordinates associated with the deformed (Eulerian) configuration and $t$ is time. The Eulerian thickness of the film is denoted by $h_f(x_1, x_2,t)$ . Hence, the Eulerian position of the free surface of the film is given by $z = H_p/ 2 + w(x_1, x_2,t) + h_f(x_1, x_2,t)$ .

2.1. Notation

Occasionally, vector and tensor quantities will be decomposed into in-plane and vertical components that are parallel and perpendicular to the undeformed plate, respectively. We let $\boldsymbol{{e}}_1$ , $\boldsymbol{{e}}_2$ , and $\boldsymbol{{e}}_3 = \boldsymbol{{e}}_z$ denote Cartesian basis vectors for the $X_1$ , $X_2$ , and $X_3 = Z$ directions, as well as the $x_1$ , $x_2$ , and $z = x_3$ directions. If $\boldsymbol{{a}} = a_i \boldsymbol{{e}}_i$ denotes an arbitrary vector, then we can write $\boldsymbol{{a}} = \boldsymbol{{a}}_{\parallel } + a_z \boldsymbol{{e}}_z$ , where $\boldsymbol{{a}}_{\parallel } = a_\alpha \boldsymbol{{e}}_\alpha$ is defined as the vector of in-plane components and $a_z = a_3$ is the vertical component. Einstein summation notation is used, and we adopt the convention that Greek indices are equal to 1 or 2. We let $\nabla$ denote the material gradient taken with respect to the Lagrangian coordinates $\boldsymbol{{X}} = X_i \boldsymbol{{e}}_i$ . The in-plane gradient operator is defined as $\nabla _{\parallel } = \nabla - \boldsymbol{{e}}_z\,\partial/\partial Z$ . Tensors as written as ${\sf\boldsymbol T} = {\sf\boldsymbol T}_{\parallel } + {\sf\boldsymbol T}_{\perp } \otimes \boldsymbol{{e}}_z +{\textsf{T}}_{z\alpha } \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_\alpha +{\textsf{T}}_{zz} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z$ where ${\sf\boldsymbol T}_{\parallel } ={\textsf{T}}_{\alpha \beta } \boldsymbol{{e}}_\alpha \otimes \boldsymbol{{e}}_\beta$ and ${\sf\boldsymbol T}_\perp ={\textsf{T}}_{\alpha z} \boldsymbol{{e}}_\alpha$ .

2.2. Kinematics

The governing equations are formulated in terms of Lagrangian coordinates $\boldsymbol{{X}} = X_i \boldsymbol{{e}}_i$ associated with the initial (undeformed) configuration of the film and plate. We let $\boldsymbol{{x}} = x_i(\boldsymbol{{X}},t) \boldsymbol{{e}}_i$ denote Eulerian coordinates associated with the current (deformed) configuration. During drying, the solid element originally located at $\boldsymbol{{X}}$ is displaced to $\boldsymbol{{x}}$ , thereby generating a displacement $\boldsymbol{{u}} = \boldsymbol{{x}}(\boldsymbol{{X}},t) - \boldsymbol{{X}}$ . The deformation gradient tensor ${\sf\boldsymbol F}$ describes the distortion of material elements and is given by

(2.1) \begin{align} {\sf\boldsymbol F} = \nabla \boldsymbol{{x}} = {\sf\boldsymbol I} + \nabla \boldsymbol{{u}}, \end{align}

where ${\sf\boldsymbol I}$ is the identity tensor. We adopt the convention that the gradient of a vector $\boldsymbol{{a}} = a_i \boldsymbol{{e}}_i$ is given by $\nabla \boldsymbol{{a}} = (\partial a_{i}/\partial X_j)\,\boldsymbol{{e}}_i\otimes \boldsymbol{{e}}_j$ . The determinant of ${\sf\boldsymbol F}$ , denoted by $J = \det {\sf\boldsymbol F}$ , accounts for volumetric changes of material elements.

2.3. A model for a poroelastic film

Poroelastic materials consist of a porous and deformable solid matrix that is filled with fluid. The theory of poroelasticity was first developed by Biot [3] in the context of soil mechanics. Coussy [7] provides a comprehensive overview of the theory. We treat the drying film as a poroelastic material and describe it using a Lagrangian version of the model proposed by MacMinn et al. [29]. In essence, the model couples the equations of nonlinear elasticity to those for flow in a porous medium.

Conservation of liquid in the porous film leads to

(2.2a) \begin{align} \frac{\partial \Phi }{\partial t} + \nabla \cdot \boldsymbol{{Q}} &= 0, \end{align}

where $\Phi$ is the nominal volume fraction of fluid and $\boldsymbol{{Q}}$ is the nominal fluid flux. The Eulerian volume fraction of fluid, which is equivalent to the film porosity, is given by $\phi = \Phi/ J$ . The transport of fluid within the pore space is governed by Darcy’s law, which can be written in terms of the reference configuration as

(2.3) \begin{align} \boldsymbol{{Q}} &= -{{\sf\boldsymbol {K}}} \nabla p, \end{align}

where ${\sf\boldsymbol {K}}$ is a permeability tensor and $p$ is the fluid pressure. The permeability tensor can be written in terms of the scalar permeability $k(\phi )$ as ${\sf\boldsymbol {K}} = (k(\phi )/ \mu _f) J {\sf\boldsymbol {C}}^{-1}$ , where $\mu _f$ is the fluid viscosity and ${\sf\boldsymbol {C}} = {\sf\boldsymbol F}^T {\sf\boldsymbol F}$ is the right Cauchy–Green tensor. The factor of $J {\sf\boldsymbol {C}}^{-1}$ in ${\sf\boldsymbol {K}}$ is the result of mapping Darcy’s law in the current configuration to the reference configuration. We let $k_0 = k(\phi _0)$ denote the initial permeability, where $\phi _0$ is the initial volume fraction (or porosity) of the film. For simplicity, we assume that the initial fluid fraction is spatially uniform.

At the microscopic level, the solid and fluid phases are assumed to be incompressible. However, at the macroscopic level, the film is compressible, with volumetric changes being accommodated by rearrangements of the pore geometry and the removal (or addition) of fluid from material elements. The connection between macroscopic volume changes and the amount of fluid in the pore space of a material element is captured through the incompressibility condition

(2.4) \begin{align} J = 1 + \Phi - \phi _{0}. \end{align}

Evaluating (2.4) at $t = 0$ gives $J = 1$ . Since drying leads to the removal of fluid from the pore space, $\Phi$ decreases relative to $\phi _0$ , and we expect that $J \leq 1$ . If the composition of the film remains uniform during drying, then $J$ can be written in terms of the total film volume $V(t)$ as $J = V(t)/ V(0)$ . In light of this, we will refer to $J$ as the contraction ratio.

Conservation of linear and angular momentum in the film leads to

(2.5) \begin{align} \nabla \cdot {\sf\boldsymbol S} &= \textbf{0},\\[-9pt]\nonumber \end{align}

(2.6) \begin{align} {\sf\boldsymbol S} {\sf\boldsymbol F}^T &= {\sf\boldsymbol F} {\sf\boldsymbol S}^T. \end{align}

The first Piola–Kirchhoff (PK1) stress tensor ${\sf\boldsymbol S}$ is decomposed as ${\sf\boldsymbol S} = \boldsymbol{\mathsf{\Sigma }} - p J {\sf\boldsymbol F}^{-T}$ , where the first component, $\boldsymbol{\mathsf{\Sigma }}$ , represents the effective (Terzaghi) elastic stress of the solid. The second contribution to ${\sf\boldsymbol S}$ accounts for the stress exerted by the fluid. The solid matrix is assumed to be isotropic and obey a neo-Hookean equation of state. The elastic component of the stress tensor can be written as

(2.7) \begin{align} \boldsymbol{\mathsf{\Sigma }} = \frac{\nu E_f}{(1 + \nu _f)(1-2\nu _f)} J (J - 1){\sf\boldsymbol F}^{-T} + \frac{E_f}{2(1+\nu _f)}\left({\sf\boldsymbol F} - {\sf\boldsymbol F}^{-T}\right), \end{align}

where $E_f$ and $\nu _f$ are the Young’s modulus and Poisson’s ratio of the film, respectively. Both of these parameters are assumed to remain constant during the drying process. In the limit of small deformations, $\nabla \boldsymbol{{u}} \ll 1$ , we find that ${\sf\boldsymbol F} \sim {\sf\boldsymbol I} + \nabla \boldsymbol{{u}}$ , ${\sf\boldsymbol F}^{-T} \sim {\sf\boldsymbol I} - (\nabla \boldsymbol{{u}})^T$ , and $J = \det {\sf\boldsymbol F} \sim 1 + \nabla \cdot \boldsymbol{{u}}$ . Hence, the stress-strain relation (2.7) reduces to

(2.8) \begin{align} \boldsymbol{\mathsf{\Sigma }} \sim \frac{\nu E_f}{(1 + \nu _f)(1-2\nu _f)}(\nabla \cdot \boldsymbol{{u}}){\sf\boldsymbol I} + \frac{E_f}{2(1+\nu _f)}\,\left (\nabla \boldsymbol{{u}} + (\nabla \boldsymbol{{u}})^T\right ), \end{align}

thus recovering linear elasticity. When carrying out a scaling analysis in Section 3, it will be convenient to write the stress balance (2.5) in component form as

(2.9a) \begin{align} \nabla _{\parallel } \cdot {\sf\boldsymbol S}_{\parallel } + \frac{\partial {\sf\boldsymbol S}_{\perp }}{\partial Z} = 0,\\[-9pt]\nonumber \end{align}

(2.9b) \begin{align} \frac{\partial \textsf{S}_{z \alpha }}{\partial X_\alpha } + \frac{\partial \textsf{S}_{zz}}{\partial Z} = 0, \end{align}

where the definitions of ${\sf\boldsymbol S}_{\parallel }$ and ${\sf\boldsymbol S}_{\perp }$ can be found in Section 2.1.

2.4. A model for a deformable plate

Conservation of linear and angular momentum for the plate leads to

(2.10a) \begin{align} \nabla \cdot {\sf\boldsymbol S} &= \textbf{0}, \end{align}

(2.10b) \begin{align} {\sf\boldsymbol S} {\sf\boldsymbol F}^T &= {\sf\boldsymbol F} {\sf\boldsymbol S}^T, \end{align}

where ${\sf\boldsymbol S}$ is the PK1 stress tensor, ${\sf\boldsymbol F} = {\sf\boldsymbol I} + \nabla \boldsymbol{{u}}$ is the deformation gradient tensor, and $\boldsymbol{{u}}$ is the plate displacement. The mechanical response of the plate is described using the Saint-Venant–Kirchhoff constitutive relation given by

(2.11) \begin{align} {\sf\boldsymbol S} = {\sf\boldsymbol F} {\sf\boldsymbol {P}}, \quad {\sf\boldsymbol {P}} = \frac{\nu _p E_p}{(1+\nu _p)(1-2\nu _p)}\textrm{tr}\,( {\sf\boldsymbol {E}}) {\sf\boldsymbol I} + \frac{E_p}{1+\nu _p}{\sf\boldsymbol {E}}, \end{align}

where $\nu _p$ and $E_p$ are the Poisson’s ratio and Young’s modulus of the plate, respectively; ${\sf\boldsymbol {P}}$ is the second Piola–Kirchhoff stress tensor, and ${\sf\boldsymbol {E}} = (1/2) ({\sf\boldsymbol F}^T {\sf\boldsymbol F} - {\sf\boldsymbol I})$ is the strain tensor.

2.5. Boundary and initial conditions

In the reference configuration, the film thickness is constant in time and given by $H_f(\boldsymbol{{X}}_{\parallel })$ . The position of the film surface is then $Z =H_f + H_p/ 2$ . At the film surface, we assume that fluid evaporates with volumetric flux $V_e$ that depends on the Eulerian fluid fraction. We therefore impose

(2.12) \begin{align} \boldsymbol{{Q}} \cdot \boldsymbol{{N}} = V_e(\phi ) \mathcal{S}({\sf\boldsymbol F}), \quad Z =H_f(\boldsymbol{{X}}_{\parallel }) + H_p/ 2, \end{align}

where $\mathcal{S}({\sf\boldsymbol F}) = J |{\sf\boldsymbol F}^{-T} \cdot \boldsymbol{{N}}|$ is a dimensionless function that accounts for how the differential area element differs between the reference and deformed configurations. The unit normal vector can be written as $\boldsymbol{{N}} = \mathcal{N}^{-1}({-}\nabla _{\parallel } H_f + \boldsymbol{{e}}_z)$ where $\mathcal{N} = (1 + |\nabla _{\parallel } H_f|^2)^{1/2}$ . In addition, we assume that the film surface is stress-free:

(2.13) \begin{align} {\sf\boldsymbol S} \cdot \boldsymbol{{N}} = \textbf{0}, \quad Z = H_f(\boldsymbol{{X}}_{\parallel }) + H_p/ 2. \end{align}

The adhesion between the film and the plate is assumed to be perfect; consequently, neither slip nor delamination can occur. The assumption of perfect adhesion manifests as continuity of displacement,

(2.14) \begin{align} \left .\boldsymbol{{u}}\right |_{-} = \left .\boldsymbol{{u}}\right |_{+}, \quad Z = H_p/2, \end{align}

where $-$ means approaching $Z = H_p/2$ from the plate and $+$ means approaching $Z = H_p/ 2$ from the film. The traction exerted by the film on the plate is denoted by $\boldsymbol{\tau }$ . Continuity of stress across the plate-film interface can therefore be expressed as

(2.15) \begin{align} \left .{\sf\boldsymbol S}\cdot \boldsymbol{{e}}_z \right |_{-} = \boldsymbol{\tau } = \left .{\sf\boldsymbol S}\cdot \boldsymbol{{e}}_z\right |_{+}, \quad Z = H_p/2. \end{align}

The plate is assumed to be impermeable, and hence, the following no-flux condition is imposed:

(2.16) \begin{align} \boldsymbol{{Q}} \cdot \boldsymbol{{e}}_z = 0, \quad Z = H_p/ 2. \end{align}

The lower surface of the plate is assumed to be stress-free, resulting in

(2.17) \begin{align} {\sf\boldsymbol S} \cdot \boldsymbol{{e}}_z = \textbf{0}, \quad Z = -H_p/ 2. \end{align}

At the contact surface between the plate and the wall, zero-displacement conditions are imposed:

(2.18) \begin{align} \boldsymbol{{u}} = \textbf{0}, \quad X_1 = 0. \end{align}

The free edges of the plate are stress-free; thus,

(2.19a) \begin{align} {\sf\boldsymbol S} \cdot \boldsymbol{{e}}_1 &= \textbf{0}, \quad X_1 = L, \end{align}

(2.19b) \begin{align} {\sf\boldsymbol S} \cdot \boldsymbol{{e}}_2 &= \textbf{0}, \quad X_2 = \pm W/2. \end{align}

The initial condition for the nominal fluid fraction is given by $\Phi (\boldsymbol{{X}}_{\parallel },Z,0) = \phi _0$ .

3. Scaling analysis and reduction of the plate equations

The amplitude of the plate deflection can be estimated by carrying out a scaling analysis of the governing equations. This estimate will guide the asymptotic reduction of the plate model. To facilitate the scaling analysis, the film and plate are assumed to behave as linear elastic materials with stress-strain relations that are, respectively, identical or analogous to (2.8). In addition, it will be assumed that the deflection of the plate is driven by the in-plane traction, in line with Stoney-like models [12, 36]. We will show below that these two assumptions can be broken in the following ways. The vertical traction can play a comparable role to the in-plane traction. The shear stresses in the film, $\textsf{S}_{\alpha z}$ and $\textsf{S}_{z \alpha }$ , can have different orders of magnitude. These cases will be specifically addressed when constructing asymptotic solutions to the governing equations. The strategy behind the scaling analysis is to estimate the size of the in-plane traction from the film equations. Then, the magnitude of the plate displacements will be estimated by examining the plate equations.

To begin, we consider the stress balance in the film (2.9a). The size of the shear stress in the film can be found by noting that vertical gradients in the shear stress $\partial {\sf\boldsymbol S}_{\perp }/\partial Z \sim {\sf\boldsymbol S}_{\perp }/ H_f^0$ must balance in-plane gradients of the in-plane stress $\nabla _{\parallel } \cdot {\sf\boldsymbol S}_{\parallel } \sim E_f/ L$ . Consequently, the shear stress and in-plane traction scale like ${\sf\boldsymbol S}_{\perp } \sim \epsilon E_f$ and $\boldsymbol{\tau }_{\parallel } \sim \epsilon E_f$ . Applying similar arguments to the vertical component of the stress balance (2.9a) and assuming $\textsf{S}_{z \alpha } \sim \textsf{S}_{\alpha z} \sim \epsilon E_f$ shows that $\textsf{S}_{zz} \sim \epsilon ^2 E_f$ and hence $\tau _{z} \sim \epsilon ^2 E_f$ . Thus, under these assumptions, the vertical traction is much smaller than the in-plane traction.

The in-plane traction exerted on the surface of the plate generates an internal shear stress; thus, ${\sf\boldsymbol S}_{\perp } \sim \boldsymbol{\tau }_{\parallel } \sim \epsilon E_f$ . As in the film, the vertical gradients in the shear stress, $\partial {\sf\boldsymbol S}_{\perp }/\partial Z \sim {\sf\boldsymbol S}_{\perp }/ H_p$ , are balanced by the in-plane gradients of the in-plane stresses, $\nabla _{\parallel } \cdot {\sf\boldsymbol S}_{\parallel } \sim {\sf\boldsymbol S}_{\parallel }/ L$ , leading to ${\sf\boldsymbol S}_{\parallel } \sim \delta ^{-1} \epsilon E_f$ . A linear stress-strain relation implies that ${\sf\boldsymbol S}_{\parallel } \sim E_p \boldsymbol{{u}}_{\parallel }/ L$ . Thus, we find that the in-plane displacement in the plate scales as $\boldsymbol{{u}}_{\parallel } \sim \delta ^{-1} (E_f/ E_p) H_f^0$ . A scale for the vertical displacement of the plate is obtained by assuming that bending is the primary mode of deformation and balancing the components of the shear strains, $\boldsymbol{{u}}_{\parallel }/ H_p \sim u_z/ L$ , which leads to $u_z \sim \delta ^{-2} (E_f/ E_p) H_f^0$ . Since the displacements in the film and plate must match at the film-plate surface, these displacement scales also apply to the film.

The scaling analysis motivates introducing the non-dimensional parameter $\mathcal{E} \equiv \delta ^{-2} E_f/ E_p$ , which compares the film modulus $E_f$ to a reduced modulus for the plate, $\delta ^2 E_p$ . Thus, $\mathcal{E}$ characterises how stiff the film is relative to the plate. The displacements will therefore scale as $\boldsymbol{{u}}_{\parallel } \sim \delta \mathcal{E} H_f^0$ and $u_z \sim \mathcal{E} H_f^0$ . The size of $\mathcal{E}$ , therefore, plays a key role in determining the magnitudes of the displacements and controls the mechanics of the film and plate. For films that are thin relative to the plate, $\epsilon \ll \delta$ , three regimes of plate mechanics can be identified:

1. (i) Soft films: $\mathcal{E} = O(1)$ . The deflection of the plate scales like the film thickness, which is much smaller than the plate thickness. Linear plate theory can be applied.

2. (ii) Stiff films: $\mathcal{E} = O(\delta \epsilon ^{-1}) \gg 1$ . The plate deflection is proportional to its thickness. Föppl-von Kármán theory applies to the plate.

3. (iii) Very stiff films: $\mathcal{E} \gg O(\delta \epsilon ^{-1}) \gg 1$ . The deflection of the plate greatly exceeds its thickness. This regime is beyond the validity of Föppl-von Kármán theory.

The model reduction presented below will focus on the first two of these regimes. An analysis of Regime (iii), in which the film is very stiff, is left as an area of future work. The case when the film and the plate have similar thicknesses, i.e. $\delta = O(\epsilon )$ as $\epsilon \to 0$ , is a distinguished limit. In this case, Regimes (i) and (ii) coincide and the plate can be modelled using a modified form of the Föppl-von Kármán (FvK) equations. The reduced equations for the plate are formulated by assuming that $\delta = O(\epsilon )$ and are then specialised to the first two regimes when $\epsilon \ll \delta$ in Section 5.

3.1. Reduced equations for the plate

The mechanics of the plate are described using a modified form of the Föppl-von Kármán (FvK) equations that accounts for non-uniform in-plane tractions. These equations are systematically derived from the equations of nonlinear elasticity in Appendix A.1. The leading-order contribution to the vertical displacement of the plate, which we denote by $w$ , is independent of $Z$ and satisfies

(3.1a) \begin{align} -B \nabla _{\parallel }^4 w + \nabla _{\parallel } \cdot (\bar{{\sf\boldsymbol S}}_{\parallel } \nabla _{\parallel } w) &= -\frac{H_p}{2}\nabla _{\parallel } \cdot \boldsymbol{\tau }_{\parallel } - \tau _{z}, \end{align}

where $B = \bar{E}_p H_p^3/ 12$ is the bending modulus and $\bar{E}_p = E_p/ (1 - \nu _p^2)$ is the effective Young’s modulus. The first term on the right-hand side of (3.1a) is absent from traditional formulations of the FvK equations; see, e.g. Landau and Lifshitz [25] or Howell et al. [19]. The mean in-plane displacement $\bar{\boldsymbol{{u}}}_{\parallel }$ can be obtained by solving the mean in-plane stress balance given by

(3.1b) \begin{align} \nabla _{\parallel } \cdot \bar{{\sf\boldsymbol S}}_{\parallel } &= -\boldsymbol{\tau }_{\parallel }, \end{align}

where the mean in-plane stress $\bar{{\sf\boldsymbol S}}_{\parallel }$ and strain $\bar{{\sf\boldsymbol {E}}}_{\parallel }$ are defined as

(3.1c) \begin{align} \bar{{\sf\boldsymbol S}}_{\parallel } &= H_p \bar{E}_p\left [(1 - \nu _p) \bar{{\sf\boldsymbol {E}}}_{\parallel } + \nu _p \textrm{tr}\,(\bar{{\sf\boldsymbol {E}}}_{\parallel }) {\sf\boldsymbol I}_{\parallel }\right ],\\[-9pt]\nonumber \end{align}

(3.1d) \begin{align} \bar{{\sf\boldsymbol {E}}}_{\parallel } &= \frac{1}{2}\left [\nabla _{\parallel } \bar{\boldsymbol{{u}}}_{\parallel } + (\nabla _{\parallel } \bar{\boldsymbol{{u}}}_{\parallel })^T + \nabla _{\parallel } w \otimes \nabla _{\parallel } w\right ]. \end{align}

The in-plane displacements are linked to the mean in-plane and vertical displacements via

(3.2) \begin{align} \boldsymbol{{u}}_\parallel = \bar{\boldsymbol{{u}}}_{\parallel } -Z \nabla _{\parallel } w. \end{align}

3.2. Boundary conditions at the edges of the plate

At the clamped edge of the plate, located at $X_1 = 0$ , we impose

(3.3a) \begin{align} \bar{\boldsymbol{{u}}}_{\parallel } = 0, \quad w = 0, \quad \frac{\partial w}{\partial X_1} = 0; \quad X_1 = 0. \end{align}

The boundary conditions at a free edge can be derived by resolving the mechanics of the plate in a thin boundary layer where the stresses become large, details of which are provided in Appendix A.2. The boundary conditions at $X_1 = L$ are

(3.4a) \begin{align} \bar{\textsf{S}}_{\alpha 1} &= 0, &\quad X_1 &= L;\\[-9pt]\nonumber \end{align}

(3.4b) \begin{align} \frac{\partial ^2 w}{\partial X_1^2} + \nu _p \frac{\partial ^2 w}{\partial X_2^2} &= 0, &\quad X_1 &= L;\\[-9pt]\nonumber \end{align}

(3.4c) \begin{align} B \left [\frac{\partial ^3w}{\partial X_1^3} + (2-\nu _p)\frac{\partial ^3w}{\partial X_1 \partial X_2^2}\right ] &= \frac{H_p}{2}\mathcal{\tau }_{1}, &\quad X_1 &= L. \end{align}

Similarly, the boundary conditions at $X_2 = \pm W/ 2$ are

(3.5a) \begin{align} \bar{\textsf{S}}_{\alpha 2} &= 0, &\quad X_2 &= \pm W/2;\\[-9pt]\nonumber \end{align}

(3.5b) \begin{align} \nu _p \frac{\partial ^2 w}{\partial X_1^2} + \frac{\partial ^2 w}{\partial X_2^2} &= 0, &\quad X_2 &= \pm W/2;\\[-9pt]\nonumber \end{align}

(3.5c) \begin{align} B\left [(2-\nu _p)\frac{\partial ^3w}{\partial X_1^2 \partial X_2} + \frac{\partial ^3w}{\partial X_2^3}\right ] &= \frac{H_p}{2}\mathcal{\tau }_{2}, &\quad X_2 &= \pm W/2. \end{align}

4. Non-dimensionalisation

The governing equations for the film and plate are non-dimensionalised by scaling the in-plane coordinates as $\boldsymbol{{x}}_{\parallel } \sim L$ . Time is scaled as $t \sim \tau _{\text{pe}}$ , where the poroelastic time scale $\tau _{\text{pe}} = \mu _f L^2/ (k_0 E_f)$ describes the time required for a pressure gradient of magnitude $E_f/ L$ to transport fluid a distance $L$ through a porous medium with permeability $k_0$ .

4.1. Film model

In the film, the vertical coordinate is scaled according to $Z \sim H_f^0$ . The components of the fluxes are scaled as $\boldsymbol{{Q}}_\parallel \sim (k_0/ \mu _f) E_f/ L$ and $Q_z \sim \epsilon (k_0/ \mu _f) E_f/ L$ . The film displacements and stresses are scaled according to the estimates from Section 3; thus, $\boldsymbol{{u}}_{\parallel } \sim \delta \mathcal{E} H_f^0$ , $u_z \sim \mathcal{E} H_f^0$ , ${\sf\boldsymbol S}_{\parallel } \sim E_f$ , ${\sf\boldsymbol S}_{\perp } \sim \epsilon E_f$ , $\textsf{S}_{z \alpha } \sim \epsilon E_f$ , and $\textsf{S}_{zz} \sim \epsilon ^2 E_f$ . The pressure and the elastic stress tensor are scaled as $p \sim E_f$ and $\boldsymbol{\mathsf{\Sigma }} \sim E_f$ .

Under this choice of scales, the conservation law for the nominal fluid fraction (2.2a) becomes

(4.1) \begin{align} \frac{\partial \Phi }{\partial t} + \nabla _{\parallel } \cdot \boldsymbol{{Q}}_\parallel + \frac{\partial Q_z}{\partial Z} = 0. \end{align}

By writing the symmetric permeability tensor as ${\sf\boldsymbol {K}} = {\sf\boldsymbol {K}}_{\parallel } + {\sf\boldsymbol {K}}_{\perp } \otimes \boldsymbol{{e}}_z + \boldsymbol{{e}}_z \otimes {\sf\boldsymbol {K}}_{\perp } +{\textsf{K}}_{zz} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z$ , the non-dimensional components of the fluid flux (2.3) are

(4.2a) \begin{align} \boldsymbol{{Q}}_\parallel &= -{\sf\boldsymbol {K}}_{\parallel } \nabla _{\parallel } p - \epsilon ^{-1} {\sf\boldsymbol {K}}_{\perp } \frac{\partial p}{\partial Z},\\[-9pt]\nonumber \end{align}

(4.2b) \begin{align} Q_z &= -\epsilon ^{-1} {\sf\boldsymbol {K}}_{\perp } \cdot \nabla _{\parallel } p - \epsilon ^{-2}{\textsf{K}}_{zz} \frac{\partial p}{\partial Z}. \end{align}

The stress balances representing conservation of linear momentum (2.9) are

(4.3a) \begin{align} \nabla _{\parallel } \cdot {\sf\boldsymbol S}_{\parallel } + \frac{\partial {\sf\boldsymbol S}_{\perp }}{\partial Z} = 0,\\[-9pt]\nonumber \end{align}

(4.3b) \begin{align} \frac{\partial S_{z\alpha }}{\partial X_\alpha } + \frac{\partial \textsf{S}_{zz}}{\partial Z} = 0, \end{align}

where the PK1 stress tensor is given by

(4.4) \begin{align} {\sf\boldsymbol S}_{\parallel } + \epsilon \left ({\sf\boldsymbol S}_{\perp } \otimes \boldsymbol{{e}}_z + \textsf{S}_{z \alpha } \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_\alpha \right ) + \epsilon ^2 \textsf{S}_{zz} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z = \boldsymbol{\mathsf{\Sigma }} - p J {\sf\boldsymbol F}^{-T}. \end{align}

The elastic stress tensor (2.7) is written as

(4.5) \begin{align} \boldsymbol{\mathsf{\Sigma }} = a(\nu _f) (J - 1) J {\sf\boldsymbol F}^{-T} + b(\nu _f)\left({\sf\boldsymbol F} - {\sf\boldsymbol F}^{-T}\right), \end{align}

where the functions $a$ and $b$ are defined, for convenience, as

(4.6) \begin{align} a(\nu ) \equiv \frac{\nu }{(1 + \nu )(1 - 2 \nu )}, \quad b(\nu ) \equiv \frac{1}{2(1 + \nu )}. \end{align}

The deformation gradient tensor (2.1) is given by

(4.7) \begin{align} {\sf\boldsymbol F} = {\sf\boldsymbol I} + \mathcal{E} \left (\delta \epsilon \nabla _{\parallel } \boldsymbol{{u}}_{\parallel } + \delta \frac{\partial \boldsymbol{{u}}_{\parallel }}{\partial Z} \otimes \boldsymbol{{e}}_z + \epsilon \boldsymbol{{e}}_z \otimes \nabla _{\parallel } u_z + \frac{\partial u_z}{\partial Z} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z\right ). \end{align}

The appearance of $\mathcal{E}$ , $\epsilon$ and $\delta$ in the deformation gradient tensor leads to distinct asymptotic regimes that are defined by the relative sizes of these parameters. Finally, the incompressibility condition is given by (2.4).

4.2. Plate model

In the plate, the displacements are scaled according to $\boldsymbol{{u}}_{\parallel } \sim \delta \mathcal{E} H_f^0$ , $\bar{\boldsymbol{{u}}}_{\parallel } \sim \delta \mathcal{E} H_f^0$ and $w \sim \mathcal{E} H_f^0$ . The mean in-plane strains and stresses are scaled as $\bar{{\sf\boldsymbol {E}}}_{\parallel } \sim \delta \epsilon \mathcal{E}$ and $\bar{{\sf\boldsymbol S}}_{\parallel } \sim \delta \epsilon \mathcal{E} H_p E_p$ . The tractions are scaled as $\boldsymbol{\tau }_{\parallel } \sim \epsilon E_f$ and $\tau _{z} \sim \epsilon ^2 E_f$ ; see Section 3.

The rescaled FvK equation (3.1) are then given by

(4.8a) \begin{align} -\mathcal{B} \nabla _{\parallel }^4 w + \epsilon \delta ^{-1} \mathcal{E} \nabla _{\parallel } \cdot (\bar{{\sf\boldsymbol S}}_{\parallel } \nabla _{\parallel } w) &= -\frac{1}{2}\nabla _{\parallel } \cdot \boldsymbol{\tau }_{\parallel } - \epsilon \delta ^{-1} \tau _{z}, \end{align}

(4.8b) \begin{align} \nabla _{\parallel } \cdot \bar{{\sf\boldsymbol S}}_{\parallel } &= -\boldsymbol{\tau }_{\parallel }, \end{align}

where $\mathcal{B} = 1/[12(1-\nu _b^2)]$ is a non-dimensional bending modulus. The rescaled mean stresses and strains, obtained from (3.1c) and (3.1d), are

(4.8c) \begin{align} \bar{{\sf\boldsymbol S}}_{\parallel } &= \frac{1}{1-\nu _p^2}\left [(1 - \nu _p) \bar{{\sf\boldsymbol {E}}}_{\parallel } + \nu _p \textrm{tr}\,(\bar{{\sf\boldsymbol {E}}}_{\parallel }) {\sf\boldsymbol I}_{\parallel }\right ],\\[-9pt]\nonumber \end{align}

(4.8d) \begin{align} \bar{{\sf\boldsymbol {E}}}_{\parallel } &= \frac{1}{2}\big[\nabla _{\parallel } \bar{\boldsymbol{{u}}}_{\parallel } + (\nabla _{\parallel } \bar{\boldsymbol{{u}}}_{\parallel })^T + \epsilon \delta ^{-1} \mathcal{E} \nabla _{\parallel } w \otimes \nabla _{\parallel } w\big]. \end{align}

By scaling $Z \sim H_p$ in the plate, the leading-order contribution to the in-plane displacement (3.2) is given by

(4.9) \begin{align} \boldsymbol{{u}}_{\parallel } = \bar{\boldsymbol{{u}}}_{\parallel } - Z \nabla _{\parallel } w. \end{align}

4.3. Boundary conditions

In the context of the film model, in which $H_f^0$ has been used to non-dimensionalise the vertical coordinate $Z$ , the plate has a non-dimensional thickness of $\mathcal{H} = \delta \epsilon ^{-1}$ . The upper surface of the plate is therefore located at $Z = \mathcal{H}/ 2$ ; the free surface of the film is located at $Z = H_f(\boldsymbol{{X}}_{\parallel }) + \mathcal{H}/2$ .

The boundary conditions at the free surface of the film, obtained from (2.12) to (2.13), are given by

(4.10a) \begin{align} -\nabla _{\parallel } H_f \cdot \boldsymbol{{Q}}_\parallel + Q_z &= \textrm{Pe}\, \mathcal{V}(\phi )\, \mathcal{S}({\sf\boldsymbol F})\,\mathcal{N},\\[-9pt]\nonumber \end{align}

(4.10b) \begin{align} {\sf\boldsymbol S}_{\perp } - {\sf\boldsymbol S}_{\parallel } \nabla _{\parallel } H_f &= 0,\\[-9pt]\nonumber \end{align}

(4.10c) \begin{align} \textsf{S}_{zz} - \textsf{S}_{z \alpha } \frac{\partial H_f}{\partial X_\alpha } &= 0, \end{align}

where the Péclet number is defined as $\textrm{Pe} = V_e^0 H_f^0 \mu _f/ (\epsilon ^2 k_0 E_f)$ , where $V_e^0$ is a characteristic evaporative flux. Moreover, $\mathcal{V} = V_e/ V_e^0$ denotes the composition-dependent non-dimensional evaporative flux and $\mathcal{N} = (1 + \epsilon ^2 |\nabla _{\parallel } H_f|^2)^{1/2}$ . At the upper plate surface, $Z = \mathcal{H}/ 2$ , the no-flux condition on the fluid (2.16) can be written as

(4.11) \begin{align} Q_z|_{+} = 0. \end{align}

Using (4.9), continuity of displacement (2.14) can be written as

(4.12a) \begin{align} \left .\boldsymbol{{u}}_{\parallel }\right |_{+} &= \bar{\boldsymbol{{u}}}_{\parallel } - (1/2) \nabla _{\parallel } w,\\[-9pt]\nonumber \end{align}

(4.12b) \begin{align} \left .u_z\right |_{+} &= w. \end{align}

The non-dimensionalised boundary conditions at the clamped and free edges of the plate are identical to those written in Section 3.2. However, the free edges are now located at $x_1 = 1$ and $x_2 = \pm \mathcal{W}/2$ , where $\mathcal{W} = W/ L$ ; $B$ can be replaced with $\mathcal{B}$ , and $H_p$ can be set to one.

The tractions can be found by evaluating the normal components of the PK1 stress in the film at $Z = \mathcal{H}/ 2$ to give $\boldsymbol{\tau }_{\parallel } = {\sf\boldsymbol S}_{\perp } \cdot \boldsymbol{{e}}_z|_{+}$ and $\tau _{z} = \textsf{S}_{zz}|_{+}$ . However, equivalent expressions that are more convenient to use can be found by integrating the stress balances (4.3) across the film thickness, resulting in

(4.13a) \begin{align} \boldsymbol{\tau }_{\parallel } &= \nabla _{\parallel } \cdot \int _{\mathcal{H}/2}^{H_f + \mathcal{H}/ 2} {\sf\boldsymbol S}_{\parallel }\, \textrm{d} Z,\\[-9pt]\nonumber \end{align}

(4.13b) \begin{align} \tau _{z} &= \frac{\partial }{\partial X_\alpha }\int _{\mathcal{H}/2}^{H_f + \mathcal{H}/ 2} \textsf{S}_{z \alpha } \, \textrm{d} Z. \end{align}

5. Asymptotic reductions of the coupled film-plate model

The coupled system of equations for the film and the plate can be asymptotically reduced by taking the limit $\epsilon \to 0$ . Substantial simplifications occur through a decoupling of the models. In particular, the equations for the film can be solved independently of those for the plate. The problem for the film also simplifies through a decoupling of the mechanical and fluid-transport problems. This decoupling occurs because the deformation gradient tensor, to leading order, is simply ${\sf\boldsymbol F} =\textrm{diag}(1, 1, J)$ . Thus, the mechanical problem can be readily solved in terms of $J = 1 + \Phi - \phi _0$ and used to formulate a closed-form problem for $\Phi$ . Using the asymptotically reduced model in practice therefore involves three steps:

1. 1. Solve a reduced model for the nominal fluid fraction $\Phi$ and compute the local volumetric contraction $J$ .

2. 2. Compute the tractions, which can be expressed in terms of $J$ and $\nabla _{\parallel } w$ , and substitute them into the FvK equations.

3. 3. Solve the FvK equations for the plate displacements.

Despite the simplified pathway to obtaining a solution, there are several asymptotic regimes in both the mechanical and fluid-transport problems for the film that must be delicately handled. These regimes arise from the various sizes that $\mathcal{E}$ and $\textrm{Pe}$ can take. The former alters the structure of the mechanical problem whereas the latter alters the fluid-transport problem. In terms of the mechanical problem, we consider three cases that can be summarised as follows:

• $\delta = O(\epsilon )$ with $\mathcal{E} = O(1)$ : the deformation of the film is driven by a combination of volumetric contraction and plate bending. Both the in-plane and vertical tractions drive bending.

• $\delta \gg \epsilon$ with $\mathcal{E} = O(1)$ : in-plane film deformation is dominated by plate bending. The vertical film deformation driven by volumetric contraction and plate bending. Only the in-plane traction drives bending.

• $\delta \gg \epsilon$ with $\mathcal{E} = O(\delta \epsilon ^{-1})$ : the deformation of the film is dominated by plate bending. Both the in-plane and vertical traction drive bending.

In terms of the drying dynamics, there are three cases to consider:

• $\textrm{Pe} = O(1)$ : the nominal fluid fraction is uniform in the vertical direction. The flow of fluid predominantly occurs in the in-plane directions. The case $\textrm{Pe} \ll 1$ is a sub-limit in which the film remains homogeneous during drying; that is, the fluid fraction is spatially uniform.

• $\textrm{Pe} = O(\epsilon ^{-1})$ : the nominal fluid fraction is uniform in the vertical direction. Fluid flow is two-dimensional, with the in-plane and vertical flux components having similar orders of magnitude.

• $\textrm{Pe} = O(\epsilon ^{-2})$ : the nominal fluid fraction is non-uniform in the vertical direction. Fluid flow mainly occurs in the vertical direction.

In Section 5.1, asymptotic solutions to the mechanical problems are detailed for the three cases described above. In Section 5.2, the fluid-transport problem is asymptotically reduced according to the size of the Péclet number.

5.1. Reduction of the mechanical problems

Asymptotic expansions are used to solve the governing equations for the film mechanics. The goal is to obtain expressions for the traction. As the plate model has already been reduced, we do not asymptotically expand the solutions to the FvK equations. Instead, we comment on which terms can be neglected depending on the sizes of $\delta$ and $\mathcal{E}$ relative to $\epsilon$ .

5.1.1. Films and plates of similar thicknesses

We first consider the distinguished limit in which the film and the plate have similar thicknesses by letting $\delta = \delta _0 \epsilon$ , where $\delta _0 = O(1)$ as $\epsilon \to 0$ . To ensure that the FvK limit remains valid, we also assume that $\mathcal{E} = O(1)$ as $\epsilon \to 0$ . With this choice of $\delta$ , the deformation gradient tensor for the film (4.7) can be asymptotically expanded as ${\sf\boldsymbol F} = {\sf\boldsymbol F}^{(0)} + \epsilon {\sf\boldsymbol F}^{(1)} + O(\epsilon ^2)$ , where

(5.1a) \begin{align} {\sf\boldsymbol F}^{(0)} &= {\sf\boldsymbol I}_{\parallel } + \left(1 + \mathcal{E} \frac{\partial u_z^{(0)}}{\partial Z}\right) \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z,\qquad\qquad\qquad\quad\end{align}

(5.1b) \begin{align} {\sf\boldsymbol F}^{(1)} &= \mathcal{E} \left(\delta _0 \frac{\partial \boldsymbol{{u}}_{\parallel }^{(0)}}{\partial Z}\otimes \boldsymbol{{e}}_z + \boldsymbol{{e}}_z \otimes \nabla _{\parallel } u_z^{(0)}\right) +{\textsf{F}}_{zz}^{(1)} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z. \end{align}

Importantly, the leading-order contribution to ${\sf\boldsymbol F}$ is diagonal and can be expressed as ${\sf\boldsymbol F}^{(0)} = \textrm{diag}(1,1, J^{(0)})$ , where the local contraction ratio is given by $J^{(0)} = 1 + \Phi ^{(0)} - \phi _0$ from the incompressibility condition (2.4). Thus, removal of fluid from the pore space drives a uniaxial contraction of the solid matrix along the vertical direction.

The components of the PK1 stress tensor for the film are expanded as ${\sf\boldsymbol S}_{\parallel } = {\sf\boldsymbol S}_{\parallel }^{(0)} + O(\epsilon )$ , ${\sf\boldsymbol S}_{\perp } = {\sf\boldsymbol S}_{\perp }^{(1)} + O(\epsilon )$ , $\textsf{S}_{z \alpha } = \textsf{S}_{z \alpha }^{(1)} + O(\epsilon )$ and $\textsf{S}_{zz} = \textsf{S}_{zz}^{(2)} + O(\epsilon )$ . Non-zero superscripts are used when a component has already been scaled by a power of $\epsilon$ when non-dimensionalising the model. The elastic stress tensor and the pressure are written as $\boldsymbol{\mathsf{\Sigma }} = \boldsymbol{\mathsf{\Sigma }}^{(0)} + \epsilon \boldsymbol{\mathsf{\Sigma }}^{(1)} + O(\epsilon ^2)$ and $p = p^{(0)} + O(\epsilon )$ . From (5.1a) and (4.5), we can deduce that the leading-order contribution to the elastic stress tensor is diagonal and given by $\boldsymbol{\mathsf{\Sigma }}^{(0)} = \boldsymbol{\mathsf{\Sigma }}_{\parallel }^{(0)} + \boldsymbol{\mathsf{\Sigma }}_{zz}^{(0)} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z$ , where

(5.2a) \begin{align} \boldsymbol{\mathsf{\Sigma }}_{\parallel }^{(0)} &= a(\nu _f)(J^{(0)} - 1) J^{(0)} {\sf\boldsymbol I}_{\parallel },\qquad\qquad\ \qquad\end{align}

(5.2b) \begin{align}\boldsymbol{\mathsf{\Sigma }}_{zz}^{(0)} &= a(\nu _f)(J^{(0)} - 1) + b(\nu _f)\left(J^{(0)} - \frac{1}{J^{(0)}}\right). \end{align}

In addition, by considering the $O(1)$ contributions of the PK1 stress tensor (4.4), we see that

(5.3a) \begin{align} {\sf\boldsymbol S}_{\parallel }^{(0)} &= \boldsymbol{\mathsf{\Sigma }}_{\parallel }^{(0)} - p^{(0)} J^{(0)} {\sf\boldsymbol I}_{\parallel },\\[-9pt]\nonumber \end{align}

(5.3b) \begin{align} 0 &=\boldsymbol{\mathsf{\Sigma }}_{zz}^{(0)} - p^{(0)}. \end{align}

Therefore, by substituting (5.2) into (5.3), we can deduce that the pressure and the in-plane PK1 stress are given by

(5.4a) \begin{align} p^{(0)} &=a(\nu _f)(J^{(0)} - 1) + b(\nu _f)\left(J^{(0)} - \frac{1}{J^{(0)}}\right), \end{align}

(5.4b) \begin{align} {\sf\boldsymbol S}_{\parallel }^{(0)} &= b(\nu _f)\left [1 - (J^{(0)})^2\right ]{\sf\boldsymbol I}_{\parallel }.\qquad\qquad\qquad\ \end{align}

Since $J^{(0)} \lt 1$ , the pressure will be negative. In dimensional terms, this means the fluid pressure will be below atmospheric pressure (which has been set to zero). The in-plane elastic stresses will be compressive, while the total (PK1) in-plane stresses will be tensile. For convenience, we define

(5.5) \begin{align} \textsf{S}_{\parallel }(J) \equiv b(\nu _f)\left (1 - J^2\right ) \end{align}

so that the in-plane PK1 stress tensor can be written as ${\sf\boldsymbol S}_{\parallel }^{(0)} = \textsf{S}_{\parallel }(J^{(0)}) {\sf\boldsymbol I}_{\parallel }$ . From the $O(\epsilon )$ contributions to the PK1 stress tensor (4.4), we find, after simplification, that

(5.6a) \begin{align} {\sf\boldsymbol S}_{\perp }^{(1)} = \mathcal{E} b(\nu _f) \left(\delta _0 \frac{\partial \boldsymbol{{u}}_{\parallel }^{(0)}}{\partial Z} + J^{(0)} \nabla _{\parallel } u_z^{(0)}\right), \end{align}

(5.6b) \begin{align} \textsf{S}_{z \alpha }^{(1)} \boldsymbol{{e}}_\alpha = \mathcal{E} b(\nu _f) \left(\delta _0 J^{(0)} \frac{\partial \boldsymbol{{u}}_{\parallel }^{(0)}}{\partial Z} + \nabla _{\parallel }u_z^{(0)}\right). \end{align}

Equation (5.6a) will enable the in-plane displacement to be determined, whereas (5.6b) will enable the vertical traction to be computed via (4.13b).

The shear stress ${\sf\boldsymbol S}_{\perp }^{(1)}$ in the film can be determined by integrating the in-plane stress balance (4.3a) and imposing the stress-free boundary condition (4.10b) to obtain, after simplification,

(5.7) \begin{align} {\sf\boldsymbol S}_{\perp }^{(1)} = \nabla _{\parallel } \int _{Z}^{H_f + \mathcal{H}/ 2}\textsf{S}_{\parallel }(J^{(0)})\, \textrm{d} Z. \end{align}

A differential equation for the vertical component of the film displacement $u_z^{(0)}$ can be obtained by calculating $J^{(0)} = \det {\sf\boldsymbol F}^{(0)}$ using (5.1a) to obtain

(5.8) \begin{align} 1 + \mathcal{E} \frac{\partial u_z^{(0)}}{\partial Z} = J^{(0)}. \end{align}

Having determined ${\sf\boldsymbol S}_{\perp }^{(1)}$ and $\partial u_z^{(0)}/\partial Z$ , we can use (5.6a) to formulate a differential equation for $\boldsymbol{{u}}_{\parallel }^{(0)}$ . After solving (5.6a) and (5.8) and imposing continuity of displacement (4.12), we find that the film displacements can be written as

(5.9a) \begin{align} \boldsymbol{{u}}_{\parallel }^{(0)} &= \bar{\boldsymbol{{u}}}_{\parallel } - (1/ 2)\nabla _{\parallel }w \nonumber \\ & \quad + (\delta _0 \mathcal{E} b(\nu _f))^{-1}\Bigg \{ \int _{\mathcal{H}/2}^{Z}\bigg [{\sf\boldsymbol S}_{\perp }^{(1)} - \mathcal{E} b(\nu _f) \nabla _{\parallel } w - b(\nu _f) J^{(0)} \int _{\mathcal{H}/2}^{Z} \nabla _{\parallel } J^{(0)}\,\textrm{d} Z'\bigg ]\,\textrm{d} Z \Bigg \},\\[-9pt]\nonumber \end{align}

(5.9b) \begin{align} u_z^{(0)} &= \mathcal{E}^{-1} \int _{\mathcal{H}/ 2}^{Z}(J^{(0)} - 1)\,\textrm{d} Z + w. \end{align}

The Eulerian film thickness can be obtained from (5.9b) as

(5.10) \begin{align} h_f = \int _{\mathcal{H}/2}^{H_f + \mathcal{H}/2} J^{(0)}(\boldsymbol{{X}}_{\parallel }, Z, t)\,\textrm{d} Z. \end{align}

We are now in a position to calculate the tractions. Substituting (5.4b) into (4.13a) and using (5.5) leads to an expression for the in-plane traction given by

(5.11a) \begin{align} \boldsymbol{\tau }_{\parallel } = \nabla _{\parallel } \int _{\mathcal{H}/2}^{H_f + \mathcal{H}/ 2}\textsf{S}_{\parallel }(J^{(0)})\, \textrm{d} Z. \end{align}

The vertical traction can be obtained by first substituting the expressions for the displacements (5.9) into (5.6b) and then inserting the result into (4.13b) to find

(5.11b) \begin{align} \tau _{z} &= \nabla _{\parallel }^2 \int _{\mathcal{H}/2}^{H_f + \mathcal{H}/2} J^{(0)} \left (\int _{Z}^{H_f + \mathcal{H}/2}{\textsf{S}_{\parallel }(J^{(0)})}\,\textrm{d} Z'\right )\,\textrm{d} Z \nonumber \\[4pt] &- \nabla _{\parallel } \cdot \int _{\mathcal{H}/ 2}^{H_f + \mathcal{H}/ 2}\left (\int _{Z}^{H_f + \mathcal{H}/2}{\textsf{S}_{\parallel }(J^{(0)})}\, \textrm{d} Z'\right )\nabla _{\parallel } J^{(0)}\,\textrm{d} Z \nonumber \\[4pt] &+ \nabla _{\parallel } \cdot \int _{\mathcal{H}/ 2}^{H_f + \mathcal{H}/ 2}{\textsf{S}_{\parallel }(J^{(0)})} \left (\int _{\mathcal{H}/2}^{Z} \nabla _{\parallel } J^{(0)}\, \textrm{d} Z'\right )\, \textrm{d} Z\ \nonumber \\[4pt] &+ \nabla _{\parallel } \cdot \int _{\mathcal{H}/2}^{H_f + \mathcal{H}/2} \mathcal{E}{\textsf{S}_{\parallel }(J^{(0)})} \nabla _{\parallel } w\, \textrm{d} Z. \end{align}

At this point, the leading-order mechanical problem for the film has been solved in terms $J^{(0)} = 1 + \Phi ^{(0)} - \phi _0$ . Thus, all that remains is to solve the transport problem for the nominal solvent fraction $\Phi ^{(0)}$ .

The expressions for the traction (5.11) can be inserted into the FvK equation (4.8). In the distinguished limit $\delta = O(\epsilon )$ with $\mathcal{E} = O(1)$ as $\epsilon \to 0$ , the full set of FvK equations must be solved; no further simplifications to the plate model are possible. In particular, both the in-plane and vertical tractions drive plate bending despite the differences in their orders of magnitude.

5.1.2. Soft films on thick plates

We now consider the limit in which the film is thin relative to the plate, $\epsilon \ll \delta$ . We further assume that the film is soft so that $\mathcal{E} = O(1)$ . In this case, the deflection of the plate will be proportional to the film thickness. The film displacements are expanded as $\boldsymbol{{u}}_{\parallel } = \boldsymbol{{u}}_{\parallel }^{(0)} + o(1)$ and $u_z = u_z^{(0)} + o(1)$ . The deformation gradient tensor ${\sf\boldsymbol F}$ remains diagonal to leading order; see (4.7). Therefore, many of the results obtained in Section 5.1.1 for the case $\delta = O(\epsilon )$ still apply. In particular, the leading-order contributions to pressure, in-plane stresses, shear stresses and vertical displacement are given by (5.4a), (5.4b), (5.7), and (5.9b), respectively. The main difference occurs in the form of the in-plane displacement, which can be found from examining the $\alpha z$ component of the stress-strain relation (4.4), which reads as

(5.12) \begin{align} \epsilon \textsf{S}_{\alpha z} = \boldsymbol{{e}}_\alpha \boldsymbol{\mathsf{\Sigma }} \boldsymbol{{e}}_z - p J \boldsymbol{{e}}_\alpha {\sf\boldsymbol F}^{-T} \boldsymbol{{e}}_z = \mathcal{E} b(\nu _f)\left (\delta \frac{\partial u_\alpha ^{(0)}}{\partial Z} + \epsilon J^{(0)} \frac{\partial u_z^{(0)}}{\partial X_\alpha } \right ) + o(\epsilon ), \end{align}

where $J^{(0)} = 1 + \mathcal{E} \partial u_z^{(0)}/\partial Z = 1 + \Phi ^{(0)} - \phi _0$ . The leading-order contributions to (5.12), which are $O(\delta )$ in size, imply that $\partial u_\alpha ^{(0)}/\partial Z = 0$ . Solving this equation and imposing continuity of in-plane displacements (4.12a) then shows that $\boldsymbol{{u}}_{\parallel }^{(0)} = \bar{\boldsymbol{{u}}}_{\parallel } - (1/2) \nabla _{\parallel } w$ . Therefore, the in-plane displacement of the film simply matches that of the upper plate surface. In-plane deformation due to volumetric contraction is a higher-order effect.

The leading-order problem for the plate can be obtained by taking $\epsilon \to 0$ with $\epsilon \ll \delta$ and $\mathcal{E} = O(1)$ in the FvK equations given by (4.8). The vertical traction $\tau _{z}$ drops out of (4.8a) and is not required. Thus, plate bending is dominated by the in-plane traction $\boldsymbol{\tau }_{\parallel }$ , which can be obtained from (5.11a). The nonlinear terms from $\nabla _{\parallel } \cdot (\bar{{\sf\boldsymbol S}}_{\parallel } \nabla _{\parallel } w)$ also drop out of (4.8a); thus, the vertical plate displacement satisfies a linear plate equation. The mean in-plane plate displacements $\bar{\boldsymbol{{u}}}_{\parallel }$ can be obtained by solving (4.8b)–(4.8d) after neglecting the nonlinear terms in the strain tensor (4.8d).

5.1.3. Stiff films on thick plates

The final case to consider is when the film is thin relative to the plate, $\epsilon \ll \delta$ , and stiff, so that $\mathcal{E} = O(\delta \epsilon ^{-1})$ . We therefore write the effective film stiffness as $\mathcal{E} = \delta \epsilon ^{-1} \mathcal{E}_1$ . The deformation gradient tensor for the film (4.7) is given by

(5.13) \begin{align} {\sf\boldsymbol F} = {\sf\boldsymbol I} + \mathcal{E}_1 \left (\delta ^2 \nabla _{\parallel } \boldsymbol{{u}}_{\parallel } + \delta ^2 \epsilon ^{-1} \frac{\partial \boldsymbol{{u}}_{\parallel }}{\partial Z} \otimes \boldsymbol{{e}}_z + \delta \boldsymbol{{e}}_z \otimes \nabla _{\parallel } u_z + \delta \epsilon ^{-1}\frac{\partial u_z}{\partial Z} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z\right ). \end{align}

The requirement that ${\sf\boldsymbol F}$ is invertible implies that the large $O(\delta \epsilon ^{-1})$ term in (5.13) must be smaller in magnitude; consequently, the vertical displacement must have the expansion

(5.14) \begin{align} u_z(\boldsymbol{{X}}_{\parallel }, Z, t) = u_z^{(0)}(\boldsymbol{{X}}_{\parallel },t) + \epsilon \delta ^{-1}\tilde{u}_z(\boldsymbol{{X}}_{\parallel }, Z, t) + o(\epsilon \delta ^{-1}). \end{align}

The expansion for the in-plane film displacement can be found using a similar argument to that in Section 5.1.2, leading to

(5.15) \begin{align} \boldsymbol{{u}}_{\parallel }(\boldsymbol{{X}}_{\parallel }, Z, t) = \boldsymbol{{u}}_{\parallel }^{(0)}(\boldsymbol{{X}}_{\parallel },t) + \epsilon \delta ^{-1}\tilde{\boldsymbol{{u}}}_{\parallel }(\boldsymbol{{X}}_{\parallel }, Z, t) + o(\epsilon \delta ^{-1}). \end{align}

By imposing continuity of displacement (4.12), we conclude that $\boldsymbol{{u}}_{\parallel }^{(0)} = \bar{\boldsymbol{{u}}}_{\parallel } - (1/2) \nabla _{\parallel } w$ and $u_z^{(0)} = w$ . Thus, the leading-order film displacement matches the plate displacement, i.e. the film simply bends with the plate. Boundary conditions for $\tilde{\boldsymbol{{u}}}_{\parallel }$ and $\tilde{u}_z$ could be obtained by calculating a higher-order solution to the plate model and imposing continuity of displacement; however, we do not do this here. The deformation gradient tensor for the film therefore reduces to

(5.16) \begin{align} {\sf\boldsymbol F} = {\sf\boldsymbol I} + \mathcal{E}_1 \frac{\partial \tilde{u}_z}{\partial Z} \boldsymbol{{e}}_z \otimes \boldsymbol{{e}}_z + \delta \mathcal{E}_1 \left (\frac{\partial \tilde{\boldsymbol{{u}}}_{\parallel }}{\partial Z} \otimes \boldsymbol{{e}}_z + \boldsymbol{{e}}_z \otimes \nabla _{\parallel } w\right ) + O(\epsilon ). \end{align}

Although ${\sf\boldsymbol F}$ is diagonal to leading order, the $O(\delta )$ corrections result in substantially different film mechanics compared to when the film is soft ( $\mathcal{E} = O(1)$ ). To see this, the deformation gradient tensor (5.16) can be substituted into the stress-strain relation (4.4). Collecting the $O(\delta )$ contributions to the $\alpha z$ and $z \alpha$ components leads to the following equations, respectively,

(5.17a) \begin{align} \frac{\partial \tilde{\boldsymbol{{u}}}_{\parallel }}{\partial Z} + J^{(0)} \nabla _{\parallel } w &= \textbf{0},\\[-9pt]\nonumber \end{align}

(5.17b) \begin{align} J^{(0)}\frac{\partial \tilde{\boldsymbol{{u}}}_{\parallel }}{\partial Z} + \nabla _{\parallel } w &= \textbf{0}, \end{align}

where $J^{(0)} = 1 + \mathcal{E}_1 \partial \tilde{u}_z/\partial Z$ . Solving (5.17) shows that $\nabla _{\parallel } w = \textbf{0}$ , which implies that the plate remains flat. This unphysical result stems from the assumption that the shear stresses $\textsf{S}_{\alpha z}$ and $\textsf{S}_{z \alpha }$ in the film have the same order of magnitude, which would be the case in linear elasticity. Thus, a resolution can be found by elevating the asymptotic order of $\textsf{S}_{z \alpha }$ and writing $\textsf{S}_{z \alpha } = \delta \epsilon ^{-1} \tilde{\textsf{S}}_{z \alpha }$ . From the vertical stress balance (4.3b), we must also elevate the order of $\textsf{S}_{zz}$ by writing $\textsf{S}_{zz} = \delta \epsilon ^{-1} \tilde{\textsf{S}}_{zz}$ . Importantly, the elevated order of the vertical stress implies that the non-dimensional vertical traction must also be rescaled according to $\tau _{z} = \delta \epsilon ^{-1} \tilde{\tau }_{z}$ . In terms of dimensional quantities, these rescalings imply that $\textsf{S}_{z \alpha } \sim \delta E_f$ , $\textsf{S}_{zz} \sim \epsilon \delta E_f$ , and $\tau _{z} \sim \epsilon \delta E_f$ .

With these rescalings, we proceed by expanding ${\sf\boldsymbol S}_{\perp } = {\sf\boldsymbol S}_{\perp }^{(1)} + o(1)$ and $\tilde{\textsf{S}}_{z \alpha }$ as $\tilde{\textsf{S}}_{z \alpha } = \tilde{\textsf{S}}_{z \alpha }^{(1,1)} + o(1)$ . The $\alpha z$ component of the stress-strain relation still leads to (5.17a); however, the $z \alpha$ component now gives

(5.18) \begin{align} \tilde{\textsf{S}}_{z \alpha }^{(1,1)} \boldsymbol{{e}}_\alpha = \mathcal{E}_1 b(\nu _f) \left (J^{(0)}\frac{\partial \tilde{\boldsymbol{{u}}}_{\parallel }}{\partial Z} + \nabla _{\parallel } w \right ) = \mathcal{E}_1 \textsf{S}_{\parallel }(J^{(0)}) \nabla _{\parallel } w, \end{align}

which replaces (5.17b). The second equality arises from using (5.17a) and the definition of $\textsf{S}_{\parallel }$ in (5.5). The vertical traction can now be obtained by substituting (5.18) into (4.13b) to find that

(5.19) \begin{align} \tilde{\tau }_{z} = \nabla _{\parallel } \cdot \int _{\mathcal{H}/ 2}^{H_f + \mathcal{H}/2} \mathcal{E}_1 \textsf{S}_{\parallel }(J^{(0)}) \nabla _{\parallel } w\,\textrm{d} Z. \end{align}

Despite the elevated asymptotic order of the vertical stress $\textsf{S}_{zz}$ , the $zz$ component of the stress-strain relation (4.4) leads to an expression for the pressure that is identical to (5.4a). Consequently, the in-plane stresses and traction are given by (5.4b) and (5.11a), respectively.

In the regime of stiff films on thick plates, $\epsilon \ll \delta$ with $\mathcal{E} = O(\delta \epsilon ^{-1})$ , the full set of FvK equations given by (4.8) must be solved. The small $O(\epsilon \delta ^{-1})$ prefactor of $\tau _{z}$ in (4.8a) is exactly offset by the large $O(\delta \epsilon ^{-1})$ magnitude of $\tau _{z}$ . Consequently, the vertical traction must be considered despite the film being much thinner than the plate.

5.2. Film drying

Having solved the mechanical problem in the film, the fluid-transport problem is asymptotically reduced by taking the limit $\epsilon \to 0$ . The goal is to derive a hierarchy of evolution equations for the nominal solvent fraction that capture different drying regimes depending on the magnitude of the Péclet number, $\textrm{Pe}$ .

Physically, the different drying regimes arise from differences in the relative rates at which fluid is removed from the pore space by evaporation and replenished by bulk transport. When $\textrm{Pe} = O(1)$ , the time scales of evaporation and in-plane fluid transport are commensurate. If $\textrm{Pe} \ll 1$ , evaporation is slow and a homogeneous drying process can be expected. When $\textrm{Pe} \gg 1$ , evaporation is fast relative to in-plane transport and a highly non-uniform film composition is likely to arise.

In reducing the fluid-transport problem, we work within the distinguished limit of $\delta = O(\epsilon )$ as $\epsilon \to 0$ by assuming the film and plate have comparable thicknesses. In this case, the deformation gradient tensor has the asymptotic expansion ${\sf\boldsymbol F} = \textrm{diag}(1, 1, J) + O(\epsilon )$ . Thus, the off-diagonal components of the permeability tensor are small, specifically, ${\sf\boldsymbol {K}}_{\perp } = O(\epsilon )$ . We therefore expand the pressure, nominal fluid fraction, and components of the permeability tensor as $p = p^{(0)} + O(\epsilon )$ , $\Phi = \Phi ^{(0)} + O(\epsilon )$ , ${\sf\boldsymbol {K}}_{\parallel } = {\sf\boldsymbol {K}}_{\parallel }^{(0)} + O(\epsilon )$ , ${\sf\boldsymbol {K}}_{\perp } = \epsilon {\sf\boldsymbol {K}}_{\perp }^{(1)} + O(\epsilon ^2)$ , and ${\textsf{K}}_{zz} ={\textsf{K}}_{zz}^{(0)} + O(\epsilon )$ . The expansions of the flux will depend on the size of $\textrm{Pe}$ .

Although the reduction of the fluid-transport problem will be carried out in the limit that $\delta = O(\epsilon )$ , the results also apply when the film is thin relative to the plate and $\delta \gg \epsilon$ .

5.2.1. Slow evaporation: $\textrm{Pe} = O(1)$

We first consider the drying problem when the evaporation is slow, as characterised by $\textrm{Pe} = O(1)$ as $\epsilon \to 0$ . The components of the flux are expanded as $\boldsymbol{{Q}}_\parallel = \boldsymbol{{Q}}_\parallel ^{(0)} + O(\epsilon )$ and $Q_z = Q_z^{(0)} + O(\epsilon )$ . Collecting the $O(\epsilon ^{-2})$ terms in (4.2b) lead to the conclusion that the pressure $p^{(0)}$ , and hence $J^{(0)}$ and $\Phi ^{(0)}$ , are independent of $Z$ . Therefore, the film composition is uniform along the vertical direction. The $O(1)$ terms in (4.2a) provide an expression for the in-plane flux, $\boldsymbol{{Q}}_\parallel ^{(0)} = -{\sf\boldsymbol {K}}_{\parallel }^{(0)} \nabla _{\parallel } p^{(0)}$ , where ${\sf\boldsymbol {K}}_{\parallel }^{(0)} = k(\phi ^{(0)}) J^{(0)} {\sf\boldsymbol I}_{\parallel }$ . The Eulerian fluid fraction is given by $\phi ^{(0)} = \Phi ^{(0)}/ J^{(0)}$ .

By integrating the solvent conservation equation (4.1) from $Z = \mathcal{H}/2$ to $Z = H_f(\boldsymbol{{X}}_{\parallel }) + \mathcal{H}/2$ and using the boundary conditions (4.10a) and (4.11), we find that the nominal solvent fraction evolves according to

(5.20) \begin{align} H_f \frac{\partial \Phi ^{(0)}}{\partial t} = \nabla _{\parallel } \cdot \left (H_f J^{(0)} k(\phi ^{(0)}) \nabla _{\parallel } p^{(0)}\right ) - \textrm{Pe}\,\mathcal{V}(\phi ^{(0)}), \end{align}

where $J^{(0)} = 1 + \Phi ^{(0)} - \phi _0$ . Upon substituting the expression for the pressure $p^{(0)}$ given by (5.4a) into (5.20), a nonlinear diffusion equation for $\Phi ^{(0)}$ is obtained. Equation (5.20) can be solved with the boundary conditions $\boldsymbol{{Q}}_\parallel ^{(0)} \cdot \boldsymbol{{N}} = 0$ at the contact line, $X_1=0$ , $X_1 = 1$ and $X_2 = \pm \mathcal{W}/2$ , along with the initial condition $\Phi ^{(0)} = \phi _0$ .

Upon solving for $\Phi ^{(0)}$ and computing $J^{(0)}$ , the Eulerian film thickness can be obtained from (5.10) as $h_f = J^{(0)} H_f$ . Moreover, the traction can be determined from (5.11) as

(5.21a) \begin{align} \boldsymbol{\tau }_{\parallel } &= \nabla _{\parallel } (H_f \textsf{S}_{\parallel }),\qquad\qquad\qquad\qquad\qquad\ \ \end{align}

(5.21b) \begin{align} \tau _{z} &= \frac{1}{2}\nabla _{\parallel }^2 \left(J^{(0)} H_f^2 \textsf{S}_{\parallel }\right) + \mathcal{E} \nabla _{\parallel } \cdot \left(H_f \textsf{S}_{\parallel } \nabla _{\parallel } w\right), \end{align}

where $\textsf{S}_{\parallel } = \textsf{S}_{\parallel }(J^{(0)})$ is given by (5.5). In the thin-film regime $\delta \gg \epsilon$ with $\mathcal{E} = O(1)$ , the in-plane traction is given by (5.21a) and the vertical traction (5.21b) is not required. When $\delta \gg \epsilon$ with $\mathcal{E} = O(\delta \epsilon ^{-1})$ , the in-plane traction is given by (5.21a) and the vertical traction is

(5.22) \begin{align} \tau _{z} = \mathcal{E} \nabla _{\parallel } \cdot \left (H_f \textsf{S}_{\parallel } \nabla _{\parallel } w\right ). \end{align}

If evaporation is very slow, $\textrm{Pe} \ll 1$ , then time can be rescaled as $t \sim \textrm{Pe}^{-1}$ . The leading-order (in $\textrm{Pe}$ ) part of (5.20) then implies that $p^{(0)}$ and hence $\Phi ^{(0)}$ are uniform in space. Integrating (5.20) over the plate surface leads to a simple differential equation for the nominal volume fraction, which can be expressed in terms of the original non-dimensional time as

(5.23) \begin{align} \frac{\textrm{d} \Phi ^{(0)}}{\textrm{d} t} = -\frac{\textrm{Pe} \mathcal{W} \mathcal{V}(\phi ^{(0)})}{V_0}, \end{align}

where $V_0$ is the initial volume of the film. By taking $J^{(0)}$ and hence $\textsf{S}_{\parallel }$ to be constant in (5.21a), the in-plane traction reduces to $\boldsymbol{\tau }_{\parallel } = \textsf{S}_{\parallel } \nabla _{\parallel } H_f$ . Given that $\textsf{S}_{\parallel } \gt 0$ , the direction of the traction follows the gradient of the film thickness. Near the contact line, where $H_f$ tends to zero, the film will always pull the plate inwards, that is, towards the centre of the plate.

5.2.2. Moderate evaporation: $\textrm{Pe} = O(\epsilon ^{-1})$

The limit of moderate evaporation is characterised by $\textrm{Pe} = O(\epsilon ^{-1})$ . Thus, we write $\textrm{Pe} = \epsilon ^{-1} \textrm{Pe}_0$ . In order to obtain a balance in the boundary condition (4.10a), we must rescale the vertical component of the flux as $Q_z = \epsilon ^{-1} \tilde{Q}_z$ . Both flux components now have the same order of magnitude. Despite the elevated asymptotic order of $Q_z$ , the vertical component of Darcy’s law (4.2b) implies that $p^{(0)}$ , $J^{(0)}$ and $\Phi ^{(0)}$ are still independent of $Z$ . Balancing terms in the conservation law for the nominal fluid fraction (4.1) requires rescaling time as $t = \epsilon \tilde{t}$ . Integrating the $O(\epsilon ^{-1})$ terms in the conservation law across the film thickness leads to

(5.24) \begin{align} H_f \frac{\partial \Phi ^{(0)}}{\partial \tilde{t}} = - \textrm{Pe}_0 \,\mathcal{V}(\phi ^{(0)}). \end{align}

Equation (5.24) can be treated as an ordinary differential equation for $\Phi ^{(0)}$ , with the in-plane coordinates $\boldsymbol{{X}}_{\parallel }$ acting as parameters, and solved with the initial condition $\Phi ^{(0)} = \phi _0$ at $\tilde{t} = 0$ . The expressions for the traction that are provided in Section 5.2.1 also apply to this drying regime.

The reduced transport problem defined by (5.24) does not capture the in-plane components of the fluid flux. However, near the contact line, where $H_f$ is small, (5.24) predicts that the nominal fluid fraction will undergo a rapid decrease. Consequently, large in-plane gradients in the fluid fraction are expected near the contact line, and these will locally amplify the in-plane flux. From an asymptotics perspective, near the contact line, the order of the in-plane flux must be elevated, resulting in additional terms that must be considered in the reduced model. The dynamics near the contact line can therefore be resolved by carrying out a boundary-layer analysis using matched asymptotic expansions. However, we do not resolve these boundary layers because they are unimportant in terms of plate bending; see Section 6.2 for more details.

5.2.3. Fast evaporation: $\textrm{Pe} = O(\epsilon ^{-2})$

The limit of fast evaporation occurs when $\textrm{Pe} = O(\epsilon ^{-2})$ . Thus, we write $\textrm{Pe} = \epsilon ^{-2} \textrm{Pe}_1$ . By following the same strategy as in Section 5.2.2, we deduce that capturing the dynamics requires rescaling time as $t = \epsilon ^2 \hat{t}$ and the vertical component of the flux as $Q_z = \epsilon ^{-2} \hat{Q}_z$ . By writing $\hat{Q}_z = \hat{Q}_z^{({-}2)} + O(\epsilon )$ and equating the $O(\epsilon ^{-2})$ components of (4.2b), we find that $\hat{Q}_z^{({-}2)} = -{\textsf{K}}_{zz}^{(0)} \partial _Z p^{(0)}$ , with ${\textsf{K}}_{zz}^{(0)} = k(\phi ^{(0)})/ J^{(0)}$ . Therefore, the rate of evaporation is now sufficiently large that vertical gradients will occur in the pressure $p^{(0)}$ , the volumetric contraction ratio $J^{(0)}$ , and the nominal fluid fraction $\Phi ^{(0)}$ . Equating the $O(\epsilon ^{-2})$ contributions of the conservation law (4.1) leads to a nonlinear diffusion equation for the fluid fraction given by

(5.25) \begin{align} \frac{\partial \Phi ^{(0)}}{\partial \hat{t}} + \frac{\partial \hat{Q}_z^{({-}2)}}{\partial Z} = 0. \end{align}

The boundary conditions and initial for (5.25) are

(5.26a) \begin{align} \hat{Q}_z^{({-}2)} &= 0, &\quad Z &= \mathcal{H}/2;\\[-9pt]\nonumber \end{align}

(5.26b) \begin{align} \hat{Q}_z^{({-}2)} &= \textrm{Pe}_1\,\mathcal{V}(\phi ^{(0)}), &\quad Z &= H_f + \mathcal{H}/2;\\[-9pt]\nonumber \end{align}

(5.26c) \begin{align} \Phi ^{(0)} &= \phi _0, &\quad \hat{t} &= 0. \end{align}

Again, the in-plane coordinates $\boldsymbol{{X}}_{\parallel }$ play the role of parameters. After solving for $\Phi ^{(0)}$ and hence $J^{(0)}$ , the tractions can be evaluated by computing the integrals in (5.11a) and (5.11b) or (5.19).

The reduced transport problem (5.25) does not capture in-plane fluid transport. However, in contrast to the case $\textrm{Pe} = O(\epsilon ^{-1})$ , the contributions from the in-plane fluid flux are always asymptotically smaller than those from the large vertical flux. Therefore, there is no need to resolve the dynamics near the contact line using a separate boundary-layer analysis.

6. Validation and parametric studies

The asymptotic reductions are validated by comparing them against finite element (FE) solutions of the fully coupled film-plate model proposed in Section 2. The FE method is implemented in Python with FEniCS [1, 28]. The multiphenics [2] package is used to couple the film and plate problems. The traction $\boldsymbol{\tau }$ is treated as a Lagrange multiplier that enables the continuity of displacement at the film-plate interface to be imposed. The film and plate displacements are represented using P2 elements. The pressure and nominal fluid fraction are represented using P1 elements. Finally, DGT1 elements are used for the traction. Time is discretised using a fully implicit Euler method with a fixed time step. The film and plate problems are simultaneously solved as a monolithic nonlinear system. We employ the PETSc SNES nonlinear solver, using MUMPS to solve the linear systems during each Newton iteration.

The FE simulations are carried out in a two-dimensional Cartesian geometry under the assumption of plane strain ( $\mathcal{W} \gg 1$ ). In this case, the plate reduces to a beam and the FvK equations (4.8) simplify considerably. The case of plane stress, valid for plates with small widths ( $\mathcal{W} \ll 1$ ), can be considered by replacing $\mathcal{B}$ with $(1 - \nu _p^2) \mathcal{B}$ in the expressions below. The film is assumed to have a parabolic initial profile that can be written in dimensional form as $H_f(X) = 4 \epsilon X (1 - X/ L)$ , where we define $X \equiv X_1$ for simplicity. The plate is assumed to be made of steel with a Poisson’s ratio of $\nu _p = 0.27$ . The Poisson’s ratio of the film is set to $\nu _f = 0.2$ , corresponding to a porous matrix formed from a dried colloidal dispersion [37].

Snapshots from a typical FE simulation are shown in Figure 2. The Eulerian fluid fraction in the film $\phi = \Phi/ J$ is plotted as a heat map. It has been assumed that the film and plate have the same thickness, $\epsilon = \delta = 0.1$ . The Péclet number has been set to $\textrm{Pe} = 100$ . The effective film stiffness has been set to $\mathcal{E} = 3$ . The other parameters are discussed in Section 6.2. Due to Péclet number being $O(\epsilon ^{-2})$ in size, film drying is highly non-uniform. Drying first occurs at the contact line, where the film is the thinnest, and then proceeds inwards. Due to the high rate of evaporation, fluid that is removed from the pore space near the free surface cannot be replenished by fluid in the bulk; hence, vertical gradients in the composition arise. Despite drying initially being confined to the contact line, the plate undergoes an appreciable deflection for small times. The deflection then monotonically increases as the film approaches its homogeneous steady state.

Figure 2. Finite element simulations of the film-plate model presented in Sec. 2. The heat map represents the Eulerian fluid fraction (porosity) of the film, $\phi = \Phi / J$ . Simulation snapshots are shown at times (a) $t = 0$ , (b) $t = 0.005$ , (c) $t = 0.010$ , (d) $t = 0.015$ , (e) $t = 0.020$ , and (f) $t = 0.025$ . The parameter values are $\epsilon = \delta = 0.1$ , $\textrm{Pe} = 100$ , $\varepsilon = 3$ with $\phi_0 = 0.68$ , and $k(\phi) \equiv 1$ . The non-dimensional evaporative flux is $\mathcal{V}(\phi) = \phi - \phi_\infty$ , where $\phi_\infty = 0.36$ .

6.1. Comparison of the steady states

Calculating the steady states of the film-plate system enables the accuracy of the asymptotic solutions of the mechanical problem to be established. In particular, the steady states satisfy a purely mechanical problem and can be obtained without solving the fluid-transport problem. Thus, in this section, it is assumed that the film has dried to a homogeneous steady state with uniform contraction ratio $J \lt 1$ . In this case, the poroelastic model for the film developed in Section 2.3 reduces to the equilibrium equations of incompressible nonlinear elasticity, $\nabla \cdot {\sf\boldsymbol S} = \textbf{0}$ with $\det {\sf\boldsymbol F} = J$ , where $J$ is now treated as a prescribed constant.

When the full problem is posed in two Cartesian dimensions under the assumption of plane strain, the FvK equations (4.8) become one-dimensional. Moreover, when $J$ is uniform in the vertical direction so that the expressions for the tractions in (5.21) apply, then the mean in-plane plate stress can be obtained from (4.8b) as $\bar{\textsf{S}}_{11} = -H_f \textsf{S}_{\parallel }(J)$ , where $\textsf{S}_{\parallel }(J)$ is given by (5.5). Inserting this result along with the tractions (5.21) into (4.8a) leads to cancellations that allow two integrations to be performed. After imposing the boundary conditions at the free edge $X = 1$ , the equation for the vertical plate displacement (4.8a) simplifies to

(6.1) \begin{align} -\mathcal{B} \frac{\textrm{d}^2 w}{\textrm{d} X^2} = -\frac{1}{2} H_f \textsf{S}_{\parallel }(J) - \frac{1}{2}\epsilon \delta ^{-1} J H_f^2 \textsf{S}_{\parallel }(J). \end{align}

The boundary conditions for (6.1) are $w(0) = 0$ and $\partial _X w(0) = 0$ . The first term on the right-hand side of (6.1) captures the influence of the in-plane traction, whereas the second term captures the vertical traction. When $J$ is a constant, (6.1) can be solved to find

(6.2) \begin{align} w(X) = \frac{1 }{30 \mathcal{B}}\, \textsf{S}_{\parallel }(J) X^3 \left [10 - 5X + 4 \epsilon \delta ^{-1} J (2X^2 - 6X + 5) X \right ]. \end{align}

Equation (6.1) and its solution (6.2) apply to all of the mechanical regimes considered in Section 5.1; however, for asymptotic consistency, the terms that are proportional to $\epsilon \delta ^{-1}$ should be dropped when $\epsilon \ll \delta$ .

We first compare asymptotic and FE solutions when the film and plate have similar thicknesses. In particular, we set $\delta = \epsilon$ and $\mathcal{E} = 1$ . A comparison of the deflection of the end of the plate, $w(X = 1)$ , shows that the asymptotic and FE solutions are in good agreement across a wide range of $J$ values, especially when $\epsilon = \delta \leq 0.05$ ; see Figure 3(a). The deflection is seen to increase as $J$ is decreased from one, corresponding to greater volumetric contraction of the film. However, the deflection reaches a maximum at $J \simeq 0.30$ , after which it decreases with further decreases in $J$ . The non-monotonic behaviour results from a competition between the in-plane stress in the film $\textsf{S}_{\parallel }$ and the contribution to (6.1) from the vertical traction. The former increases as $J$ decreases, whereas the latter decreases.

Figure 3. Steady-state solutions of the vertical plate displacement when the film and plate have similar thicknesses, $\delta = O(\epsilon)$ with $\mathcal{E} = O(1)$ . Asymptotic solutions are shown as lines and obtained from (6.2); finite element solutions are shown as symbols. (a) The deflection at the end of the plate as a function of the imposed (uniform) contraction ratio J in the film. (b) The vertical displacement of the plate as a function of space when $J = 0.2$ .

Profiles of the vertical plate displacement are found to compare favourably in Figure 3(b). As expected, the plate deflection monotonically increases as the distance from the wall increases.

To validate the asymptotic reduction of the mechanical problem when the film is thin relative to the plate, we set $\epsilon = 0.01$ and $\delta = \epsilon ^{1/2} = 0.1$ . Solutions are computed when $\mathcal{E} = 1$ and $\mathcal{E} = \delta \epsilon ^{-1} = 10$ across a range of $J$ values. By computing the deflection at the end of the plate, we again find excellent agreement between the asymptotic and FE solutions; see Figure 4(a). The strong agreement when $\mathcal{E} = 10$ is remarkable given that the film undergoes extremely large deformations, with the vertical displacement being much greater than the film thickness.

Figure 4. Steady-state solutions when the film is thin relative to the plate, $\epsilon \ll \delta$ . Asymptotic solutions are shown as lines; finite element solutions are shown as symbols. (a) The deflection at the end of the plate as a function of the imposed (uniform) contraction ratio J in the film. The asymptotic solution for w is obtained from (6.2) with $\epsilon \delta^{-1} = 0$ . (b) The vertical traction $\tau_z$ as a function of space. The asymptotic solutions for the traction are given by (5.21b) when $\mathcal{E} = 1$ and (5.22) when $\mathcal{E} = 10$ .

A distinguishing feature of the thin-film regime ( $\epsilon \ll \delta$ ) is that the asymptotic order of the vertical traction increases when the non-dimensional stiffness $\mathcal{E}$ increases. To observe this increase, the vertical traction is computed from the FE and asymptotic solutions by fixing $\epsilon = 0.01$ , $\delta = 0.1$ , and $J = 0.2$ while considering $\mathcal{E} = 1$ and $\mathcal{E} = 10$ . When $J$ is independent of $Z$ , the asymptotic solution for the vertical traction is obtained from (5.21b) when $\mathcal{E} = 1$ and (5.22) when $\mathcal{E} = 10$ . The solutions for the traction are plotted as a function of space in Figure 4(b). There is clearly a large difference in the magnitude of the traction in the two cases, with the traction remaining $O(1)$ in size when $\mathcal{E} = 1$ and becoming $O(\delta \epsilon ^{-1})$ in size when $\mathcal{E} = 10$ . Again, excellent agreement between the asymptotic and FE solutions is found.

6.2. Dynamic simulations

Time-dependent simulations are used to explore the dynamics of drying and bending. The non-dimensional evaporative flux is written in terms of a simple phenomenological law given by $\mathcal{V}(\phi ) = \phi - \phi _\infty$ . The parameter $\phi _\infty$ can have different interpretations depending on the modelling context. For instance, $\phi _\infty$ can represent a constant amount of fluid vapour in the surrounding environment. In any case, with this model for $\mathcal{V}$ , the film will dry until its porosity (fluid fraction) reaches $\phi _\infty$ . For simplicity, we set the scalar permeability of the film, $k(\phi )$ , to be a constant. Hennessy et al. [17] showed that a porosity-dependent permeability can lead to drying fronts propagating into the bulk of the film from the contact line. These fronts separate wet and completely dry solid and hence require the porosity to become very small. By taking $\phi _\infty$ to be sufficiently large, the formation of drying fronts will be suppressed and the role of a porosity-dependent permeability will be negligible.

Our time-dependent simulations are based on the work of Bouchaudy and Salmon [4], who measured drying-induced stresses in poroelastic discs formed from nanoparticle suspensions. Gelation was estimated to occur at a solid fraction of 0.32. Hence, we take $\phi _0 = 1 - 0.32 = 0.68$ . We imagine that $1 - \phi _\infty$ represents the closest packing fraction of the colloidal dispersion. For a poroelastic film composed of monodisperse nanoparticles, $\phi _\infty \simeq 1 - 0.64 = 0.36$ . The drying that occurs after the closest packing fraction is reached can lead to pore invasion by air [27], a phenomenon that our model does not capture. For this set of parameters, the film will shed half of its volume during drying, which can be seen by calculating $J = (1 - \phi _0)/ (1 - \phi _\infty ) = 0.5$ .

6.2.1. Films drying on plates with similar thickness

The first set of time-dependent simulations captures the dynamics when the film and plate have similar thicknesses, $\delta = O(\epsilon )$ with $\mathcal{E} = O(1)$ . The evaporation is assumed to be slow, $\textrm{Pe} = O(1)$ , or moderate, $\textrm{Pe} = O(\epsilon ^{-1})$ . We set $\delta = \epsilon = 0.01$ with $\mathcal{E} = 1$ , and consider $\textrm{Pe} = 1$ or $\textrm{Pe} = 100$ . For both of these cases, the asymptotic solution to the plate displacement can be obtained from (6.1). The asymptotic solutions for the in-plane film stress, $\textsf{S}_{\parallel }$ , and local contraction ratio, $J$ , can be obtained by numerically solving (5.20) when $\textrm{Pe} = O(1)$ and (5.24) when $\textrm{Pe} = O(\epsilon ^{-1})$ .

When evaporation is slow, $\textrm{Pe} = 1$ , the Eulerian fluid fraction (porosity) of the film remains relatively uniform during drying. That is, the in-plane composition gradients, in addition to the vertical composition gradients, are weak. The weak in-plane gradients can be seen in Figure 5(a), where the fluid fraction at the bottom of the film ( $Z = 0$ ) is plotted as a function of the in-plane coordinate $X$ at various times. Black lines represent numerical solutions of the asymptotically reduced transport problem defined by (5.20); circles represent the FE solution. The uniform loss of fluid across the film induces a homogeneous contraction of the film that predominately occurs in the vertical direction. The in-plane stresses in the film, $\textsf{S}_{\parallel }$ , will also be approximately homogeneous. Since the initial film profile is quadratic, the in-plane traction is expected to vary linearly in space, which can be deduced from (5.21a) and is shown in Figure 5(b). The evolution of the vertical traction is more complicated, as seen from Figure 5(c). The vertical traction is typically positive near the contact line, implying the film pulls upwards on the plate. However, near the centre of the plate, the traction becomes negative, implying the film pushes downwards. The vertical displacement of the plate is shown in Figure 5(d); the inset depicts the time evolution of the deflection at the end of the plate, $w(1,t)$ . The plate deflection monotonically increases in space and time. The deflection at the end of the beam initially grows linearly with time and then saturates as the film approaches its steady state. The initial linear increase in the beam deflection is a result of the evaporation rate being approximately constant for small times, with an evaporative flux given by $\mathcal{V} \sim \phi _0 - \phi _\infty$ . In all of the panels of Figure 5, excellent agreement between the asymptotic and FE solutions can be seen.

Figure 5. Drying dynamics for small evaporation rates, $\textrm{Pe} = O(1)$ , when the film and plate have similar thicknesses, $\delta = O(\epsilon)$ . Asymptotic solutions are shown as lines; finite element solutions are shown as circles. (a) The Eulerian fluid fraction at the bottom of the film. (b) and (c) The in-plane and vertical traction. (d) The vertical plate displacement. The inset shows the evolution of the deflection of the end of the plate. Solutions are shown at times $t = 0.1$ , 0.4, 0.7, 1, and 3. The arrows show the direction of increasing time. The parameter values are $\epsilon = \delta = 0.01$ , $\textrm{Pe} = 1$ , and $\mathcal{E} = 1$ .

Increasing the Péclet number to $\textrm{Pe} = 100$ , corresponding to moderate evaporation, leads to a non-uniform drying process with significant in-plane gradients; see Figure 6. The film rapidly dries near the contact line, as seen in Figure 6(a), which leads to a localisation of both the in-plane and vertical components of the traction; see Figure 6(b) and (c), respectively. Despite the strongly non-uniform tractions, the evolution of the plate deflection, shown in Figure 6(d), is remarkably similar to the case when $\textrm{Pe} = O(1)$ . A more in-depth comparison of the plate deflection for different evaporation rates will be provided in Section 6.3. The asymptotic solutions (lines) are in good agreement with the FE simulations (circles). The largest discrepancy occurs in the fluid fraction near the contact line; see Figure 6(a). The asymptotically reduced model (5.24) does not capture the influence of in-plane pressure gradients, which provide a mechanism to replenish the fluid that evaporates at the contact line. Hence, the asymptotically reduced model leads to faster drying compared to the FE solutions. Despite the asymptotically reduced model not correctly capturing the dynamics near the contact line, it is still able to provide a highly accurate prediction of the plate deflection; see Figure 6(d).

Figure 6. Drying dynamics for moderate evaporation rates, $\textrm{Pe} = O(\epsilon^{-1})$ , when the film and plate have similar thicknesses, $\delta = O(\epsilon)$ . Asymptotic solutions are shown as lines; finite element solutions are shown as circles. (a) The Eulerian fluid fraction at the bottom of the film. (b) and (c) The in-plane and vertical traction. (d) The vertical plate displacement. The inset shows the evolution of the deflection of the end of the plate. Solutions are shown at times $t = 0.001$ , 0.005, 0.01, 0.02, and 0.05. The arrows show the direction of increasing time. The parameter values are $\epsilon = \delta = 0.01$ , $\textrm{Pe} = 100$ , and $\mathcal{E} = 1$ .

6.2.2. Films drying on thick plates

The final comparison that we make between the asymptotic and FE solutions considers the case of films that are thin relative to the plate, $\epsilon \ll \delta$ . Moreover, the film is assumed to be soft, $\mathcal{E} = O(1)$ , and evaporation is fast, $\textrm{Pe} = O(\epsilon ^{-2})$ . In particular, we set $\epsilon = 10^{-2}$ , $\delta = \epsilon ^{1/2} = 10^{-1}$ , $\mathcal{E} = 1$ , and $\textrm{Pe} = 5\,\epsilon ^{-2} = 5\times 10^{4}$ . In this regime, the vertical displacement of the plate can be obtained by solving

(6.3) \begin{align} -\mathcal{B} \frac{\partial ^2 w}{\partial X^2} = -\frac{1}{2} H_f \langle \textsf{S}_{\parallel }(J)\rangle, \quad \langle \textsf{S}_{\parallel }(J)\rangle = \frac{1}{H_f} \int _{\mathcal{H}/2}^{H_f + \mathcal{H}/2} \textsf{S}_{\parallel }(J)\, \textrm{d} Z, \end{align}

subject to $w(0,t) = 0$ and $\partial _X w(0,t) = 0$ , where $\langle \textsf{S}_{\parallel } \rangle$ denotes the vertically averaged in-plane film stress. The contraction ratio $J$ and the in-plane film stress $\textsf{S}_{\parallel }$ can be obtained by numerically solving (5.25).

The evolution of the fluid fraction at the bottom of the film and the plate deflection are shown in Figures 7(a) and (b), respectively. Qualitatively, the solutions appear very similar to the $\textrm{Pe} = O(\epsilon ^{-1})$ case shown in Figure 6. There is a rapid depletion of fluid near the contact line, which leads to large in-plane gradients. However, a key difference is that vertical composition gradients occur when $\textrm{Pe} = O(\epsilon ^{-2})$ , which can be seen in Figure 7(c). As expected, the fluid fraction is the greatest at the bottom of the film and the smallest at the free surface. Consequently, the in-plane film stress $\textsf{S}_{\parallel }$ monotonically increases with $Z$ , as shown in Figure 7(d).

Figure 7. Drying dynamics for large evaporation rates, $\textrm{Pe} = O(\epsilon^{-2})$ , when the film is much thinner than the plate, $\epsilon \ll \delta$ . Asymptotic solutions are shown as lines; finite element solutions are shown as circles. (a) The Eulerian fluid fraction at the bottom of the film. (b) The vertical plate displacement. The inset shows the evolution of the deflection of the end of the plate. (c) and (d) The Eulerian fluid fraction and in-plane film stress at the centre of the film. Solutions are shown at times $t = 2 \times 10^{-6}$ , $10^{-5}$ , $2 \times 10^{-5}$ , $4 \times 10^{-5}$ , and $10^{-4}$ . The arrows show the direction of increasing time. The parameter values are $\epsilon = 0.01$ , $\delta = 0.1$ , $\textrm{Pe} = 5 \times 10^{4}$ , and $\mathcal{E} = 1$ .

The asymptotically reduced model for the fluid-transport problem (5.25) provides an excellent approximation to the FE solutions of the nominal fluid fraction and the in-plane film stress; see Figures 7(a), (c), and (d). In particular, the fluid fraction near the contact line is well captured. The reduced model for the plate deflection provides a close approximation to the FE solution, as shown in Figure 7(b). Since the reduced model does not account for the vertical traction when $\epsilon \ll \delta$ and $\mathcal{E} = O(1)$ , the deflection is slightly smaller than that obtained using the FE method.

We have also compared the asymptotic and FE solutions when $\mathcal{E}$ is increased to $\mathcal{E} = \delta \epsilon ^{-1} = 10$ , keeping the other parameters the same (not shown). The solutions are virtually the same as those in Figure 7. The agreement between the fluid fraction and in-plane stresses remains excellent. However, there is much stronger agreement in the plate deflection. This is because the asymptotically reduced model for the plate now captures terms associated with the in-plane (plate) stress and the vertical traction. However, these terms serendipitously cancel out so that (6.3) still applies.

6.3. Impact of the evaporation rate on bending

The insets of Figures 5(d) and 6(d) show that the deflection of the plate can undergo similar evolutions despite the evaporation rate differing by two orders of magnitude. We now explore the dependence of the deflection on the evaporation rate in more detail. The asymptotically reduced models for the fluid-transport problem are numerically solved when $\textrm{Pe} \ll 1$ , $\textrm{Pe} = 1$ , $\textrm{Pe} = \epsilon ^{-1}$ and $\textrm{Pe} = \epsilon ^{-2}$ . The equations for the reduced models are given by (5.23), (5.20), (5.24), and (5.25), respectively. The plate deflection is then computed by assuming that the film is thin and solving (6.3). In doing so, we take $\mathcal{E} = 1$ .

By plotting the deflection at the end of the beam as a function of $\textrm{Pe}\, t$ , we find that the curves obtained from different values of $\textrm{Pe}$ nearly overlap; see Figure 8(a). Thus, once the time scale of evaporation is accounted for, the evolution of the beam deflection is relatively insensitive to the evaporation rate for this range of Péclet numbers. The small differences between the curves shown in Figure 8(a) could be attributed to two physical mechanisms. The first is that larger Péclet numbers lead to non-uniform in-plane tractions that become increasingly localised at the contact line. The second is that a larger Péclet number will lead to a more rapid depletion of fluid from the free surface, which, in turn, will reduce the evaporative flux $\mathcal{V}(\phi )$ and hence the rate of deflection. The first mechanism can be ruled out by noticing that the curves in Figure 8(a) do overlap for $\textrm{Pe}\, t \lt 0.5$ . However, the in-plane tractions during this time are very different, as shown in Figure 8(b) when $\textrm{Pe} = 1$ and $\textrm{Pe} = \epsilon ^{-2}$ . Thus, the differences in the curves can be attributed to the depletion of fluid near the film surface as $\textrm{Pe}$ increases.

Figure 8. (a) Asymptotic solutions for plate deflection plotted in terms of $\textrm{Pe}\, t$ for five different drying regimes. (b) Asymptotic solutions for the in-plane traction when $\textrm{Pe} = 1$ (top) and $\textrm{Pe} = \epsilon ^{-2}$ (bottom). The solutions are shown when $\textrm{Pe}\, t = 0.1$ , $0.3$ and $0.6$ . The arrows show the direction of increasing time. See text for full details.

To further demonstrate the impact of fluid depletion at the free surface on the rate of deflection, we consider an extreme case and numerically solve (5.25) when $\textrm{Pe}_1 = 100$ , corresponding to $\textrm{Pe} \gg O(\epsilon ^{-2})$ . In this case, a compositional boundary layer of width $O(\textrm{Pe}_1^{-1})$ develops at the free surface, where the Eulerian fluid fraction decreases to approximately $\phi _\infty$ . As a result, the evaporation flux $\mathcal{V}(\phi )$ undergoes a considerable decrease. When the plate deflection is obtained from (6.3) and plotted in terms of $\textrm{Pe}\,t$ , the corresponding curve is markedly different from those obtained from smaller Péclet numbers; see Figure 8(a). When $\textrm{Pe} = O(\epsilon ^{-2})$ or smaller, the deflection initially grows linearly with time. However, when $\textrm{Pe} \gg O(\epsilon ^{-2})$ , the kinetics are qualitatively different and the deflection grows approximately with $t^{1/2}$ .

7. Discussion and conclusion

Stoney’s equation relating the in-plane film stress, $\textsf{S}_{\parallel }$ , to the radius of curvature of the cantilever beam, $R$ , can be written in dimensional form as

(7.1) \begin{align} \textsf{S}_{\parallel } = \frac{E_p H_p^2}{6 H_f^0 R}, \end{align}

where $H_f^0$ denotes the constant film thickness. To fairly compare against the one-dimensional models derived in Section 6, we must consider the thin-film limit $\epsilon \ll \delta$ . By re-dimensionalising (6.3) and using $\partial ^2 w/ \partial X^2 \sim 1/ R$ , we obtain

(7.2) \begin{align} \langle \textsf{S}_{\parallel } \rangle = \frac{\bar{E}_p H_p^2}{6 H_f R}, \end{align}

where $\bar{E}_p = E_p/ (1 - \nu _p^2)$ . Thus, the equations presented here are consistent with Stoney’s. The differing factor of $1 - \nu _p^2$ can be attributed to Stoney working under the assumption of plane stress, whereas (6.3) has been obtained under the assumption of plane strain. The agreement between Stoney’s result and our own suggests that Stoney’s choice of neutral axis is not incorrect after all, as was suggested by Chiu [6]. That being said, Stoney’s derivation is based on the assumption that bending is driven by the in-plane film traction, whereas the asymptotic analysis reveals that the vertical traction can be just as relevant when the film is thin and sufficiently stiff. However, when reducing the plate model to a beam model, the terms associated with the vertical traction cancel out. Although Stoney did not recognise the importance of the vertical traction, their formula remains correct due to this serendipitous cancellation. A direct comparison of (7.1) and (7.2) shows that the in-plane stress in Stoney’s formula should be interpreted as the mean stress across the film thickness. However, given that significant in-plane stress gradients can also arise during drying, it would be more accurate to interpret Stoney’s $\textsf{S}_{\parallel }$ as the mean stress across and along the film.

Tomar et al. [39] conducted cantilever experiments using drops of polymer solutions. By analysing their data using an equation that is similar to (7.1), they found that the film stress increases like $t$ for thin drops and like $t^{1/2}$ for thicker drops. These observations are consistent with the results in Section 6.3, in which the small-time evolution of the beam deflection is seen to change from $t$ to $t^{1/2}$ as the Péclet number, $\textrm{Pe} = V_e^0 H_f^0 \mu _f/ (\epsilon ^2 k_0 E_f)$ increases. For a given polymer solution, the aspect ratio $\epsilon$ of the drop will be fixed and determined by the equilibrium contact angle. Thus, increasing the drop thickness is equivalent to increasing the Péclet number.

Using asymptotic methods, we have systematically derived a hierarchy of simplified models of the cantilever experiment. Our models extend Stoney-like formulae by accounting for the time dependence of the drying process and the generation of non-uniform in-plane and transverse stresses in the film. Fitting the time-dependent solutions to experimental data could lead to new insights into the solid mechanics of drying films and enable difficult-to-measure quantities, such as the Young’s modulus of the film, to be inferred.

Due to the complex rheology of drying films, a number of extensions to this work are possible. For instance, the Young’s modulus of a colloidal mixture can increase by several orders of magnitude during drying [4]. As a result, the system will pass through many of the asymptotic regimes that we have identified. Matched asymptotic expansions could be used to obtain a reduced drying model that spans all of regimes and captures the evolution of the material properties of the film. The asymptotic analysis could also be adapted to films with viscoelastic or elasto-viscoplastic rheologies and applied to other contexts in which the film undergoes thermal expansion, growth or swelling. Finally, the asymptotic framework that we have laid out could also be applied to drying-induced delamination of thin films.