Numerical Simulation of the Spin Coating of the Interior of Metal Beverage Cans

Abstract

Using scaling arguments and perturbation theory, we derive the lubrication form of the fluid mechanical equations governing the motion of a thin liquid film on an arbitrarily curved, rotating, axisymmetric substrate. The resulting equations are discretized and then solved numerically using an efficient implicit finite difference algorithm. The primary application for this work is to model the spin coating of the interior of two-piece metal beverage cans, and we consider this problem in some depth. Specifically, we show how adjusting several parameters can eliminate one possible defect in the spin coating process: the tendency for droplets to detach from the substrate when the can is spun at high rotation rates.

1. Introduction

When a liquid is applied to a spinning substrate, or if a pre-wetted substrate is spun, centrifugal forces act to drive any irregularities in the film thickness outward, away from the axis of rotation. The result is that the film becomes thinner and more uniform as the rotation proceeds. Consequently spin coating is used in such applications as coating magnetic storage discs, optical devices, and semiconductor wafers to obtain very thin but uniform films on flat substrates.

When the substrate is curved, centrifugal forces will act to produce a uniform layer only in horizontal regions where the substrate is perpendicular to the axis of rotation. But the coating layer may be very irregular in regions where the substrate is highly curved and the normal vector from the surface is not parallel to the axis of rotation.

In this work we will derive the lubrication form of the fluid mechanical equations for a thin liquid film sprayed on an arbitrarily curved, rotating, axisymmetric substrate. Our goal is to predict how the coating thickness changes with time as a function of the substrate geometry, the rotation rate, the rheological properties of the coating liquid, and the geometry and flux of the spray gun. We then discretize the equations and solve the partial differential equations governing the flow. Using an implicit method of solving the finite difference representation of the partial differential equations, we require a minimum of computer resources.

The theory we develop in this work is used to analyze one specific application: the spray/spin coating of the interior of aluminum beverage containers. When a spinning can is spray painted, centrifugal forces help to cause a more uniform coating layer on the can substrate. Indeed this is the purpose of rotating the can a high spin rates while spray coating the interior. But centrifugal forces can also cause such coating irregularities as drop formation. In this work we will demonstrate how such parameters as the spray gun placement and the rotation rate contribute to how long the beverage can may remain in the spin phase of its coating process before this potential defect occurs.

2. Derivation of the evolution equation

Consider an arbitrarily curved, axisymmetric substrate with the parameter s representing arc length along the substrate. The parameter $\varphi$ represents the circumferential angle and n is the distance normal to the substrate (see Figure 1). u is the velocity parallel to the substrate, in the s direction; v is the velocity normal to the substrate, in the n direction; and w the velocity in the circumferential direction $\varphi$. $\omega$ is the rotation rate around the z axis. Note that while the substrate is axisymmetric, the coating need not be. r is the distance from the can centerline to a given point $\left(s,\varphi \right)$ on the can substrate. For this geometry, the metric coefficients are listed in Table 1. The Appendix of a text book on fluid mechanics by Kee et al. [1] lists the momentum and continuity equations governing the motion of Newtonian fluids on orthogonal, curved substrates. From here we find that for our geometry the continuity equation is given by Figure 1.

Figure 1.

Curvilinear coordinate system. $\varphi$ is the circumferential angle, ${\mathrm{e}}^{\left(s\right)}$ is the unit vector parallel to the substrate, ${\mathrm{e}}^{\left(n\right)}$ the unit vector perpendicular to the substrate, and $r\left(s\right)$ is the distance from the centerline to a given point on the substrate. $\theta$ is the angle of the substrate with respect to the horizontal, $\omega$ is the rotation rate, and $g$ is the acceleration of gravity.

Metric coefficient	Coordinate	Velocity
$a_1=1$	$x_1=s$	$u$
$a_2=1$	$x_2=n$	$v$
$a_3=r$	$x_3=\varphi$	$w$

Table 1.

Coordinates and their metric coefficients.

$$\frac{\partial \left(\mathit{ru}\right)}{\partial s}+\frac{\partial \left(\mathit{rv}\right)}{\partial n}+\frac{\partial \left(w\right)}{\partial \varphi }=0.$$E1

If $\theta$ is the angle the substrate makes with the horizontal, and the gravity vector points in the negative $z$ direction, and then if the can is set upright during the coating process, $\mathit{\rho g}sin\theta$ is the force parallel to the substrate, and $\mathit{\rho g}cos\theta$ is the force perpendicular to the substrate. The centrifugal force points in the $r$ direction due to a rotation about the $z$ axis; consequently $\mathit{\rho \omega r}sin\theta$ is the centrifugal force on the coating perpendicular to the substrate and $\mathit{\rho \omega r}cos\theta$ the force parallel to the substrate. When the coating layer is sufficiently thin and the Reynolds number is sufficiently small that we can neglect the inertial terms, the momentum equations governing the motion of the liquid coating reduce to

$$\begin{array}{l}\frac{\partial p}{\partial s}=-\mathit{\rho g}sin\theta +{\mathit{\rho \omega }}^{2}rcos\theta +\mu \left[2u_{\mathit{ss}}+v_{\mathit{sn}}+u_{\mathit{nn}}\right]+\frac{\mu }{r}\left[2r_s{u}_s+w_{\mathit{s\varphi }}\right]+\frac{\mu }{r^2}\left[w_{\mathit{\varphi \varphi }}-w_{\varphi }-2r_s\left(w_{\varphi }+u\right)\right],\\ \frac{\partial p}{\partial n}=\mathit{\rho g}cos\theta +{\mathit{\rho \omega }}^{2}rsin\theta +\mu \left\{\left[v_{\mathit{ss}}+u_{\mathit{ns}}+2{\mathit{rv}}_{\mathit{nn}}+w_{\mathit{n\varphi }}\right]+\frac{\left[2r_s{v}_s+2r_s{u}_n+2w_n-v_{\mathit{\varphi \varphi }}/r\right]}{r}\right\},\\ \frac{1}{r}\frac{\partial p}{\partial \varphi }=\mu \left\{\left[w_{\mathit{ss}}+v_{\mathit{\varphi n}}+w_s{r}_s\right]+\frac{\left[u_{\mathit{\varphi s}}+w_s+{\mathit{rw}}_{\mathit{nn}}+2w_{\mathit{\varphi \varphi }}+2u_{\varphi }{r}_s^2\right]}{r}\right\}.\end{array}$$

Here $p$ is the pressure, $\rho$ is the density of the liquid, and $\mu$ is the viscosity of the liquid. Meyers et al. [2] found that the Coriolis forces were of the same magnitude as the inertial forces. Because we have neglected the inertial forces on the grounds that the Reynolds number is small, we shall also neglect Coriolis forces in our analysis. We assume that the coating layer obeys the no-slip boundary condition at the surface of the substrate, $n=0$

$$u\left(s,n=0,\varphi \right)=v\left(s,n=0,\varphi \right)=w\left(s,n=0,\varphi \right)=0.$$E2

If $F$, a function of $s$, $n$, $\varphi$, and time $t$, is always zero at the free surface, then it must be given by the coating thickness $h$ minus the value of the coordinate $n$

$$F\left(s,n,\varphi ,t\right)=h\left(s,\varphi ,t\right)-n=0.$$E3

We define $e^{\left(s\right)}$ as the unit vector in the direction parallel to the substrate $s$, $e^{\left(n\right)}$ as the unit vector in the direction perpendicular to the substrate $n$, and $e^{\left(\varphi \right)}$ as the unit vector in the circumferential direction $\varphi$. From Eq. (3) we find that the unit vector perpendicular to the coating surface, denoted by $\widehat{\mathbf{f}}$, is defined by

$$\widehat{\mathbf{f}}=\frac{\nabla F}{\mid \nabla F\mid }=\frac{h_s{\mathbf{e}}^{\left(\mathbf{s}\right)}-{\mathbf{e}}^{\left(\mathbf{n}\right)}+\left(h_{\varphi }/r\right){\mathbf{e}}^{\left(\varphi \right)}}{\sqrt{{\left(h_s\right)}^2+1+{\left(h_{\varphi }/r\right)}^2}}.$$

At the free surface, the substantial derivative of F must be zero. Consequently the kinematic condition on the free surface is given by

$$\frac{\mathit{DF}}{\mathit{Dt}}=0=h_t-{\mathit{uh}}_s+v+w\left(\frac{h_{\varphi }}{r}\right).$$

Employing the substitution

$$\mathcal{L}=\sqrt{{\left(h_s\right)}^2+1+{\left(\frac{h_{\varphi }}{r}\right)}^2},$$

the mean curvature of the free surface for this geometry is given by

$$\kappa =-\nabla \cdot \widehat{\mathbf{f}}+\frac{1}{R^{\left(\mathrm{sub}\right)}}=-\frac{r_s{h}_s+{\mathit{rh}}_{\mathit{ss}}}{r\mathcal{L}}+\frac{h_s^2{h}_{\mathit{ss}}}{{\left(\mathcal{L}\right)}^3}-\frac{h_{\mathit{\varphi \varphi }}}{r^2\mathcal{L}}+\frac{h_{\varphi }^2{h}_{\mathit{\varphi \varphi }}}{r^4{\left(\mathcal{L}\right)}^3}+\frac{1}{R^{\left(\mathrm{sub}\right)}},$$E4

at $n=h$. Here $R^{\left(\mathrm{sub}\right)}$ is the radius of curvature of the substrate. Adding the curvature of the free surface relative to the substrate with the curvature of the substrate, given by $1/R^{\left(\mathrm{sub}\right)}$, to obtain the total curvature of the free surface has been shown to be valid by Schwartz and Weidner [3] as long as $h_s\ll 1$ and $\mathit{h\kappa }\ll 1$.

If the atmospheric pressure is zero, then the tensor equation relating the change in pressure across the free surface due to surface tension is given by

$$\left[{\tau }_{\mathit{ij}}+p{\delta }_{\mathit{ij}}\right]{\widehat{f}}_j=-\mathit{\sigma \kappa }{\widehat{f}}_i.$$E5

Here $\sigma$ is the surface tension of the liquid, and ${\tau }_{\mathit{ij}}$ is the stress tensor in the liquid at the free surface. For our geometry, the components of the stress tensor are given by

$$\begin{array}{l}{\tau }_{\mathit{ss}}=2\mu u_s,{\tau }_{\mathit{nn}}=2\mu v_n,{\tau }_{\mathit{\varphi \varphi }}=\frac{2\mu }{r}\left[w_{\varphi }+{\mathit{ur}}_s\right],\\ {\tau }_{\mathit{sn}}={\tau }_{\mathit{ns}}=\mu \left[v_s+u_n\right],{\tau }_{\mathit{n\varphi }}={\tau }_{\mathit{\varphi n}}=\mu \left[w_n+\frac{v_{\varphi }}{r}\right],\\ {\tau }_{\mathit{s\varphi }}={\tau }_{\mathit{\varphi s}}=\mu \left[\frac{u_{\varphi }}{r}+w_s-\frac{{\mathit{wr}}_s}{r}\right].\end{array}$$

If we sum the three equations found in Eq. (5) over $j=s,n,\varphi$ for each $i=s$, $i=n$, and $i=\varphi$, we have the following three equations:

$$-\mathit{\sigma \kappa }{h}_s={\mathit{ph}}_s+{\tau }_{\mathit{ss}}{h}_s+{\tau }_{\mathit{sn}}+{\tau }_{\mathit{s\varphi }}\frac{h_{\varphi }}{r},$$E6

$$-\mathit{\sigma \kappa }=p+{\tau }_{\mathit{ns}}{h}_s+{\tau }_{\mathit{nn}}+{\tau }_{\mathit{n\varphi }}\frac{h_{\varphi }}{r},$$E7

$$-\mathit{\sigma \kappa }\frac{h_{\varphi }}{r}=p\frac{h_{\varphi }}{r}+{\tau }_{\mathit{\varphi s}}{h}_s+{\tau }_{\mathit{\varphi n}}+{\tau }_{\mathit{\varphi \varphi }}\frac{h_{\varphi }}{r},$$E8

all evaluated at the free surface $n=h$. Because we have assumed that the pressure in the air is zero, the pressure in the liquid is positive if $\kappa$ is negative. If we eliminate the pressure $p$ and the curvature of the free surface $\kappa$ from the three Eqs. (6)–(8), we find the following equations are valid at the free surface $n=h$:

$${\tau }_{\mathit{sn}}\left(1-h_s^2\right)+\left({\tau }_{\mathit{ss}}-{\tau }_{\mathit{nn}}\right)h_s-{\tau }_{\mathit{n\varphi }}\frac{h_s{h}_{\varphi }}{r}+{\tau }_{\mathit{s\varphi }}\frac{h_{\varphi }}{r}=0,$$E9

$${\tau }_{\mathit{n\varphi }}\left(1-{\left[\frac{h_{\varphi }}{r}\right]}^2\right)+{\tau }_{\mathit{\varphi \varphi }}\frac{h_{\varphi }}{r}+{\tau }_{\mathit{s\varphi }}{h}_s-{\tau }_{\mathit{ss}}{h}_s\frac{h_{\varphi }}{r}-{\tau }_{\mathit{ns}}\frac{h_s{h}_{\varphi }}{r}=0.$$E10

We now scale the dependent and independent variables with various characteristic lengths of the substrate geometry. These include $h_0$, a characteristic coating thickness, $R_c$ the maximum radius of the substrate, $u_0$ a characteristic velocity, and $\alpha$ a small parameter representing the ratio of a typical length perpendicular to the substrate to a typical length parallel to the substrate. Denoting dimensionless variables by bars, we have

$$\begin{array}{l}\overline{\varphi }=\frac{\alpha R_c\varphi }{h_0},\overline{s}=\frac{\mathit{\alpha s}}{h_0},\overline{n}=\frac{n}{h_0},\overline{w}=\frac{w}{u_0},\overline{u}=\frac{u}{u_0},\overline{v}=\frac{v}{\alpha u_0},\\ \overline{r}=\frac{r}{R_c},\overline{t}=\frac{\alpha u_{0}t}{h_0},\overline{\kappa }=R_c\kappa ,\overline{p}=\frac{{\mathit{ph}}_0}{u_0\mu },\overline{Q}=\frac{Q}{u_0{h}_0}.\end{array}$$

If $Q$ is the flux, with units of length cubed over time, then the flux in the s and $\varphi$ directions are defined by

$${\overline{Q}}^{\left(\overline{s}\right)}={\int }_0^{\overline{h}}\overline{u}d\overline{n},{\overline{Q}}^{\left(\overline{\varphi }\right)}={\int }_0^{\overline{h}}\overline{w}d\overline{n}.$$

We can eliminate the pressure from the $s$, $n$, and $\varphi$ momentum equations using identities for partial differentiation of functions:

$$\begin{array}{cc}{\left({\overline{p}}_{\overline{s}}\right)}_{\overline{n}}-{\left({\overline{p}}_{\overline{n}}\right)}_{\overline{s}}& =0=\mu {\overline{u}}_{\overline{n}\overline{n}\overline{n}}+\mathit{\alpha \mu }\left(\frac{-{\overline{w}}_{\overline{\varphi }\overline{n}}}{{\overline{r}}^2}-\frac{2{\overline{u}}_{\overline{n}}{\overline{r}}_{\overline{s}}}{{\overline{r}}^2}+\frac{{\overline{w}}_{\overline{n}\overline{s}}}{\overline{r}}\right)\\ & +\alpha \left[-\frac{\mathit{\rho g}sin{\mathit{\theta \theta }}_{\overline{s}}}{R_c}+{\mathit{\rho \omega }}^2{\overline{r}}_{\overline{s}}sin\theta +\frac{{\mathit{\rho \omega }}^{2}cos{\mathit{\theta \theta }}_{\overline{s}}}{R_c}\right]+O\left({\alpha }^2\right),\end{array}$$E11

$${\left({\overline{p}}_{\overline{\varphi }}\right)}_{\overline{n}}-{\left({\overline{p}}_{\overline{n}}\right)}_{\overline{\varphi }}=0=\mu {\overline{w}}_{\overline{n}\overline{n}\overline{n}}\overline{r}+O\left({\alpha }^2\right).$$E12

The scaled form of the continuity equation is given by

$$\frac{\partial }{\partial \overline{s}}\left(\overline{r}\overline{u}\right)+\frac{\partial }{\partial \overline{n}}\left(\overline{r}\overline{v}\right)+\frac{\partial }{\partial \overline{\varphi }}\left(\overline{w}\right)=0.$$E13

At the substrate, $\overline{n}=0$, we have the scaled representation of the no-slip boundary conditions:

$$\overline{u}\left(\overline{s},\overline{n}=0,\overline{\varphi },\overline{t}\right)=\alpha \overline{v}\left(\overline{s},\overline{n}=0,\overline{\varphi },\overline{t}\right)=\overline{w}\left(\overline{s},\overline{n}=0,\overline{\varphi },\overline{t}\right)=0.$$E14

The scaled form of the kinematic condition is shown to be

$${\overline{h}}_{\overline{t}}={\overline{u}\overline{h}}_{\overline{s}}-\overline{v}+\overline{w}\frac{{\overline{h}}_{\overline{\varphi }}}{\overline{r}},$$E15

at $\overline{n}=\overline{h}$. The boundary conditions, Eq. (9) and Eq. (10), are given by

$${\overline{u}}_{\overline{n}}+O\left({\alpha }^2\right)=0,{\overline{w}}_{\overline{n}}+O\left({\alpha }^2\right)=0,$$E16

at $\overline{n}=\overline{h}$.

Instead of transforming the equations representing the pressure discontinuity across the liquid interface, given by Eq. (9) and Eq. (10), we will eliminate the pressure term in Eqs. (9) and (10) using the momentum equations and the identities for the partial differentiation of implicit functions. Eq. (7) gives the pressure in the liquid at the free surface $\overline{p}\left(\overline{s},\overline{n}=\overline{h},\overline{\varphi },\overline{t}\right)$. The momentum equations give us ${\overline{p}}_{\overline{n}}$, ${\overline{p}}_{\overline{s}}$, and ${\overline{p}}_{\overline{\varphi }}$ in the liquid. The directional derivatives of $\overline{p}$, along with $\overline{n}=\overline{h}\left(\overline{s},\overline{\varphi },\overline{t}\right)$, in the $\overline{s}$ and $\overline{\varphi }$ directions, respectively, are

$$\begin{array}{l}{\left.\left[{\overline{p}}_{\overline{s}}+{\overline{p}}_{\overline{n}}{\overline{h}}_{\overline{s}}\right]\right|}_{\overline{n}=\overline{h}}=\frac{\partial }{\partial \overline{s}}\overline{p}\left(\overline{s},\overline{n}=\overline{h},\overline{\varphi },\overline{t}\right),\\ {\left.\left[{\overline{p}}_{\overline{\varphi }}+{\overline{p}}_{\overline{n}}{\overline{h}}_{\overline{\varphi }}\right]\right|}_{\overline{n}=\overline{h}}=\frac{\partial }{\partial \overline{\varphi }}\overline{p}\left(\overline{s},\overline{n}=\overline{h},\overline{\varphi },\overline{t}\right).\end{array}$$

Plugging the pressure at the interface, found by Eq. (7), into the momentum equations yields

$$\begin{array}{ll}{\left.\mu {\overline{u}}_{\overline{n}\overline{n}}\right|}_{n=h}=& \lambda sin\theta -\mathrm{\Gamma }\left(\overline{r}cos\theta \right)+\alpha \left[\frac{\mathrm{\Lambda }}{H}{\overline{\kappa }}_{\overline{s}}+\lambda sin\theta {\overline{h}}_{\overline{s}}+\mathrm{\Gamma }\overline{r}cos\theta {\overline{h}}_{\overline{s}}\right]\\ & +\mathit{\alpha \mu }\left[\frac{2{\overline{r}}_{\overline{s}}\overline{u}}{{\overline{r}}^2}+\frac{{\overline{w}}_{\varphi }}{\overline{r}}-\frac{2{\overline{w}}_{\overline{n}}}{{\overline{r}}^2}{\overline{h}}_{\overline{s}}\right]+O\left({\alpha }^2\right),\end{array}$$E17

$${\left.\mu {\overline{w}}_{\overline{n}\overline{n}}\right|}_{\overline{n}=h}=\alpha \left[\frac{\mathrm{\Lambda }}{H}\frac{{\overline{\kappa }}_{\overline{\varphi }}}{\overline{r}}+\lambda cos\theta {\overline{h}}_{\overline{\varphi }}+\mathrm{\Gamma }\overline{r}sin\theta {\overline{h}}_{\overline{\varphi }}\right]+\mathit{\alpha \mu }\frac{{\overline{w}}_{\overline{s}}}{r}+O\left({\alpha }^2\right),$$E18

where

$$\lambda =\frac{\rho {\mathit{gR}}_c^2}{\mu u_0},\mathrm{\Gamma }=\frac{{\mathit{\rho \omega }}^2{R}_c^3}{\mu u_0},\mathrm{\Lambda }=\frac{\mathit{\sigma \alpha }}{\mu u_0},H=\frac{h_0}{R_c}.$$

Here the $\kappa$ terms have been assumed to be of order $O\left(1\right)$, while Eq. (4) would indicate that they are $O\left({\alpha }^2\right)$. In order to bring the surface tension terms into the $O\left(\alpha \right)$ problem, we are in effect assuming that $\mathrm{\Lambda }$ is $O\left(1/{\alpha }^2\right)$. Others, including Benney [4] and Atherton and Homsy [5], have used a similar “specific promotion” for similar thin film problems.

We now expand the scaled velocities in a regular perturbation series expansion in powers of the small parameter $\alpha$:

$$\begin{array}{l}\overline{u}={\overline{u}}^{\left(0\right)}+\alpha {\overline{u}}^{\left(1\right)}+\dots {\alpha }^m{\overline{u}}^{\left(m\right)},\\ \overline{v}={\overline{v}}^{\left(0\right)}+\alpha {\overline{v}}^{\left(1\right)}+\dots {\alpha }^m{\overline{v}}^{\left(m\right)},\\ \overline{w}={\overline{w}}^{\left(0\right)}+\alpha {\overline{w}}^{\left(1\right)}+\dots {\alpha }^m{\overline{w}}^{\left(m\right)}.\end{array}$$

We then use boundary conditions to determine the appropriate constants of integration and solve for ${\overline{u}}^{\left(0\right)}$, ${\overline{v}}^{\left(0\right)}$, and ${\overline{w}}^{\left(0\right)}$ and ${\overline{u}}^{\left(1\right)}$, ${\overline{v}}^{\left(1\right)}$, and ${\overline{w}}^{\left(1\right)}$. The $O\left(1\right)$ velocities are given by

$${\overline{u}}^{\left(0\right)}\left(\overline{s},\overline{n},\overline{\varphi }\right)=A_0\left[\frac{{\overline{n}}^2}{2}-\overline{h}\overline{n}\right],{\overline{w}}^{\left(0\right)}=0,$$E19

where $A_0=\lambda sin\theta -\mathrm{\Gamma }\overline{r}cos\theta$. We now determine the order $\alpha$ velocities

$${\overline{u}}^{\left(1\right)}=A_0\frac{2{\overline{r}}_{\overline{s}}}{{\overline{r}}^2}\left[\frac{{\overline{n}}^4}{24}-\frac{{\overline{h}\overline{n}}^3}{6}+\frac{{\overline{h}}^2{\overline{n}}^4}{4}-\frac{{\overline{h}}^3\overline{n}}{6}\right]+B_0\left[\frac{{\overline{n}}^2}{2}-\overline{h}\overline{n}\right]+C_0\left[{\overline{h}}^2\overline{n}-\frac{{\overline{h}\overline{n}}^2}{2}\right],$$E20

$${\overline{w}}^{\left(1\right)}=D_0\left[\frac{{\overline{n}}^2}{2}-\overline{h}\overline{n}\right],$$E21

with the constants $B_0$, $C_0$, and $D_0$ given by

$$\begin{array}{l}{B}_0=\left[\frac{\mathrm{\Lambda }}{H}{\overline{\kappa }}_{\overline{s}}+\lambda cos\theta {\overline{h}}_{\overline{s}}+\mathrm{\Gamma }\overline{r}sin\theta {\overline{h}}_{\overline{s}}\right],\\ C_0=\left[-\lambda sin{\mathit{\theta \theta }}_{\overline{s}}+\mathrm{\Gamma }{\overline{r}}_{\overline{s}}sin\theta +\mathrm{\Gamma }\overline{r}cos{\mathit{\theta \theta }}_{\overline{s}}\right],\\ D_0=\left[\frac{\mathrm{\Lambda }}{H}\frac{{\overline{\kappa }}_{\overline{\varphi }}}{\overline{r}}+\lambda cos\theta \frac{{\overline{h}}_{\overline{\varphi }}}{\overline{r}}+\mathrm{\Gamma }\overline{r}sin\theta \frac{{\overline{h}}_{\overline{\varphi }}}{\overline{r}}\right].\end{array}$$

Using the appropriate boundary conditions, the order $\alpha$ velocity in the normal direction $n$ may be determined by integrating the continuity equation. But if we limit our consideration to the order $\alpha$ problem, we can write the evolution equation in divergence form, negating the necessity of determining ${\overline{v}}^{\left(1\right)}$. Integrating Eq. (13) between the boundaries $\overline{n}=0$ and $\overline{n}=\overline{h}$, and employing Leibniz’s rule, we have the kinematic condition

$${\overline{h}}_{\overline{t}}=-\frac{1}{\overline{r}}\left[\frac{\partial \left({\overline{r}\overline{Q}}^{\left(\overline{s}\right)}\right)}{\partial \overline{s}}+\frac{\partial {\overline{Q}}^{\left(\overline{\varphi }\right)}}{\partial \overline{\varphi }}\right].$$

The evolution equation, to order $\alpha$, is then determined to be

$$\begin{array}{ll}{\overline{h}}_{\overline{t}}=-& \frac{1}{\overline{r}}\frac{\partial }{\partial \overline{s}}\left[\left[\left(-\lambda sin\theta +\mathrm{\Gamma }\overline{r}cos\theta \right)-\alpha \left(+\frac{\mathrm{\Lambda }}{H}{\overline{\kappa }}_{\overline{s}}+\lambda cos\theta {\overline{h}}_{\overline{s}}-\mathrm{\Gamma }\overline{r}sin\theta {\overline{h}}_{\overline{s}}\right)\right]\overline{r}\frac{{\overline{h}}^3}{3}\right.\\ & \left.+\alpha \left(-\lambda sin{\mathit{\theta \theta }}_{\overline{s}}+\mathrm{\Gamma }{\overline{r}}_{\overline{s}}sin\theta +\mathrm{\Gamma }\overline{r}cos{\mathit{\theta \theta }}_{\overline{s}}\right)\overline{r}\frac{{\overline{h}}^4}{8}+\alpha \left(-\lambda sin\theta +\mathrm{\Gamma }\overline{r}cos\theta \right)\frac{{\overline{r}}_{\overline{s}}}{\overline{r}}\frac{{\overline{h}}^5}{3}\right]\\ & -\frac{1}{\overline{r}}\frac{\partial }{\partial \overline{\varphi }}\left[\alpha \left(-\frac{\mathrm{\Lambda }}{H}\frac{{\overline{\kappa }}_{\overline{\varphi }}}{\overline{r}}-\lambda cos\theta \frac{{\overline{h}}_{\overline{\varphi }}}{\overline{r}}-\mathrm{\Gamma }\overline{r}sin\theta \frac{{\overline{h}}_{\overline{\varphi }}}{\overline{r}}\right)\frac{{\overline{h}}^3}{3}\right].\end{array}$$

The dimensional evolution equation is found to be given by

$$\begin{array}{ll}{h}_t& =-\frac{1}{\mathit{\mu r}}\frac{\partial }{\partial s}\left[\left(-\mathit{\rho g}sin\theta +{\mathit{\rho \omega }}^{2}rcos\theta -{\mathit{\sigma \kappa }}_s-\mathit{\rho g}cos\theta h_s-{\mathit{\rho \omega }}^{2}rsin\theta h_s\right)r\left(\frac{h^3}{3}\right)\right]\\ & -\frac{1}{\mathit{\mu r}}\frac{\partial }{\partial s}[\left(-\mathit{\rho g}sin{\mathit{\theta \theta }}_s+\sigma r_{s}sin\theta +{\mathit{\rho \omega }}^{2}rcos{\mathit{\theta \theta }}_s\right)r\left(\frac{h^4}{8}\right)\\ & +\frac{r_s}{r}\left(-\mathit{\rho g}sin\theta +{\mathit{\rho \omega }}^{2}rcos\theta \right)\left(\frac{h^5}{3}\right)]\\ & -\frac{1}{\mathit{\mu r}}\frac{\partial }{\partial \varphi }\left[\left(-\sigma \frac{{\kappa }_{\varphi }}{r}-\mathit{\rho g}cos\theta \frac{h_{\varphi }}{r}-{\mathit{\rho \omega }}^{2}rsin\theta \frac{h_{\varphi }}{r}\right)\left(\frac{h^3}{3}\right)\right].\end{array}$$E22

3. Nondimensionalization

We shall employ the maximum radius of the substrate as our length scale, $r_{\max }=R_c$, and scale time with the quantity $t^{\ast }=3\mu R_c/\sigma$. We include the following substitutions:

$$\widehat{h}=\frac{h}{R_c},\widehat{r}=\frac{r}{R_c},\widehat{s}=\frac{s}{R_c},\widehat{\kappa }=\kappa R_c,\widehat{\varphi }=\varphi ,\widehat{t}=\frac{t}{t^{\ast }}.$$

Substituting these scale factors in the dimensional evolution equation, Eq. (22), leads to the nondimensional evolution equation

$$\begin{array}{ll}{\widehat{h}}_{\widehat{t}}& =-\frac{1}{\widehat{r}}\frac{\partial }{\partial \widehat{s}}\left[\left(-{\lambda }_{2}sin\theta +{\mathrm{\Gamma }}_2\widehat{r}cos\theta -{\widehat{\kappa }}_{\widehat{s}}-{\lambda }_{2}cos\theta {\widehat{h}}_{\widehat{s}}-{\mathrm{\Gamma }}_2\widehat{r}sin\theta {\widehat{h}}_{\widehat{s}}\right){\widehat{r}\widehat{h}}^3\right]\\ & -\frac{1}{\widehat{r}}\frac{\partial }{\partial \widehat{s}}\left[\left(+{\lambda }_{2}sin{\mathit{\theta \theta }}_{\widehat{s}}-{\mathrm{\Gamma }}_2{\widehat{r}}_{\widehat{s}}sin\theta +{\mathrm{\Gamma }}_2\widehat{r}cos{\mathit{\theta \theta }}_{\widehat{s}}\right)\widehat{r}\left(\frac{3{\widehat{h}}^4}{8}\right)\right]\\ & -\frac{1}{\widehat{r}}\frac{\partial }{\partial \widehat{s}}\left[\left[-{\lambda }_{2}sin\theta +{\mathrm{\Gamma }}_{2}cos\theta \right]\frac{{\widehat{r}}_{\widehat{s}}}{\widehat{r}}\left({\widehat{h}}^5\right)\right],\\ & -\frac{1}{\widehat{r}}\frac{\partial }{\partial \widehat{\varphi }}\left[\left[-\frac{{\widehat{\kappa }}_{\widehat{\varphi }}}{\widehat{r}}-{\lambda }_{2}cos\theta \frac{{\widehat{h}}_{\widehat{\varphi }}}{\widehat{r}}-{\mathrm{\Gamma }}_2\widehat{r}sin\theta \frac{{\widehat{h}}_{\widehat{\varphi }}}{\widehat{r}}\right]{\widehat{h}}^3\right],\end{array}$$E23

in nondimensional units. Here ${\lambda }_2$ and ${\mathrm{\Gamma }}_2$ are given by

$${\lambda }_2=\frac{\rho {\mathit{gR}}_c^2}{\sigma },{\mathrm{\Gamma }}_2=\frac{{\mathit{\rho \omega }}^2{R}_c^3}{\sigma }.$$

The geometry of the substrate is delineated by a schematic drawing which gives us the radius of each circular arc, $R_{\mathrm{i}}$, the subtended angle ${\phi }_{\mathrm{i}}$ of each element, and the length of each straight segment (see Figure 2). Consequently there is a discontinuity in the curvature between adjacent elements, but the slope is always continuous. In order to use the lubrication approximations, the coating thickness must always be much smaller than a typical substrate length, $R_{\mathrm{i}}{\phi }_{\mathrm{i}}$ or $l_{\mathrm{i}}$. At $t=0$ we assume that there is already a very thin paint layer on the substrate and always use a sufficiently small time step such that $\widehat{h}\left(\widehat{s},\widehat{t}\right)$ is never zero. Consequently we do not have to introduce the effects of a contact angle in our analysis. We not only assume that there is no $\varphi$ variation in the substrate but that the coating is also axisymmetric. Consequently we only have to solve for the coating height as a function of s and t.

Figure 2.

Substrate made up of straight segments and circular arcs of radius $R_{\mathrm{i}}$ and subtended angle ${\phi }_{\mathrm{i}}$.

4. Application of coating layer

We assume that the coating layer is laid down over time by a spray gun that emits a fan of gas that is directed toward the can substrate. Thus the accumulation of coating with time due to this fan is a function of s and t: $\partial h^{\mathrm{fan}}\left(s,t\right)/\partial t$. The spray sector is bounded by the angles $\mathrm{\Delta }\phi$ and $\mathrm{\Delta }\psi$ as shown in Figure 3. We assume that the subtended angle $\mathrm{\Delta }\psi$ is significantly smaller than the angle $\mathrm{\Delta }\phi$ so that if the spray fan geometry is an ellipse, it will have an eccentricity near 1.0.

Figure 3.

Assumed spray pattern.

If $Q_{\mathrm{tot}}^{\mathrm{fan}}$ is the total flux from the fan, and $q$ the flux in an elemental area of the fan, then

$$Q_{\mathrm{tot}}^{\mathrm{fan}}=\underset{{\phi }_1}{\overset{{\phi }_2}{\int }}\underset{{\psi }_1}{\overset{{\psi }_2}{\int }}q\left(\phi ,\psi \right)\mathit{\text{dφdψ}}.$$

If the flux does not vary significantly in the $\psi$ direction, then the flux is constant in the elemental shaded sector in Figure 3. Consequently the flux per unit distance across the thin dimension of the fan is given by $q\left(\phi \right)/\left(\mathrm{ϵ}\upsilon \right)$, where $\mathrm{ϵ}=\mathrm{\Delta }\psi$. When the width of the spray fan in the $\psi$ direction is less than one-half the circumference $2\mathit{\pi r}\left(s\right)$ wetted by an element of the fan (see Figure 4), the average flux deposited per revolution is given by

Figure 4.

Beverage can and spray pattern geometry. Here $\upsilon$ is the distance from the orifice of the spray gun to a given point on the can substrate.

$$\left(q/\mathrm{ϵ}\upsilon \right)\left(\mathrm{ϵ}\upsilon /\left(2\mathit{\pi r}\left(s\right)\right)\right)=q/\left(2\mathit{\pi r}\left(s\right)\right).$$E24

When expressed as a function of $\partial h^{\mathrm{fan}}/\partial t$, Eq. (24) becomes

$$\mathit{qd\phi }=\frac{\partial h^{\mathrm{fan}}}{\partial t}\mathit{ds}2\mathit{\pi r}\left(s\right).$$

From Figure 4 we have the relation

$$sin\gamma =\mathit{\upsilon d\phi }/\mathit{ds},$$

For a sufficiently thin fan, where the point on the substrate is sufficiently far from the centerline $r=0$, the circumferential-average rate of increase of height is given by

$$\frac{\partial h^{\mathrm{fan}}}{\partial t}=q\left(\phi \right)sin\gamma /\left(2\mathit{\pi r}\left(s\right)\upsilon \right).$$E25

When we are near the centerline, this formula must be modified to

$$\frac{\partial h^{\mathrm{fan}}}{\partial t}=\left[1/\upsilon \left(s\right)\right]\left[q\left(s\right)/\left(\mathrm{ϵ}\upsilon \left(s\right)\right)\right]sin\gamma \left(s\right)\cdot Min\left[(\mathrm{ϵ}\upsilon \left(s\right)/2\mathit{\pi r}\left(s\right),1/2\right],$$E26

to account for the fact that here the substrate is constantly being reached by the spray fan.

We also assume that there is a secondary “gas” which is uniform in the interior of the can and results in a constant $\partial h_{\mathrm{gas}}/\partial t$. If $0<\zeta <100$ represents the percentage of secondary gas, then $\partial h_{\mathrm{gas}}/\partial t$ from this source is given by

$$\frac{\partial h_{\mathrm{gas}}^{\mathrm{fan}}}{\partial t}=\frac{\zeta }{100}\frac{Q_{\mathrm{tot}}}{{\mathrm{A}}_{\mathrm{tot}}},$$

where ${\mathrm{A}}_{\mathrm{tot}}$ is the total area of the interior of the can. We also must include the effects of geometric shading. For example, the dome blocks some flux from the spray gun from reaching the inner wall of the moat.

5. Pendant drop of maximum volume

In an effort to determine the maximum value a droplet can form without detaching from the substrate, we will consider an axisymmetric droplet forming on the underside of a ceiling, as illustrated in Figure 5. Here $z$ is the distance perpendicular to the substrate, $x$ is the distance from the centerline of the drop parallel to the substrate, and $h$ is the distance from the substrate to the free surface of the drop. $g$ is the acceleration of gravity, and the contact angle with the substrate $z=0$ is given by ${\theta }_c$. The pressure jump across the free surface due to surface tension is given by $\mathit{\sigma \kappa }$. At equilibrium this must equal the hydrostatic pressure due to gravity, or

Figure 5.

Profile view of axisymmetric pendant drop.

$$\mathit{\sigma \kappa }-\mathit{\rho gz}=\text{constant}.$$

If we nondimensionalize the problem by scaling $z$, $x$, and $h$ with the capillary length, $l_c=\sqrt{\sigma /\left(\mathit{\rho g}\right)}$, then the drop shape is defined by the following differential equation:

$$\frac{x_{\mathit{zz}}}{{\left(1+x_z^2\right)}^{3/2}}-\frac{1}{x{\left(1+x_z^2\right)}^{1/2}}+z={\kappa }_0.$$E27

Here ${\kappa }_0$ is the mean curvature at $x=0$ and $z=h$.

Using a standard Runge–Kutta method, we can numerically integrate Eq. (27) assuming ${\theta }_c=0$. We find that the maximum volume is given by $V_{\max }\approx 19.0l_c^3$ and the maximum height is given by $h_{\max }\approx 2.2l_c$. For a curved, rotating substrate, the body force g is replaced by a centrifugal force that is a function the rotation rate $\omega$, radius from the axis of rotation $r\left(s\right)$, and the inclination angle of the substrate $\theta \left(s\right)$. For this case $l_c$ is modified to

$$l_c=\sqrt{\frac{\sigma }{\left(\mathit{\rho r}sin{\mathit{\theta \omega }}^2\right)}}.$$

If $h\left(s,t\right)$ is greater than $2.2l_c$, we assume that the drop has detached from the substrate and name this time as $t_{\max }$.

6. Applications

Approximately 400 billion two-piece, all-aluminum cans are produced annually for the purpose of storing beverages for distribution worldwide. The interior of each of these cans must be coated to protect the aluminum from onslaught due to corrosive elements in the contained beverage, and the beverage must be protected from picking up metal ions or other off-flavors from the aluminum substrate. Consequently the coating must be as uniform as possible for thick regions may slough off and thin regions may not offer adequate protection. To achieve a uniform film thickness, spin coating is employed using a spray fan to distribute the coating on the can substrate. But because the can is highly curved due to structural considerations, achieving a uniform final film thickness is much more complicated than for a flat substrate. In this section we will apply the analytical and numerical model we developed in previous sections to determine how the many parameters are influencing the flow of the paint coating. These parameters include the rotation rate, the shape of the can, the coating fluids physiochemical properties, and the geometry and flux of the spray fan, all affecting the final film thickness distribution.

The can body is initially punched from sheet aluminum and then goes through a washing process to produce a substrate suited for the spray coating. The can is then spun at between 2500 and 3500 rotations per minute, and one or two spray guns spray the interior with the liquid paint film. Centrifugal and gravitational forces redistribute this liquid layer as the can continues to spin after the initial spray process. The can is then placed in an oven where the solvent is allowed to evaporate leaving only the hardened resin. Then the can is filled with the beverage and the top of the can is attached in place. This conveyer process can produce as many as 1700 filled cans per minute.

In practice, one or two spray guns are used to coat the interior of the spinning cans. These are oriented at between 5° and 30° with respect to the vertical axis of the can and placed between 0.5 and 1.5 cm vertically from the top of the sidewall [6, 7]. Typically the can is sprayed for between 0.05 and 0.2 s and spun for an additional 0.1–0.5 s [6] so that centrifugal forces can act to redistribute the coating layer.

The industry uses schematic drawings which plot the substrate as a function of circular arcs of radius R_i, and subtended angles ${\phi }_{\mathrm{i}}$, as shown in Figure 2. Consequently the slope is always continuous, but there is a discontinuity in the radius of curvature between adjacent segments. Figure 6 shows a profile view of a typical can bottom, as well as a possible coating distribution, during the course of the fast spin with spray phase of the coating process. The functions $\kappa \left(s\right)$, $\theta \left(s\right)$, and $r\left(s\right)$ were derived from a published schematic drawing of a typical beverage can from a can manufacturer [8]. We define the “dome” of the can as the region from the centerline to the first inflection point, as labeled by the symbol A. The “moat” is defined as the region from A to D, with the inner moat defined as the region from A to B, and the outer portion of the moat from C to D. The intersection of the bottom of the can with the side wall occurs at E. In this study we will consider one phenomenon in the can coating process that can lead to a defect in the final coating layer. As centrifugal forces drive the coating outward along the dome, they can form a “front” which if it reaches the moat can form pendant drops on the inner wall of the moat. If centrifugal forces are sufficiently large, these drops can detach from the substrate and flung horizontally outward to the outer side of the moat. This can create a coating that is too thin on the inner wall of the moat and too thick on the opposite side (Figures 7 and 8).

Figure 6.

Profile of the can substrate in red lines and the coating layer in blue lines at $\widehat{t}=.10$, in the middle of fast spin with spray phase. Dimensionless variables are used.

Figure 7.

Profile of the can substrate in red lines and coating layer in blue at $\widehat{t}=.20$, in the middle of the fast spin phase. Dimensionless variables are used.

Figure 8.

Profile of the can substrate in red lines and coating layer in blue at $\widehat{t}={\widehat{t}}_{\max }=.45$, near the end of the fast spin phase. Dimensionless variables are used.

We will assume that the spray fan has an elliptical cross section with a ratio of the major axis to the length of the minor axis of 10. The parameters determining the placement and orientation of the spray gun are illustrated in Figure 9. For the simulation considered in this work, the nondimensional parameters are listed in Table 1, the dimensional parameters in Table 2, and typical properties of the coating liquid in Table 3.

Figure 9.

Profile view of can and spray gun placement and orientation parameters.

Parameter	Symbol	Value
Distance of spray gun from centerline: fan #1, fan #2	$\widehat{A}$	−0.15, −0.45
Distance of spray gun above can: fan #1, fan #2	$\widehat{B}$	0.15, 0.45
Angle of spray gun wrt vertical: fan #1, fan #2	$\vartheta$	28^°, 15^°
Subtended angle of spray fan	$\beta$	100°
Time spray gun acts in fast spin phase	${\widehat{t}}_{\text{spray}}$	0.15
Percent of secondary spray	$\zeta$	5%

Table 2.

Nondimensional parameters for this simulation.

Parameter	Symbol	Value
Can radius	R	3.33 cm
Rotation rate	$\omega$	2500–3500 RPM
Distance of spray gun from centerline: fan #1, fan #2	A	−0.5 cm, −1.5 cm
Distance of spray gun above can: fan #1, fan #2	B	0.5 cm, 1.5 cm
Average wet coating thickness	$h_{\mathrm{avg}}$	0.0028 cm
Time spray gun acts in fast spin phase	$t_{\text{spray}}$	0.05 s

Table 3.

Dimensional parameters for this simulation.

From Eq. (26), we see that the coating applied by the spray gun is an inverse function of the radius of the can substrate $r$. Consequently, the coating thickness near the centerline will be significantly greater than that far away from the centerline. This effect is illustrated in Figure 6 which shows a profile view of the can substrate and the coating thickness in the middle of the spray phase at time $\widehat{t}=0.10$. Centrifugal forces are mostly parallel to the substrate in the dome region and act to drive liquid outward toward the moat region, as shown in Figure 7. At the start of the inner wall of the moat, near symbol A in Figure 8, the can substrate becomes perpendicular to the centrifugal force, and the coating liquid tends to accumulate at this point. If the coating becomes too thick in this region, centrifugal forces can overcome surface tension forces and cause a pendant drop of liquid to form, which can detach from the inside of the dome edge and be thrown to the region between D and E. Surface tension forces act against the centrifugal forces in this region and reduce the chance that a pendant drop will detach from the coating. But if a droplet does detach, it will cause a starving of coating liquid in the inner wall region and an excess in the outer wall region, possibly leading to coating irregularities such as blisters after the final bake.

The dimensionless parameter $\rho R^2{\omega }^2/\left(\sigma /R\right)$ expresses the ratio of the pressure drop across the interface due to centrifugal forces versus surface tension forces. The dimensionless characteristic time for droplet detachment is plotted versus this dimensionless parameter in Figure 10. This graph estimates the maximum value of $\widehat{t}$, termed ${\widehat{t}}_{\max }$, before droplet ejection occurs for a Newtonian liquid. Here we estimate the value of $\widehat{t}$ when droplet detachment occurs using the criteria developed in Eq. (5) in Section 5. One can then determine how long the can may remain in the fast spin phase to avoid this coating irregularity.

Figure 10.

Critical time for droplet formation, ${\widehat{t}}_{\max }$, versus centrifugal force parameter. Here dimensionless variables are used.

The surface tension, density, and viscosity of the coating liquid are difficult to significantly alter as they depend on the required organic solvent content and surfactant levels in the paint formula. Similarly, R, the beverage can radius, is fixed by industry production standards. This leaves the rate of rotation as the only significant production parameter for changing the nondimensional (centrifugal force)/(surface tension force) parameter.

The centrifugal ejection of coating liquid from the inner wall of the moat is also predicted to be a strong function of the position and orientation of the spray gun. In the above example, the gun is placed 0.5 to the left of the centerline, is 0.5 cm above the top of the can, and is angled at ${28}^{\circ }$ with respect to the vertical. For our second simulation, the gun is placed 1.5 cm to the left of the centerline, 1.5 cm above the top of the can, and angled at ${15}^{\circ }$ with respect to the vertical. All other parameters are the same as in the previous example. This second spray gun position results in approximately 40% less spray falling on the dome and approximately 10% more coating on the sidewall than the first case. The total liquid volume emitted from the two spray guns is identical. Because the flux parallel to the substrate $Q^{\left(s\right)}$ is a function of the coating thickness cubed, we would expect that it would take longer before the coating is driven off the dome due to centrifugal forces for the second spray gun placement than the first spray gun placement. This is illustrated in Figure 11 where we plot ${\widehat{t}}_{\max }$ versus the dimensionless centrifugal force parameter for the two spray gun positions. With the new gun placement and orientation, one may keep the can in the fast spin phase for between 150 and 190% longer than the original spray gun configuration.

Figure 11.

Critical time for droplet formation, ${\widehat{\mathit{t}}}_{\mathit{max}}$, versus centrifugal force parameter for two different spray gun placements. For fan #1, the angle of inclination is ${\mathit{28}}^{\circ }$, ${\mathit{A}}_{\mathit{1}}=-\mathit{0.15}$, ${\mathit{B}}_{\mathit{1}}=\mathit{0.45}$. For fan #2, the angle of inclination is ${\mathit{15}}^{\circ }$, ${\mathit{A}}_{\mathit{2}}=-\mathit{0.45}$, ${\mathit{B}}_{\mathit{2}}=\mathit{0.45}$. Dimensionless variables are used.

7. Conclusion

In this work we have used scaling arguments and perturbation theory to derive the lubrication form of the governing fluid mechanical equations for a thin liquid film coating an arbitrarily curved, axisymmetric, rotating substrate. Our main purpose has been to develop mathematical model that can be employed to numerically simulate the application of a paint film to the interior of beverage cans, though the basic algorithm may be useful in other applications. We have used our algorithm to predict the time of centrifugal ejection of coating from the inner moat wall as a function of several input parameters: the physiochemical properties of the coating liquid, the rotation rate, and the spray gun placement. The model can also be used to predict other coating defects and how the input parameters can be used to avoid them. The effect of solvent evaporation during the drying phase, when gravity and surface tension forces affect the coating distribution as the viscosity increases until only a final, hard, film remains, may also be modeled. With so many parameters regulating the final coating thickness, using experiments to model the coating evolution and measure the final, dry, film thickness is an almost impossible task. Instead we can utilize the power and versatility of computer simulation to predict the coating profile as a function of input parameters. This understanding will be useful in optimizing the current application process. It may also be essential in acquiring a satisfactory coating when environmental regulations require a change to high solids and latex paints.