Weakly Nonlinear Hydrodynamic Stability of the Thin Newtonian Fluid Flowing on a Rotating Circular Disk

The main object of this paper is to study the weakly nonlinear hydrodynamic stability of the thin Newtonian fluid flowing on a rotating circular disk. A long-wave perturbation method is used to derive the nonlinear evolution equation for the film flow. The linear behaviors of the spreading wave are investigated by normal mode approach, and its weakly nonlinear behaviors are explored by the method of multiple scales. The Ginzburg-Landau equation is determined to discuss the necessary condition for the existence of such flow pattern. The results indicate that the superctitical instability region increases, and the subcritical stability region decreases with the increase of the rotation number or the radius of circular disk. It is found that the rotation number and the radius of circular disk not only play the significant roles in destabilizing the flow in the linear stability analysis but also shrink the area of supercritical stability region at high Reynolds number in the weakly nonlinear stability analysis.


Introduction
The study of the hydrodynamic stability of a thin liquid film is important for a wide range of situations, varying from engineering science and chemical science.Due to various applications, attention is gradually focused on this subject.For instance, a process for coating a surface on a spinning substrate, wherein a liquid coating material is dispensed radially from the center to the edge or from the edge to the center above the surface, and after application of the coating material the coating is cured.This process usually referred to as spin coating.To control the uniform and stable thin film is an interesting subject in the technological development for photolithography in wafer manufacturing 1, 2 .
The linear stability theories for various film flows have been clearly presented by Lin 3 and Chandrasekhar 4 .Kapitza  over inclined planes.Various stability behaviors of the layer flows were analyzed.Benney 6 studied the nonlinear evolution equation for free surfaces of the film flows by using the method of small parameters.The theory of Kapitza was later found not so compatible with the nonlinear theory of Landau 7 .The same film flow stability problem was studied by Yih 8 using a numerical approach.The transition mechanism from laminar flow to turbulent flow was elegantly explained by the Landau equation.This study sheds a light later for further development on the theory of nonlinear film stability.The Landau equation was rederived by Stuart 9 using the disturbed energy balance equation and Reynolds stresses.Pumir et al. 10 further included the effect of surface tension on the film flow model and solved for the solitary wave solutions.Hwang and Weng 11 showed that the conditions of both supercritical stability and subcritical instability are possible to occur for a film flow system.Atherton and Homsy 12 discussed the derivation of complicated nonlinear partial equations, evolution equations that describe the movement of a fluid-fluid interface.Ruyer-Quil and Manneville 13-15 derived several models of the film flows down an inclined plane by the numerical simulation.The results indicate that it allows to better capture the viscous effects which are dominant at small Reynolds numbers.Amaouche et al. 16 showed an accurate modeling of a wavy film flow down an inclined plane in the inertia-dominated regimes by using the weighted residual technique.It is shown that the model follows quite closely, for a suitable choice of α, the Orr-Sommerfeld equation for all Weber and Kapitza numbers, in linear stability analysis.Samanta 17 derived the wave solution of a viscous film flowing down on a vertical nonuniformly heated wall.The results indicate that supercritical unstable region increases, and the subcritical stable region decreases with the increase of Peclet number.
Emslie et al. 18 were the pioneers who analyzed a Newtonian liquid flowing on rotating disk.The flow is governed by a balance between the centrifugal force and the viscous resisting force.It was shown that the nonuniform distribution in the initial film profile tends to become uniform during spinning.This model has been widely employed in the subsequent investigations.Higgins 19 analyzed the flow of a Newtonian liquid placed on an impulsively started rotating disk.In this work, a uniform film thickness is assumed, and the method of matched asymptotic expansions is adopted.It is showed that for films that maintains a planar interface there exists a self-similar form for the velocity field that allows the radial dependence to be factored out of the Navier-Stokes equations and boundary conditions.Kitamura et al. 20 solved the unsteady thin liquid film flow of nonuniform thickness on a rotating disk by asymptotic methods.
In this paper the authors are interested in investigating the weakly nonlinear hydrodynamic stability of the thin Newtonian liquid flowing on a rotating circular disk.It is assumed that the disk radius is much larger than film thickness.Therefore, the peripheral effects are neglected by comparing with total film area 20 .It is focused that the effects of instability due to inertia and centrifugal forces were revealed in the region near the rotating axis.The influence of the rotating motion and the disk size effect on the equilibrium finite amplitude is studied and characterized mathematically.In an attempt to verify computational results and to illustrate the effectiveness of the proposed modeling approach, several numerical examples are also presented.

Mathematical Formulation
Consider a two-dimensional incompressible, viscous liquid film flowing on a rotating circular disk which rotates with constant velocity Ω * see Figure 1 .The variable with a superscript " * stands for a dimensional quantity.Here the cylindrical polar coordinate axes r * , θ * , and z * are chosen as the radial direction, the circumferential direction, and the axial direction, respectively.All associated physical properties and the rate of film flow are assumed to be constant i.e., time-invariant .Let u * and w * be the velocity components in the radial direction r * and the axial direction z * .The governing equations of motion are where v * is the tangential velocity, ρ is the constant fluid density, p * is the fluid pressure, g is the acceleration due to gravity, and μ is the fluid dynamic viscosity.On the disk surface z * 0, the no-slip boundary conditions for the velocity fields are u * 0,

2.2
On the free surface z * h * , the boundary condition approximated by the vanishing of shear stress is expressed as The normal stress condition obtained by solving the balance equation in the direction normal to the free surface is given as

2.4
where h * is the local film thickness, S * is surface tension, and p * a is the atmosphere pressure.The kinematic condition that the flow does not travel across a free surface can be given as By introducing a stream function ϕ * , the dimensional velocity components can be expressed as Now the following variables are used to form the dimensionless governing equations and boundary conditions:

2.7
where h * 0 , α, u * 0 , Re, ν, Fr, and λ are the average film thickness, perturbed wavelength, scale of velocity, Reynolds number, the kinematic viscosity, Froude number, and dimensionless wavenumber, respectively.For simplification, it is assumed that Newtonian liquid film is very thin h * r * .In consequence, it is reasonable to assume that the tangential velocity is constant throughout the radial direction in the thin film, that is, v * r * Ω * .In order to investigate the effect of angular velocity, Ω * , on the stability of the flow field, the dimensionless rotation number is introduced as follows: In terms of these nondimensional variables, the equations of motion become

2.9
Using the nondimensional variables, the boundary conditions at the surface of disk z 0 are reduced to ϕ ϕ r ϕ z 0.

2.10
And the boundary conditions at the free surface of disk z h become 12 Hence the term α 2 S can be treated as a quantity of zeroth order 21 .Since the modes of longwavelength that give the smallest wave number are most likely to induce flow instability for the film flow, this can be done by expanding the stream function and flow pressure in terms of some small wave number α 1 as

2.14
Following procedure described in 6 , the complicated nonlinear system of 2.9 and boundary conditions 2.10 -2.13 is reduced to a single nonlinear evolution equation for the film thickness h r, t .After inserting the above two equations into the system of equations and boundary conditions, the zeroth order α 0 terms in the governing equations can be expressed as

2.15
The boundary conditions associated with the equations of zeroth order are given as z 0,

Mathematical Problems in Engineering
The solutions of the zeroth-order equations are

2.17
After considering all terms of first order α 1 from the system of equations and boundary conditions, the first-order equations can be achieved as and boundary conditions can be achieved as z 0,

2.19
The solutions of the first-order equations are given as

2.20
By substituting the solutions of the zeroth-order and first-order equations into the dimensionless-free surface kinematic equation of 2.13 , the nonlinear evolution equation is derived

and expressed as h t A h h r B h h rr C h h rrr D h h rrrr E h h 2 r F h h r h rrr 0, 2.21
where

Stability Analysis
The dimensionless film thickness when expressed in perturbed state can be given as h r, t 1 η r, t , where η is a perturbed quantity to the stationary film thickness.Substituting the value of h r, t into the evolution equation 2.21 and all terms up to the order of η 3 are collected, the evolution equation of η becomes where the values ofA, B, C, D, E, F, and their derivatives are all evaluated at the dimensionless height of the film h 1.

Linear Stability Analysis
As the nonlinear terms of 3.2 are neglected, the linearized equation is given as η t Aη r Bη rr Cη rrr Dη rrrr 0.

Mathematical Problems in Engineering
In order to use the normal mode analysis, we assume that η a exp i r − dt c.c., 3.4 where a is the perturbation amplitude, and c.c. is the complex conjugate counterpart.The complex wave celerity, d, is given as where d r and d i are regarded as the linear wave speed and linear growth rate of the disturbance, respectively.The solution of the disturbance about h r, t 1 is asymptotically stable or unstable according to d i < 0 or d i > 0. This is equivalent to the inequality of B < D or B > D.

Weakly Nonlinear Stability Analysis
Nonlinear effects, when they are weak enough, do not fundamentally alter the nature of the motion.A weakly nonlinear solution can still be usefully expressed as a superposition of plane waves, but the amplitudes of these waves do not remain constant; they are modulated by nonlinear interactions.In order to characterize the weakly nonlinear behaviors of thin film flows, the method of multiple scales 22 is employed here, and the resulting Ginburg-Landau equation 23 can be derived as where

3.7
Equation 3.6 can be used to investigate the weak nonlinear behavior of the fluid film flow.
In order to solve for 3.6 , solution is taken for a filtered wave in which spatial modulation does not exist.So for a filtered wave a can be given as a a 0 exp − ib t 2 t 2 .Subcritical instability

3.8
After substituting 3.8 into 3.6 , one can obtain The threshold amplitude εa 0 in the supercritical stable region is given as and the nonlinear wave speed Nc r is given as The detail derivation of the earlier mentioned equations can refer to Cheng and Lai 24 .
If E 1 0, 3.9 is reduced to a linear equation.The second term on the right-hand side of 3.9 is due to the nonlinearity and may moderate or accelerate the exponential growth of the linear disturbance according to the signs of d i and E 1 .Equation 3.10 is used to modify the perturbed wave speed caused by infinitesimal disturbances appearing in the nonlinear system.The condition for the film flow to present the behavior of subcritical instability in the linearly stable region d i < 0 is given as E 1 < 0, and the threshold amplitude of the wave is given as εa 0 .The subcritical stable region can only be found as E 1 > 0. The neutral stability curve can only be derived and plotted for the condition of E 1 0. On the basis of earlier mentioned discussion, it is obvious that the Ginzburg-Landau equation can be used to characterize various flow states.The Landau equation can be summarized and presented in Table 1.Re Re Re an increase of the radius of circular disk.Hence one can say that in linear stability analysis rotation number and the radius of circular disk give the same destabilizing effects.

Weakly Nonlinear Stability Analysis
The main purpose of the nonlinear stability analysis is to study the weakly nonlinear analysis of the evolution equation.Figures 3 a to 3 c indicates that the area of shaded subcritical instability region and subcritical stability region decreases, and the area of shaded supercritical instability region increases with an increase of rotation number or the radius of circular disk.Also, from these figures, it is interesting to find, in high Reynolds regime, that the area of supercritical stability region decreases with an increase of rotation number or an increase of circular disk radius.It is found that the rotation number and the radius of circular disk not only play the significant roles in destabilizing the flow in the linear stability analysis but also decrease the ranges of supercritical stability region at large Reynolds number in  the weakly nonlinear stability analysis.The reason for this phenomenon is the existence of the centrifugal force term, which is a radius-related force in the governing equation.As the size of radius and Reynolds number gradually increased, the centrifugal force is enhanced significantly.Due to inertia forces, that may accelerate the growth of the linear disturbance, the trend of instability for the flow with larger radius and larger Reynolds number is higher than those with smaller ones.Figure 4 a shows the threshold amplitude in subcritical instability region for various wave numbers with different Ro values at Re 6 and r 50.The results indicate that the threshold amplitude εa 0 becomes smaller as the value of rotation number increases.Figure 4 b shows the threshold amplitude in subcritical instability region for various wave numbers with different r values at Re 6 and Ro 0.1.The results indicate that the threshold amplitude εa 0 becomes smaller as the value of the radius increases.In such situations, the film flow which holds the higher threshold amplitude value will become more stable than that which holds smaller one.If the initial finite amplitude disturbance is less than the threshold amplitude, the system will become conditionally stable.
Figure 5 a shows the threshold amplitude in the supercritical stability region for various wave numbers with different Ro values at Re 6 and r 50. Figure 5 b shows the threshold amplitude in the supercritical stability region for various wave numbers with different r values at Re 6 and Ro 0.1.It is found that the decrease of rotation number or the radius of circular disk will lower the threshold amplitude, and the flow will become relatively more stable.
The wave speed of 3.5 predicted by using the linear theory is a constant value for all wave number and rotation number.However, the wave speed of 3.12 predicted by using nonlinear theory is no longer a constant.It is actually a function of wave number, Reynolds number, rotation number, and the radius of disk.

Figure 1 :
Figure 1: Schematic diagram of thin Newtonian liquid flow upon a rotating circular disk.

Figure 3 :
Figure 3: a Neutral stability curves of Newtonian film flows for Ro 0.1 and r 50.b Neutral stability curves of Newtonian film flows for Ro 0.15 and r 50.c Neutral stability curves of Newtonian film flows for Ro 0.1 and r 100.

Figure 4 :
Figure 4: a Threshold amplitude in subcritical instability region for three different Ro values at Re 6 and r 50.b Threshold amplitude in subcritical instability region for three different r values at Re 6 and Ro 0.1.

Figure 6 aFigure 5 :Figure 6 :
Figure 5: a Threshold amplitude in supercritical stability region for two different Ro values at Re 6 and r 50.b Threshold amplitude in supercritical stability region for two different r values at Re 6 and Ro 0.1.

Table 1 :
Various states of Landau equation.