1. Introduction
Spin coating is widely used to apply functional and protective coatings in the manufacturing of electronic and optical components, such as microprocessors, light-emitting diode displays and solar panels. The spin-coating process consists of depositing a coating liquid onto a substrate surface, then rotating the substrate at high speed so that centrifugal force spreads the liquid over the surface (Cohen & Lightfoot Reference Cohen and Lightfoot2011). Once the liquid has formed a thin film over the entire substrate surface, the film is allowed to cure by solvent evaporation, photochemical or other means, to leave a uniform and highly reproducible coating. Current spin-coating techniques, however, are unable to reliably produce uniform coatings on curved substrates (Rich et al. Reference Rich, Jiang, Fukuda and Someya2021). This restricts a wide range of spin-coated products to only flat geometries or singly curved surfaces (by bending flat surfaces without stretching).
Emslie, Bonner & Peck (Reference Emslie, Bonner and Peck1958) developed a simple one-dimensional model for the dynamics of an axisymmetric thin fluid film on a rotating flat substrate, based on the lubrication approximation, considering only the effects of centrifugal force. They showed that, regardless of the initial film profile, a spin-coated film on a flat substrate will tend towards a uniform coating. Models for spin-coated films have since been extended to incorporate additional effects such as gravitational and Coriolis forces, surface tension and curved substrate geometry. The case when surface-tension and moving-contact-line effects are significant on flat substrates was later studied by Wilson, Hunt & Duffy (Reference Wilson, Hunt and Duffy2000). Chen et al. (Reference Chen, Liang, Luo, Zhang and Yang2009) and Liu et al. (Reference Liu, Fang, Meng and Wang2017) presented one-dimensional models to predict the thickness of spin-coated films on convex spherical substrates, and demonstrated close agreement with the thickness of experimentally measured films. Kang, Nadim & Chugunova (Reference Kang, Nadim and Chugunova2016) and Duruk, Boujo & Sellier (Reference Duruk, Boujo and Sellier2021) used similar models to investigate the flow of a film over the entire surface of rotating spherical and spheroidal substrates, paying particular attention to the film dynamics during the transition from gravity- to centrifugal force-driven flow as the angular velocity of the substrate is increased. Finally, these one-dimensional models were generalised by Weidner (Reference Weidner2018) to allow for an arbitrary axisymmetric substrate, such as with ridges and dips.
All of the above models have only considered the flow of axisymmetric films, and cannot capture the more complex dynamics which occurs in non-axisymmetric flows, such as fingering instabilities, as observed experimentally by Fraysse & Homsy (Reference Fraysse and Homsy1994), and angular velocity components introduced by the Coriolis force. In the case of an axisymmetric film, Myers & Charpin (Reference Myers and Charpin2001) showed that the effects of the Coriolis force would have no impact on film thickness, but would induce an angular velocity component that could affect the evolution of non-axisymmetric films. This was supported by experiments by Cho et al. (Reference Cho, Chou, Wang and Tsai2005) demonstrating the deflection of fingers formed at the contact line due to the Coriolis force. On a flat substrate, spreading from a droplet initial condition was investigated by Schwartz & Roy (Reference Schwartz and Roy2004) without assuming an axisymmetric flow and including the Coriolis force, where they were able to reproduce the fingering instability observed by Fraysse & Homsy (Reference Fraysse and Homsy1994) as well as finger deflection similar to Cho et al. (Reference Cho, Chou, Wang and Tsai2005), but only while using a viscosity 14 times smaller than used in experiments.
Better understanding the dynamics of thin liquid films on stationary curved substrates has been the topic of several studies starting with the work of Schwartz & Weidner (Reference Schwartz and Weidner1995), who included the effects of substrate curvature in the lubrication approximation for planar flows. The effects of inertia were later considered in Ruschak & Weinstein (Reference Ruschak and Weinstein2003) for a two-dimensional flow. For the three-dimensional flow of a thin liquid film on a stationary curved substrates, a general theory for thin films driven by surface tension and gravity was developed by Roy, Roberts & Simpson (Reference Roy, Roberts and Simpson2002) and Thiffeault & Kamhawi (Reference Thiffeault and Kamhawi2006), allowing for flow to be modelled on any smooth substrate geometry, and using any parameterisation of the substrate geometry. The possible inclusion of solidification in the governing equations was concurrently considered in Myers, Charpin & Chapman (Reference Myers, Charpin and Chapman2002) and the inertia effects later included in Wray, Papageorgiou & Matar (Reference Wray, Papageorgiou and Matar2017). This general framework was later used by Takagi & Huppert (Reference Takagi and Huppert2010), Balestra, Brun & Gallaire (Reference Balestra, Brun and Gallaire2016), Balestra, Nguyen & Gallaire (Reference Balestra, Nguyen and Gallaire2018), Qin, Xia & Gao (Reference Qin, Xia and Gao2021), Ledda et al. (Reference Ledda, Pezzulla, Jambon-Puillet, Brun and Gallaire2022) and McKinlay, Wray & Wilson (Reference McKinlay, Wray and Wilson2023) to study gravitational drainage and contact-line instabilities over a range of curved substrate geometries. The general theory was extended by Mayo et al. (Reference Mayo, McCue, Moroney, Forster, Kempthorne, Belward and Turner2015) to simulate the dynamics of droplets on leaves (with the addition of disjoining pressure at the contact line).
Since the pioneering work of Howell (Reference Howell2003), few studies have considered the combined effects of a non-trivial substrate kinematics with a complex substrate shape. A notable exception which builds on the large body of literature related to rimming flows on circular cylinder – see for example Evans, Schwartz & Roy (Reference Evans, Schwartz and Roy2004), Rietz et al. (Reference Rietz, Scheid, Gallaire, Kofman, Kneer and Rohlfs2017), Lopes, Thiele & Hazel (Reference Lopes, Thiele and Hazel2018) and Mitchell, Duffy & Wilson (Reference Mitchell, Duffy and Wilson2022) – is the work in Li, Carvalho & Kumar (Reference Li, Carvalho and Kumar2017) investigating the free-surface dynamics of a thin film on a rotating elliptical cylinder. Recently, Duruk et al. (Reference Duruk, Shepherd, Boujo and Sellier2023) modelled coating flows on rotating ellipsoids (with the addition of centrifugal force, but not Coriolis force).
We aim, in this work, to extend these models by presenting an extended general theory which includes all non-inertial forces and therefore can reliably simulate spreading from a droplet initial condition over a rotating, non-axisymmetric curved substrate and shed light on the combined effects of rotation and substrate curvature on the film spreading dynamics. We also aim to clarify to which extent the Coriolis force, which has commonly been assumed (and convincingly been demonstrated) to be negligible for axisymmetric thin film flow configurations, can still be ignored for non-axisymmetric surfaces. This has remained, to the best of the authors’ knowledge, an open question.
In § 2, we will derive a dimensionless general lubrication model for the evolution of a thin fluid film over the surface of an arbitrarily parameterised rotating curved substrate following a similar methodology to Roy et al. (Reference Roy, Roberts and Simpson2002) and Thiffeault & Kamhawi (Reference Thiffeault and Kamhawi2006), incorporating the effects of surface tension, disjoining pressure, gravity, centrifugal and Coriolis forces.
Section 3 gives the parameters and numerical details used to implement this model. In § 4, we present the results from a series of example simulations. In § 4.1, we first consider the spreading of a spin-coated droplet on a flat substrate in order to demonstrate the effects of the Coriolis force at different angular velocities. In § 4.2, we then show the effect of two different non-axisymmetric curved substrates (a parabolic cylinder and saddle) on the spreading of a spin-coated droplet and the onset of the fingering instability in gravity-driven, transitional and centrifugal force-driven flow regimes. Finally, in § 4.3, we present quantitative results comparing the rate of film spreading over the different substrate geometries.
2. Model
2.1. Substrate-based curvilinear coordinate system
Consider a thin film of an incompressible Newtonian fluid on a smooth substrate surface. Let the substrate surface, $\boldsymbol {s}(x^1,x^2)$, be a two-dimensional Riemannian manifold parameterised by $(x^1,x^2)\in \mathbb {R}^2$. Let
be basis vectors tangent to the substrate (for $\alpha \in \{1,2\}$, as with other Greek indices throughout), and a unit vector in the positive normal direction to the substrate, respectively. Note that $\boldsymbol {e}_\alpha$ are not necessarily orthogonal or normalised, and a suitable choice of substrate parameterisation is required to ensure the desired surface orientation. We also define cobasis vectors such that $\boldsymbol {e}_\alpha \boldsymbol {\cdot }\boldsymbol {e}^\beta = \delta _\alpha {}^\beta$ (where $\delta _\alpha {}^\beta$ is the Kronecker delta) and $\boldsymbol {e}^\alpha \boldsymbol {\cdot } \hat {\boldsymbol {n}} = 0$
Let $h(x^1,x^2,t)$ be the thickness of the fluid film above $\boldsymbol {s}(x^1,x^2)$ at time $t$, measured in the positive normal direction to the substrate surface. Each point in the fluid film can then be written as
where $0\leq n\leq h$ is the distance from the substrate surface. This coordinate system is shown in figure 1.
Let $\boldsymbol{\mathsf{G}}$ be the metric tensor on the substrate surface with components
Distance and area elements on the substrate surface can then be written as
adopting the Einstein summation convention for repeated Greek indices, where $G = \det (\boldsymbol{\mathsf{G}}) = G_{11} G_{22} - (G_{12})^2$ is the determinant of the metric. Now, let
define the symmetric substrate curvature tensor, $\boldsymbol{\mathsf{K}} = K_\alpha {}^\beta \boldsymbol {e}^\alpha \boldsymbol {e}_\beta$. We can obtain the components $K_{\alpha \beta } = G_{\beta \gamma } K_\alpha {}^\gamma$ and $K^{\alpha \beta } = G^{\gamma \beta } K_\gamma {}^\alpha$ by the usual way of raising and lowering indices with the metric. We define the mean curvature and Gaussian curvature of the substrate, respectively, as follows:
Accounting for curvature, the basis and cobasis tangent vectors can be extended to the space around the substrate so that $\boldsymbol {e}^+_\alpha \boldsymbol {\cdot }\boldsymbol {e}^{+\beta } = \delta _\alpha {}^\beta + O(n^2)$ (note that here $n^2$ refers to $n$ squared (and not the second contravariant component), since $n$ is scalar. This will be the case for powers of scalar quantities $h$, $L$ and $\varepsilon$, throughout, and should be apparent from context rather than use the notation $\{n\}^2$ or $(n)^2$):
This leads to an extended metric tensor parallel to the substrate surface, $\boldsymbol{\mathsf{G}}^{\,+}$, with components
An area element parallel to (but away from) the substrate surface is then
where the determinant of the extended metric is $G^+ = \det (\boldsymbol{\mathsf{G}}^{\,+}) = \eta ^2 G$, and the expansion or contraction of the coordinate system away from the substrate surface is characterised by
Since the extended metric becomes non-invertible at $\eta =0$, the extended coordinate system is only valid when $\eta > 0$. This condition is always satisfied on a convex or flat substrate surface, however, when the substrate is concave in any direction (that is, when either of the eigenvalues, $k_1$, $k_2$, of $\boldsymbol{\mathsf{K}}$ become strictly positive), this is only guaranteed when the normal coordinate is less than the minimum radius of curvature, $n < 1/\max (k_1,k_2)$.
2.2. Governing equations and boundary conditions
Flow within the fluid film is governed by the Navier–Stokes equations (in coordinate-free form)
where $\rho$ and $\mu$ are the density and viscosity of the fluid, $\boldsymbol {u}(x^1,x^2,n,t) = u^\alpha \boldsymbol {e}^+_\alpha + u^n\hat {\boldsymbol {n}}$ is the velocity in the fluid film (with contravariant components $u^\alpha$ tangent to the substrate, and component $u^n$ normal to the substrate), $p(x^1,x^2,n,t)$ is the pressure in the film and $\boldsymbol {f}(x^1,x^2,n,t) = f^\alpha \boldsymbol {e}^+_\alpha + f^n\hat {\boldsymbol {n}}$ is the acceleration due to the total body force acting on the fluid (referred to herein as simply the body force, which may vary with both time and space depending on the substrate kinematics), and where $\boldsymbol {\nabla }$ and $\Delta = \boldsymbol {\nabla }\boldsymbol {\cdot }\boldsymbol {\nabla }$ are the usual three-dimensional gradient and Laplacian operators in $\mathbb {R}^3$.
Let $\boldsymbol {q}(x^1,x^2,t) = q^\alpha \boldsymbol {e}_\alpha$ be the vector field of volumetric flux over the substrate surface with components
Integrating (2.12) in $n$ with the boundary conditions $u^n\rvert _{n=0}=0$ and $u^n\rvert _{n=h} = \mathrm {d}h/\mathrm {d}t$ gives the following continuity equation in terms of the volume flux:
where $\eta ^* = \eta \rvert _{n=h} = 1 - \kappa h + K h^2$ and $\boldsymbol {\nabla }_{\!S} = \boldsymbol {e}^\alpha \partial /\partial x^\alpha$ is the gradient operator over the substrate surface. The divergence over the substrate surface is given by (see Lebedev & Cloud (Reference Lebedev and Cloud2003), p. 78)
For thin film flows in spin coating, we introduce the centrifugal and Coriolis forces induced by a substrate reference frame rotating at a constant speed (Morin Reference Morin2008, p. 461). The total acceleration due to body forces acting on the fluid at the point $\boldsymbol {r}$ is then
where $g$ is the acceleration due to gravity, $\hat {\boldsymbol {g}}$ is a unit vector in the direction of gravity, $\omega$ is the angular velocity of the substrate and $\hat {\boldsymbol {\omega }}$ is a unit vector in the direction of the axis of rotation (by the right-hand rule).
At the substrate surface, $n=0$, we impose the zero-slip boundary condition $\boldsymbol {u} = \boldsymbol {0}$. At the free fluid surface, the velocity and pressure must satisfy the stress balance
where $\boldsymbol {t}_\alpha ^*$ and $\hat {\boldsymbol {n}}^*$ are tangent and unit normal vectors to the free surface, $\boldsymbol{\mathsf{T}} = \boldsymbol {\nabla }\boldsymbol {u} + \boldsymbol {\nabla }\boldsymbol {u}^\mathrm {T}$ is the strain rate tensor, $p_{a}$ is the ambient pressure, $\gamma$ is the surface tension at the fluid–air interface, $\kappa ^*$ is the mean curvature of the free surface and $\varPi$ is the disjoining pressure at the interface. The disjoining pressure can be modelled as
where $h_{p}$ is a precursor film thickness, $\theta _{e}$ is the equilibrium contact angle between the fluid and the substrate and $m_1$, $m_2$ are constants such that $m_1>m_2>1$. In this paper we will use $m_1 = 3$ and $m_2 =2$ (as used by Schwartz et al. Reference Schwartz, Roy, Eley and Petrash2001; Mayo et al. Reference Mayo, McCue, Moroney, Forster, Kempthorne, Belward and Turner2015).
2.3. Non-dimensionalisation
Let $h_{c}$ be the characteristic thickness of the fluid film, let $L$ be the characteristic length scale of the substrate and let $\varepsilon =h_{c}/L \ll 1$ be the aspect ratio of the film. We now define a dimensionless substrate coordinate system, rescaled by the characteristic lengths $h_{c}$ normal to the substrate and $L$ tangent to the substrate. This gives rise to rescaled position vectors, normal coordinate and (co)basis vectors, temporarily indicated by a tilde
and corresponding dimensionless metric and curvature tensors, $\tilde {\boldsymbol{\mathsf{G}}}$ and $\tilde {\boldsymbol{\mathsf{K}}}$, where
The dimensionless extended (co)basis vectors and metric away from the substrate surface can be written similarly to (2.8a,b) and (2.9)
The determinants of the dimensionless metrics are
where $\eta$ can be expressed in terms of dimensionless quantities as
with the dimensionless mean and Gaussian substrate curvatures
Let $f_{c}$ be the characteristic acceleration due to body forces acting on the fluid, let $\rho L^2 h_{c}$ be a characteristic mass and let $\rho L f_{c}$ be a characteristic pressure (the choice of a suitable characteristic force scale depends on the substrate kinematics, as discussed in § 2.5). Let $u_{c}$ and $t_{c}$ be characteristic velocity and time scales, as follows:
This leads to the dimensionless variables
To be consistent with the definitions, $\tilde {\boldsymbol {u}} = \tilde {u}^\alpha \tilde {\boldsymbol {e}}^+_\alpha + \tilde {u}^n\hat {\boldsymbol {n}}$ and $\tilde {\boldsymbol {q}} = \tilde {q}^\alpha \tilde {\boldsymbol {e}}_\alpha$, the components of the dimensionless velocity and flux are
2.4. Dimensionless governing equations and boundary conditions
Substituting the dimensionless variables (2.27), the Navier–Stokes momentum equation (in coordinate-free form) can be re-written as
where $Re = \rho u_{c} h_{c} / \mu$ is the Reynolds number. Furthermore, expanding and assuming that $\tilde {u}^n \sim \varepsilon \tilde {u}^\alpha$, we can express the Laplacian as (see Lebedev & Cloud (Reference Lebedev and Cloud2003), p. 79)
and the pressure gradient as (see Lebedev & Cloud (Reference Lebedev and Cloud2003), p. 63)
Substituting these into (2.29), the normal and tangential components of the NS momentum equation can be simplified to
Substituting the dimensionless variables (2.27) into (2.15) gives the dimensionless continuity equation
where $\eta ^* = 1 - \varepsilon \tilde {\kappa }\tilde {h} + \varepsilon ^2\tilde {K}\tilde {h}^2$ and $\tilde {\boldsymbol {\nabla }}_{\!S} = L \boldsymbol {\nabla }_{\!S}$ gives the dimensionless divergence over the substrate surface
Rescaled against the characteristic body force, $f_{c}$, and substituting the dimensionless position vector, $\tilde {\boldsymbol {r}} = \tilde {\boldsymbol {s}} + \varepsilon \tilde {n} \hat {\boldsymbol {n}}$, the dimensionless body force is
where $N_{grav} = g/f_{c}$ and $N_{cent} = \omega ^2 L/f_{c}$ are dimensionless groups describing the ratio of gravity and centrifugal force to the characteristic force scale, and $\mathit {Ta} = \omega \rho L h_{c} / \mu$ is the Taylor number, describing the ratio of angular momentum to viscous forces and characterising the strength of the Coriolis force.
In dimensionless form, the zero-slip condition is $\tilde {\boldsymbol {u}} = \boldsymbol {0}$ on $\tilde {n}=0$, and the stress balance simplifies to
where $N_{surf} = \gamma / \rho L^2 f_{c}$ is a dimensionless group describing the ratio of surface tension to the characteristic force scale, and $\tilde {\varPi } = \varPi /\rho L f_{c}$ is the dimensionless disjoining pressure
where $\tilde {h}_{p} = h_{p}/h_{c}$ is the dimensionless precursor film thickness. Here, $N_{surf}$ can take the form of different common dimensionless groups depending on the choice of force scale, as discussed in § 2.5. The dimensionless free-surface curvature can be approximated as
where $\tilde {\kappa }_2 = \tilde {K}_\alpha {}^\beta \tilde {K}_\beta {}^\alpha$, and $\tilde {\Delta }_{S} = \tilde {\boldsymbol {\nabla }}_{\!S} \boldsymbol {\cdot } \tilde {\boldsymbol {\nabla }}_{\!S}$ is the dimensionless Laplacian over the substrate surface. In the absence of disjoining pressure, these boundary conditions are equivalent to those used by Roy et al. (Reference Roy, Roberts and Simpson2002) and Thiffeault & Kamhawi (Reference Thiffeault and Kamhawi2006).
2.5. Choice of characteristic force scale
The form of the dimensionless groups, $N_{surf}$, $N_{grav}$ and $N_{cent}$, characterising the strength of surface tension, gravitational and centrifugal forces, is determined by the choice of the characteristic acceleration due to body forces acting on the fluid film, $f_{c}$. Depending on the fluid properties, film thickness, substrate geometry and kinematics, the dominant force acting on the film may be any of surface tension, gravity or centrifugal force. Each of these could therefore be justifiably chosen as a characteristic force scale when modelling different applications. In each case, one of $N_{surf}$, $N_{grav}$ and $N_{cent}$ reduces to 1, and the others to well-known dimensionless groups – the Bond number, rotational Weber number and rotational Froude number or their reciprocals – as summarised in table 1.
In existing literature, surface tension has typically been chosen as a force scale when modelling very thin films and droplets, including Mayo et al. (Reference Mayo, McCue, Moroney, Forster, Kempthorne, Belward and Turner2015), Roy et al. (Reference Roy, Roberts and Simpson2002) and Thiffeault & Kamhawi (Reference Thiffeault and Kamhawi2006). Gravity has been chosen as a force scale when considering the drainage of films over curved substrates, such as in Balestra et al. (Reference Balestra, Brun and Gallaire2016, Reference Balestra, Nguyen and Gallaire2018), Ledda et al. (Reference Ledda, Pezzulla, Jambon-Puillet, Brun and Gallaire2022) and Takagi & Huppert (Reference Takagi and Huppert2010). Finally, centrifugal force is a natural choice of force scale when modelling spin coating at high speeds, as used in Emslie et al. (Reference Emslie, Bonner and Peck1958) and Liu et al. (Reference Liu, Fang, Meng and Wang2017). When considering situations where several of these forces have a comparable effect on the film dynamics, there is not always a clear choice of force scale. In this case, we propose a more general choice of characteristic force
which ensures that $N_{surf} + N_{grav} + N_{cent} = 1$. This allows the dimensionless groups to be interpreted as the relative strength of each force. In the limiting cases of flow driven entirely by surface tension, gravity or centrifugal force, (2.39) reduces to one of the characteristic force scales listed in table 1. Using this generalised force scale, we can smoothly transition between appropriate scalings for gravity- and centrifugal force-driven flows in order to reconcile the differing time scales of the regimes demonstrated in § 4.
2.6. General lubrication model
Integrating the NS equations using a perturbation expansion approach similar to Roy et al. (Reference Roy, Roberts and Simpson2002) and Thiffeault & Kamhawi (Reference Thiffeault and Kamhawi2006) (see Appendix A), we obtain an expression for the components of the dimensionless volume flux over the substrate surface (omitting the tildes denoting dimensionless variables)
where $\nabla ^\alpha = G^{\alpha \beta } \partial /\partial x^\beta$ are the contravariant components of the substrate gradient, $f^\alpha _{(0)}$ and $f^n_{(0)}$ are the leading-order components of the total body force tangent and normal to the substrate surface
and $\epsilon _\beta {}^\alpha$ is a mixed component of the modified Levi-Civita tensor, defined by
Together with (2.33), this gives a partial differential equation (PDE) describing the evolution of the thickness, $h$, of a thin fluid film on an arbitrary rotating curved substrate. In the absence of substrate rotation, (2.40) is equivalent to (3.13) from Mayo et al. (Reference Mayo, McCue, Moroney, Forster, Kempthorne, Belward and Turner2015). Furthermore, in the absence of both substrate rotation and disjoining pressure, (2.40) is equivalent to (51) from Roy et al. (Reference Roy, Roberts and Simpson2002) and equations (IV.12)–(IV.14) from Thiffeault & Kamhawi (Reference Thiffeault and Kamhawi2006).
3. Numerical implementation
3.1. Simulation parameters
For the example simulations in § 4, we choose physical parameters and material properties based on spin-coating experiments using silicon oil by Wang & Chou (Reference Wang and Chou2001) and Cho et al. (Reference Cho, Chou, Wang and Tsai2005), as shown in table 2. The characteristic film thickness was chosen as $h_{c} = V/L^2$, so that the initial droplet has a dimensionless volume of $1$. For a range of substrate angular velocities from $\omega =0$ to $200\,\textrm {rad}\,\textrm {s}^{-1}$, the corresponding dimensionless groups and characteristic time scales are shown in table 3. With this choice of parameters, we see from the relative magnitudes of $N_{surf}$, $N_{grav}$ and $N_{cent}$ that the dynamics will transition from gravity dominated to centrifugal force dominated over the range of $\omega$ considered here, with surface tension always having a small effect.
3.2. Initial conditions
For the example simulations, we choose a spherical cap initial condition with dimensionless radius of curvature $r_0$, maximum height $h_0$ and surrounding precursor film height $h_{p}$, as shown in figure 2. The (dimensional) volume of the spherical cap is given by (Polyanin & Manzhirov Reference Polyanin and Manzhirov2007, p. 69)
Expressing the maximum height in terms of the contact angle, $h_{c} h_0 = L r_0 (1-\cos \theta _{e})$, and setting $V = h_{c}L^2$ to ensure a dimensionless volume of $1$ (excluding the precursor film), (3.1) can be re-written as
Finally, rearranging and recalling that $\varepsilon = h_{c}/L$ is the film aspect ratio, the radius of curvature and maximum height of a spherical cap initial condition must satisfy
For a contact angle of $\theta _{e} = 10^\circ$, this results in $r_0 = 1.77$ and $h_0 = 6.72$.
We choose a dimensionless precursor film thickness of $h_{p} = 0.1$ (equivalent to $1.5\,\%$ of the initial droplet height) as a practical consideration for efficient simulations (see Appendix B for a discussion of the effect of the precursor film thickness on the rate of film spreading, and also note that Spaid & Homsy (Reference Spaid and Homsy1996) showed that the choice of precursor film thickness does not affect the most amplified wavelength in fingering at the contact line). For all of the examples in § 4, we will always parameterise the substrate surface with $(x,y) = (x^1,x^2)$. In this case, the initial film thickness for a spherical cap with a precursor film is given by
We also consider the flow of a droplet spreading from a randomly perturbed initial condition. To introduce perturbations with a range of wavelengths, we replace $r_0$ in (3.4) with
where $\theta$ is the polar angle in the $x^1x^2$-plane and $a_n, b_n$ are normally distributed random coefficients with a mean of $0$ and standard deviation of $5\times 10^{-3}$.
3.3. Implementation using COMSOL multiphysics and MATLAB
Throughout § 4, (2.33) and (2.40) are solved using the finite-element Coefficient Form PDE Interface in COMSOL Multiphysics 5.6 with LiveLink for MATLAB R2021b. In order to implement (2.33) and (2.40), we introduce the variable $\varGamma = N_{surf}\kappa ^* + \varPi$ and express the problem as a system of second-order PDEs in $h$ and $\varGamma$. The system of PDEs is then solved by pre-computing the coefficients (including the components of the total body force, curvature tensor and modified Levi-Civita tensor) over the solution mesh in MATLAB, then using COMSOL Multiphysics with all default settings and a solver tolerance of $10^{-5}$ (to ensure that the tolerance is less than $\varepsilon ^2$ for the chosen value of $\varepsilon = 4\times 10^{-3}$) to solve for the evolution of $h$ and $\varGamma$ over time. For the example simulations, we use a $200 \times 200$ cell square mesh with linear Lagrange elements over a domain of $(x^1,x^2) \in [-1,1] \times [-1,1]$. Additionally, we choose zero-flux conditions on the boundaries of the computational domain. This does not affect the spreading of the initial droplet as we do not allow any simulations to run sufficiently long for the droplet to reach the edge of the domain.
4. Results
4.1. Effects of the coriolis force in spin coating on a flat substrate
Before investigating the flow of spin-coated films on curved substrates, we will first consider the base case of flow on a flat substrate with particular attention to the effect of the Coriolis force for the physical parameters considered here, which has not been included in previous models for flows on curved substrates (e.g. Weidner Reference Weidner2018; Duruk et al. Reference Duruk, Shepherd, Boujo and Sellier2023). Figure 3 shows the instantaneous initial volume flux in a spherical droplet (3.4) on a rotating flat substrate (given in dimensionless Cartesian coordinates by $\boldsymbol {s}(x^1,x^2) = [x^1,x^2,0]$) in two different flow regimes: transitional flow at low angular velocity ($\omega = 25\,\textrm {rad}\,\textrm {s}^{-1}$, $\mathit {Ta}=0.25$), where gravitational and centrifugal forces have a comparable effect; and centrifugal force-driven flow at high angular velocity ($\omega = 200\,\textrm {rad}\,\textrm {s}^{-1}$, $\mathit {Ta}=1.96$). In the transitional regime (figure 3a), the volume flux is entirely in the radial direction. In the centrifugal force-driven regime (figure 3b), however, the volume flux is deflected in the clockwise direction due to the onset of the Coriolis force, against the anti-clockwise direction of substrate rotation. Compared with gravity and the leading-order component of centrifugal force, which do not scale with film thickness, the Coriolis force scales with $h^2$ (as seen in (2.40) and detailed in the derivation in Appendix A). The effect will therefore be the strongest during the earliest stages of the spin-coating process (as seen in the instantaneous initial volume flux), and diminish as the droplet spreads and thins. Figure 4 shows the contact line at $t=1.5$ with and without the effect of the Coriolis force (by setting $\mathit {Ta}=0$) for the spreading of a perturbed spherical droplet (3.5) on a flat substrate in transitional and centrifugal force-driven regimes. Here, we plot the contact line as the contour where $h = 5h_{p} = 0.5$, noting that the dimensionless film thickness in the wetted area remains above 0.8 even at the end of our simulations. Again, we see that at low angular velocity there is no appreciable difference in the contact line due to the Coriolis force. With increasing angular velocity, the introduction of the Coriolis force leads to a slight deflection of the fingering at the contact line against the direction of substrate rotation with no change in the wavelength or amplitude of the fingering instability.
In the case of centrifugal force-dominated flow, where $f_{c} \approx \omega ^2 L$ (and recalling the choice of characteristic velocity from (2.25a,b)), the Reynolds number may be expressed as
This sets an upper bound of $\mathit {Ta} \lesssim 1$ in order to ensure that $\mathit {Re} \lesssim \varepsilon$, which limits the higher-order terms in (2.40) to at most $O(\varepsilon ^2)$. There is therefore a narrow range of angular velocities around $\mathit {Ta}\sim 1$ where the Coriolis force has an observable effect on the flow (such as in figures 3b and 4b), but where inertial effects can still be ignored. For the model and parameters considered here, we can conclude that while it is important to include the Coriolis force in order to accurately model the film dynamics, its effect is small enough that it is not likely to have a significant impact in practical applications.
4.2. Spin coating on non-axisymmetric curved substrates
We now consider the flow of a spin-coated film on two different non-axisymmetric curved substrates. In §§ 4.2.1 and 4.2.2, we will discuss example simulations of a perturbed spherical droplet spreading on rotating parabolic cylinder- and saddle-shaped substrates, respectively. In each case, the results from additional realisations of the randomly perturbed initial condition are reported in Appendix C.
4.2.1. Parabolic cylinder substrate
Figure 6 shows the evolution of the contact line of a perturbed spherical droplet (3.5) on a rotating parabolic cylinder substrate for angular velocities from $\omega = 0$ to $100\,\textrm {rad}\,\textrm {s}^{-1}$. The substrate is described in dimensionless Cartesian coordinates by
and shown in figure 5(a). In the gravity-driven regime on a stationary substrate ($\omega = 0\,\textrm {rad}\,\textrm {s}^{-1}$, figure 6a), the initial droplet forms rivulets flowing down either side of the substrate in the direction of ${\pm }x^2$. As the angular velocity of the substrate is increased, the film dynamics enters a transitional regime, where there are competing effects from gravitational and centrifugal forces. At $\omega = 25\,\textrm {rad}\,\textrm {s}^{-1}$ (figure 6b), the droplet begins to spread radially in all directions due to centrifugal force, but still with a preference toward the ${\pm }x^2$ directions and forming rivulets down either side of the substrate due to gravity. As the angular velocity is increased further, the dynamics becomes driven primarily by centrifugal force. At $\omega = 50$ and $100\,\textrm {rad}\,\textrm {s}^{-1}$ (figure 6c,d), there is weak gravitational influence and the droplet spreads evenly in all directions. In this regime, the spreading droplet develops a fingering instability at the contact line similarly to on a flat substrate, with the growth rate of the fingers increasing with increasing angular velocity. At the leading order, the volume flux over the substrate surface (2.40) is driven by the component of the total body force tangent to the substrate surface. In the centrifugal force-driven regime, the tangential component of the total body force is greatest along the ridge of the substrate. At high angular velocity ($\omega = 100\,\textrm {rad}\,\textrm {s}^{-1}$, figure 6d), this leads to a deflection of fingers in the ${\pm }x^1$ directions.
4.2.2. Saddle substrate
Figure 7 shows the evolution of the contact line of a perturbed spherical droplet (3.5) on a rotating saddle substrate for angular velocities from $\omega = 0$ to $100\,\textrm {rad}\,\textrm {s}^{-1}$. The substrate is described in dimensionless Cartesian coordinates by
and shown in figure 5(b). Similarly to the parabolic cylinder, in the gravity-driven regime ($\omega = 0\,\textrm {rad}\,\textrm {s}^{-1}$, figure 7a), the initial droplet forms rivulets along the downward-sloping ${\pm }x^2$ directions. At $\omega = 25\,\textrm {rad}\,\textrm {s}^{-1}$ (figure 7b), the droplet again begins to spread radially in all directions, while still forming rivulets in the downward-sloping directions. On a saddle substrate, however, the centrifugal force is not strong enough to overcome the upward slope of the substrate in the ${\pm }x^1$ directions, leading to a more pronounced elliptically shaped contact line. At high angular velocities, the film dynamics on a saddle substrate begins to differ significantly from that on the parabolic cylinder. In particular, for centrifugal force-driven flow on a saddle substrate, the tangential component of the total body force is greatest on the diagonals $x^1 = \pm x^2$, along which the substrate is horizontal. As with the parabolic cylinder substrate, we continue to observe that the thickness of the film behind the moving front is the same as on a flat substrate at high angular velocity ($\omega = 100\,\textrm {rad}\,\textrm {s}^{-1}$). At $\omega = 50$ and $100\,\textrm {rad}\,\textrm {s}^{-1}$ (figure 6c,d), the spreading droplet again develops a fingering instability at the contact line. In this case, however, the fingering does not develop in all directions, and grows in an ‘X’-shape, primarily along the diagonal $z=0$ contours ($x^1 = \pm x^2$) in the direction of the greatest tangential body force.
4.3. Film thickness and surface coverage
To quantify the film evolution, we report in figure 8(a) the film thickness at the substrate centre, $h(0,0,t)$, on the three considered substrates (flat, parabolic cylinder and saddle), for several values of the angular velocity. On all substrates and for all angular velocities, the centre film thickness decreases as $h \sim t^{-1/2}$ from $t \approx 0.1$ onwards and the variations due to the substrate and angular velocity are small compared with the overall variation over time. Furthermore, at a high angular velocity ($\omega = 100\,\textrm {rad}\,\textrm {s}^{-1}$), the centre film thickness remains the same as on a flat substrate while the droplet spreads. The same power law for film thickness evolution on cylindrical and spherical substrates without rotation (Takagi & Huppert Reference Takagi and Huppert2010) and flat substrates with rotation (Melo, Joanny & Fauve Reference Melo, Joanny and Fauve1989) was measured experimentally and predicted analytically by self-similar solutions derived from the lubrication equation with surface tension and hydrostatic pressure neglected away from the contact line.
We also define, and report in figure 8(b), the surface coverage at time $t$ as the substrate surface area within the wetted region
where $\varOmega = [-1,1]\times [-1,1]$ is the domain and
Here, we again choose $h \geq 5 h_{p}$ as the threshold for the wetted area to be consistent with the contact lines plotted in §§ 4.1 and 4.2.
As expected, larger angular velocities lead to faster initial film thinning and, therefore, to thinner films and larger surface coverage at later times, scaling as $C\sim t^{1/2}$. At smaller angular velocities (gravity-driven flows), thinner films are obtained on the parabolic cylinder substrate than on the flat and saddle substrates, which may be ascribed to the larger substrate curvature at the centre. At larger angular velocities (centrifugal force-driven flows), however, the substrate geometry does not significantly affect the overall rate of dispersal over the surface. The key difference in the results between different geometries is the shape of the droplet: although the film spreads a similar amount, its distribution is no longer axisymmetric on the parabolic cylinder and saddle substrates.
5. Conclusions
In this work, we have derived and implemented a model for the evolution of a thin fluid film on a rotating curved substrate surface, allowing for an arbitrary substrate parameterisation, and including the effects of surface tension, disjoining pressure, gravity, centrifugal and Coriolis forces. We can consider this either as an extension of existing models of spin coating (such as Schwartz & Roy Reference Schwartz and Roy2004; Weidner Reference Weidner2018) to allow for a more complex and non-axisymmetric substrate geometry, or as an extension of models for thin film flow over stationary curved substrates (such as Roy et al. Reference Roy, Roberts and Simpson2002; Thiffeault & Kamhawi Reference Thiffeault and Kamhawi2006; Mayo et al. Reference Mayo, McCue, Moroney, Forster, Kempthorne, Belward and Turner2015) to include centrifugal and Coriolis forces due to substrate rotation. We implemented our generalised model using COMSOL Multiphysics and MATLAB to simulate the evolution of spin-coated films on flat and curved substrates at a range of angular velocities.
On a flat substrate, we demonstrated that the Coriolis force has a negligible effect at low substrate angular velocity, and at high angular velocity leads to a slight deflection of the flow against the direction of substrate rotation. This deflection was observed in both the instantaneous initial volume flux at the start of spin coating, and in the fingering instability at the contact line as the spin-coated droplet spreads. Our analysis agrees with the conclusions from Myers & Charpin (Reference Myers and Charpin2001), showing that, while the Coriolis force has only a small effect, there is a range of parameters (when $\mathit {Ta} \sim 1$) where the Coriolis force cannot be neglected, but where inertial effects can still be ignored.
On parabolic cylinder and saddle substrate geometries, we showed example simulations of spin coating at a range of angular velocities to demonstrate the transition from gravity-driven to centrifugal force-driven flow. In both cases, as the angular velocity is increased, a droplet initially at the centre of the substrate shifted from draining in the downward-sloping directions to spreading in the directions where the tangential component of the centrifugal force is greatest. On a parabolic cylinder substrate at $\omega = 100\,\textrm {rad}\,\textrm {s}^{-1}$, this leads to a similar pattern of fingering at the contact line to a flat substrate, however, with the fingers deflected along the direction of the ridge (the ${\pm }x^1$ directions). On a saddle substrate at $\omega = 100\,\textrm {rad}\,\textrm {s}^{-1}$, the effect of substrate geometry was clearer, with the fingers at the contact line growing in an ‘X’-shape (along the $x^1 = \pm x^2$ diagonals). Furthermore, we showed that for both parabolic cylinder and saddle substrates, the film thickness at the substrate centre and the wetted area remain similar over time to those for a droplet spreading on a flat substrate, demonstrating that the key difference between substrate geometries is the shape of the spreading droplet and the direction of fingering.
There are several natural ways in which the present work could be extended. In particular, in order to investigate higher angular velocities where the Coriolis force would have a stronger effect (where $\mathit {Ta} > 1$, and hence $\mathit {Re} > \varepsilon$), however, a model incorporating inertial effects would be required. Additionally, a stability analysis around the contact line on a curved substrate with both gravitational and centrifugal forces may provide a more detailed explanation for the complex fingering patterns observed on a saddle substrate.
Funding
This work is part of the project ‘Development of a multi-axis spin-coating system to coat curved surfaces’ funded by the New Zealand Ministry of Business, Innovation and Employment Endeavour fund (grant number UOCX1904). This funding is gratefully acknowledged.
Declaration of interests
The authors report no conflict of interest.
Appendix A. Integration of the Navier–Stokes equations normal to the substrate
In order to integrate the simplified NS equations (2.32) in the direction normal to the substrate, we introduce a perturbation expansion in $\varepsilon$ for the velocity, pressure and body force as follows, and omit the tildes over dimensionless variables for brevity:
Substituting (2.22) for $\boldsymbol {e}^{+\alpha }$, the $O(1)$ and $O(\varepsilon )$ parts of the tangential body force can be written as
where $\boldsymbol {u}_{(0)} = u^\alpha _{(0)} \boldsymbol {e}_\alpha$ is the leading-order velocity vector. Similarly, the leading-order normal body force is
Here, we see that the $O(1)$ body force terms are constant in $n$, but the $O(\varepsilon )$ terms depend on $n$ both explicitly and implicitly due to variation in $\boldsymbol {u}_{(0)}$ in the normal direction.
Substituting the perturbation expansion (A1), together with (2.22) for $G^{+\alpha \beta }$, the lubrication equations (2.32) at $O(1)$ and $O(\varepsilon )$ are
where $\nabla ^\alpha = G^{\alpha \beta } \partial /\partial x^\beta$ is a contravariant component of the substrate gradient. The corresponding perturbation expansions of the dimensionless boundary conditions are ${u^\alpha _{(0)} = u^\alpha _{(1)} = 0}$ on $n=0$, and
Integrating (A5) with respect to $n$ with boundary conditions (A7) gives the standard leading-order lubrication model
which is equivalent to (34) from Roy et al. (Reference Roy, Roberts and Simpson2002) and (IV.7) from Thiffeault & Kamhawi (Reference Thiffeault and Kamhawi2006).
Let $\omega ^\alpha = \hat {\boldsymbol {\omega }} \boldsymbol {\cdot } \boldsymbol {e}^\alpha$ and $\omega ^n = \hat {\boldsymbol {\omega }} \boldsymbol {\cdot } \hat {\boldsymbol {n}}$ be the tangential and normal components of $\hat {\boldsymbol {\omega }}$ in the substrate coordinate system, and let $\epsilon _\beta {}^\alpha$ be a mixed component of the modified Levi-Civita tensor, $\boldsymbol{\mathsf{E}} = \epsilon _{\alpha \beta } \boldsymbol {e}^\alpha \boldsymbol {e}^\beta$, defined by
This allows the Coriolis force term in (A3) to be written as
The $O(\varepsilon )$ tangential body force can then be expressed explicitly in terms of $n$
Integrating the pressure equation from (A6) in $n$ with the boundary condition (A7) gives
Now substituting (A11) and (A12) and collecting like terms, we can re-write (A6) as
where
Integrating (A13) with boundary conditions (A7), the $O(\varepsilon )$ velocity is
Similar to (A1), we introduce a perturbation expansion for the dimensionless volume flux over the substrate surface
Substituting (A1) into (2.28a–c), the $O(1)$ and $O(\varepsilon )$ terms of the flux are
Integrating, the leading-order flux can be expressed as
and the $O(\varepsilon )$ flux as
Finally, substituting (A14) and combining the $O(1)$ and $O(\varepsilon )$ terms, the components of the dimensionless volume flux over the substrate are
Appendix B. Effect of precursor film thickness
Figure 9 shows the effect of the precursor film thickness on the evolution of the contact line on the flat substrate rotating at $\omega = 100\,\textrm {rad}\,\textrm {s}^{-1}$. Differences are observed when $h_p \leq 0.075$, which are due to a numerical instability developing when the mesh is too coarse compared with the precursor film thickness. Conversely, the results are converged and robust for $h_p \geq 0.1$: while small localised differences can be observed, the overall spreading rate and shape are unaffected. In any case, the conclusions of the present study obtained with $h_p=0.1$ remain unchanged.
Appendix C. Effect of the initial condition
Figure 10 shows the contact lines obtained at $t=1.6$ for 10 different realisations of the randomly perturbed initial condition (3.4)–(3.5) on parabolic cylinder and saddle substrates rotating at $\omega = 100\,\textrm {rad}\,\textrm {s}^{-1}$. While the location of the fingers depends on the individual realisation, this confirms that the observations of § 4.2 are robust, that is, the fingers are consistently deflected in the ${\pm }x^1$ directions on the parabolic cylinder substrate and along the $x^1 = \pm x^2$ diagonals on the saddle substrate.
We note in passing that for this specific value of the angular velocity, the ensemble-averaged contact line (shown with a thick black line) spreads non-axisymmetrically: on the parabolic cylinder substrate, the contact line moves approximately 10 % faster in the ${\pm }x^1$ directions; on the saddle substrate, the contact line develops into a square shape.