SYMMETRY REDUCTION AND NUMERICAL SOLUTION OF A THIRD-ORDER ODE FROM THIN FILM FLOW

A new approach to solving high-order ordinary differential equations numerically is presented. Instead of the usual approach of writing a high-order ordinary differential equation as a system of first-order ordinary differential equations, we write the high-order ordinary differential equation in terms of its differential invariants. The third-order ordinary differential equation y′′′ = y−k for constant k is used to illustrate this approach for the cases k = 2 and k = 3. KeywordsThin film; third-order ODE; symmetry reduction.


INTRODUCTION
We consider the ordinary differential equation (ODE) solved subject to the initial conditions (α, β, γ and λ are constants) Analytical solutions to the third-order ODE (1) are given in general in the handbook of Polyanin and Zaitsev [17] for the cases k = 0, 7/2, 5/2, 4/3, 7/6, 2, 1/2 and 5/4.The ordinary differential equation for the case k = 2 has been derived by Tanner [18] to investigate the motion of the contact line for a thin oil drop spreading on a horizontal surface.A review of thin film flow has been presented by Myers [14] in which some third-order ODEs occuring in thin film flow are discussed.The interested reader is also referred to Tuck and Schwartz [19] and the references therein.Analytical solutions admitted by the case k = 2 has also recently been investigated by Duffy and Wilson [4].The analytical solution for the case k = 2 was initially presented by Ford [5].The case k = 2, amongst others, has been the subject of intense numerical investigation by Tuck and Schwartz [19].Existence and uniqueness of solutions admitted by (1) are discussed in Troy [21] and Bernis [1].Momoniat [11] has investigated autonomous integrals of (1).A numerical investigation of boundary value problems associated with (1) have recently been considered by Momoniat [12].
In this paper we consider the equation ( 1) for arbitrary k.We demonstrate our proposed algorithm by considering the cases k = 2 and k = 3.The case k = 3, which has not been investigated before for analytical solutions, can be used to model the draining down of a dry wall of a non-Newtonian thin film.Alternatively, it can represent the steady-state solution to the nonlinear partial differential equation which describes the surface tension driven spreading of a thin viscous film (see e.g.Myers [14]).The results obtained here are general and can be applied to high-order initial value problems that admit Lie point symmetries.

Point Symmetries
We apply the Lie group method to the ODE (1) in general before considering the specific case k = 3 which has not been investigated before analytically.
The interested reader is referred to the books by Bluman and Kumei [3], Ibragimov [9], Olver [15] and Ovsiannikov [16] on the Lie group method for solving differential equations.
The operator generates a point symmetry of (1) if the determining equation X [3] ( holds, where | y ′′′ =y −k refers to an evaluation on equation (1).Also, the third extension or prolongation of the generator X is given by where with y (j) = d j y/dx j , ζ 0 = η(x, y) and d/dx is the operator of total differentiation given by By use of the approach of Lie one can calculate ξ(x, y) and η(x, y) systematically.The determining equation ( 5) can clearly be written as where ζ 3 is given in (7).The expansion of ( 9) by use of (7) (and its expansion in terms of ξ and η) results in The determining equation ( 10) can be separated by coefficients of derivatives of y.Solutions of the resultant system of linear partial differential equations are where a 0 and a 1 are arbitrary constants.Setting firstly a 1 = 0 and then a 0 = 0 we obtain the two Lie point symmetry generators Lie point symmetries admitted by differential equations are easily obtainable by the use of computer algebra packages (see e.g., Head [7], Sherring et al. [20] and Baumann [2]).

Group Invariant Solutions
Here we deduce a group invariant solution admitted by (1).We show that the initial conditions (2) are satisfied only if they are in terms of α.
We consider a linear combination of the Lie point symmetries (12) of the form where a 1 and a 2 are constants to be determined.We determine a group invariant solution y = Γ(x) admitted by ( 1) by solving the first-order ODE obtained from The solution of ( 14) yields where a 3 is a constant of integration.The substitution of ( 15) into (1) gives Without loss of generality we can choose a 2 = 1 otherwise we get the trivial solution.
We determine the constant a 1 by imposition of the initial condition to find that We then calculate the values for β and γ in terms of α.The invocation of ( 2) on (15) gives rise to For k = −1, from ( 14) we obtain y = a 3 exp(3a 2 x/a 1 ), where a 3 is a constant, which does not give rise to an invariant solution.Specifically, when k = 3 we find that To prevent the solution (20) from becoming imaginary, the following condition needs to be maintained From (19), for k = 3, we find that The group invariant solution (15) is not valid for the cases k = 2 and k = 1/2.For both the cases k = 2 and k = 1/2 we are unable to determine group invariant solutions.In the case k = 2 the exponent in (15) becomes one.On substitution into (1) we are unable to satisfy the right-hand side without the solution becoming singular.The same holds for k = 1/2 where the exponent in (15) becomes two.

Group Reduction
The algebra formed by the Lie point symmetries is two-dimensional since Hence we firstly perform a reduction on (1) using X 1 .Then X 2 is a Lie point symmetry generator of the reduced equation.An invariant corresponding to X 1 is given by We can use this invariant as follows The initial conditions (2) then imply that Equation ( 1) can be reduced to The Lie point symmetry generator X 2 is admitted by (28) and is written in terms of the new variables as The first prolongation of ( 29) is given by Invariants corresponding to (30) are given by The second-order ordinary differential equation (28) can then be reduced to the first-order ordinary differential equation For the case k = 3 the ODE (32) is solved subject to the initial conditions w ( α For the case k = 3 from (31), ( 24) and (25) we have that The first-order reduction (32) is not valid for the case k = 2.For this particular case the first prolongation of the Lie point symmetry generator (30) is Invariants corresponding to (35) are given by Equation ( 28) for k = 2 can then be written as the following first-order ODE in terms of the invariants q and h: The first-order ODE (37) is solved subject to Then from (36), ( 24) and (25) we have that We utilise the above results in the next section.

NUMERICAL STUDY
We have reduced the third-order ODE y ′′′ = y −k by successive reduction of order to a second-order ODE (28) and then to a first-order ODE (32) or (37).Unfortunately, for general k, (32) cannot be solved analytically.We can however use these reductions to determine an efficient way to solve (1) numerically.We focus on the cases k = 2 and k = 3.
A fourth-order Runge-Kutta method to solve the first-order ODE y ′ = f (x, y) gives values y i = y(x i ) where x i = λ + ih and h is the step length.An initial value y 0 = y(0) is chosen with successive values determined by where We consider the initial conditions i.e., we have taken To use a fourth-order Runge-Kutta method we write (1) as a system of three first-order equations: where The results are displayed in Table 1 for the case k = 2 and Table 2 for the case k = 3.From Duffy and Wilson [4] the exact solution for the case y ′′′ = y −2 is given in parametric form as where Ai(s) and Bi(s) denote the Airy functions and a, b and s 0 are constants to be determined.Imposing the initial conditions (42) we find that The results for the exact solution are presented in Table 1 for a comparison.The numerical solution for the cases k = 2 and k = 3 can be achieved in two ways.Firstly, the second-order ordinary differential equation ( 28) is solved numerically using a fourth-order Runge-Kutta method subject to (27).We do this by writing (28) as a system of two equations The conditions (27) imply that We form a function v(y) using cubic interpolation in MATHEMATICA [22].The first-order ODE (24) is solved numerically using a fourth-order Runge-Kutta method.These results are displayed in Table 1 for k = 2 and Table 2 for k = 3 in the column labelled Method 1.
A second method which can be employed due to the double reduction of ( 1) is to solve the equivalent first-order ordinary differential equations.For the case k = 2 we solve the first-order ODE (37) subject to (38) using the fourth-order Runge-Kutta method (41).The function z = z(q) is formed via interpolation in MATHEMATICA [22].The second-order ODE (39) is then solved using (41).The results using this second approach for k = 2 are indicated in Table 1 under the column labelled Method 2.
For k = 3 we solve firstly (32) subject to (33) using (41).We again interpolate the function w(p) using MATHEMATICA [22].The second-order ODE (34) is then solved using (41).These results are displayed in Table 2 From Table 1 we observe that the numerical results using Method 1 are correct to 5 decimal places.The numerical results from Method 2 are only correct to 2 decimal places.Applying the fourth-order Runge-Kutta method to (1) for k = 2 also yields five decimal place accuracy.
From Table 2 we observe that the numerical results using the fourth-order Runge-Kutta and Method 1 agree to 6 decimal places.Method 2 generates results correct to 2 decimal places.This is consistent with the results displayed in Table 1.

CONCLUDING REMARKS
The approach discussed in this paper is useful for validation of numerical results when no analytical results are known.The advantage in using this approach is that the same numerical algorithm is used, i.e. fourth-order Runge-Kutta.We note that for (1) we have obtained two validation approaches because of the double reduction effected by the occurence of two Lie point symmetries.The first approach, Method 1, is the most accurate of the two.The reduced second-order ODE (28) has a nonlinearity at y = 0 which is easily avoided.The second approach, Method 2, requires the solution of two highly nonlinear first-order ODEs (32) and (37).
The approach taken in this paper can easily be extended to the third-order ODEs y ′′′ = f (y) considered by Tuck and Schwartz [19] where

Table 1 :
Table comparing numerical values of the exact solution, a fourth-order

Table 2 :
Table comparing numerical values of the numerical solution obtained using a fourth-order Runge-Kutta method and our new approach at x ∈ [0, 0.2, 0.4, 0.6, 0.8, 1.0] taking h = 0.01 and k = 3 for the initial conditions y