On implementation of a semi-analytic strategy to develop an analytical solution of a steady-state isothermal tube drawing model.

In this paper, we consider an isothermal glass tube drawing model consisting of three coupled nonlinear partial differential equations. The steady-state solution of this model is required in order to investigate its stability. With the given initial and boundary conditions, it is not possible to determine an analytical solution of this model. The difficulty lies in determining the constants of integrations while solving the second order ordinary differential equation analytically appearing in the steady-state model. To overcome this difficulty, we present a numerical based approach for the first time to develop an analytical solution of the steady-state isothermal tube drawing model. We use a numerical technique called shooting method to convert the boundary value problem into a set of initial value problems. Once the model has been converted into a system of differential equations with initial values, an integrating technique is implemented to develop the analytical solution. The computed analytical solution is then compared with the numerical solution to better understand the accuracy of obtained solution with necessary discussions.

Introduction

In glass industry, tubes are drawn through various manufacturing processes used to achieve the continuous production of glass tubes having correct wall thickness and diameter. The most commonly used are the Danner process and the Vello process[1-5] having great importance in glass fabricating industry and are still in use today. In Danner process, glass is melted in a furnace to the stage where it is soft and pliable. Molten glass is then let to fall with low feeding speed on the surface of a cylindrical device called mandrel kept in a temperature controlled tank called oven. Mandrel is slightly inclined and hollow such that the air can be blown through it. By continuous rotation of mandrel about its axis of symmetry, the molten glass falling downward creates a smooth layer around the mandrel. It cools down gradually and takes the shape of a thick-walled hollow glass tube with all desired properties at just below the end of mandrel. The length of hot-forming zone is taken as L. It is then pulled out by a drawing machine with a drawing speed . This ratio is called the draw ratio. Keeping a constant temperature in the hot-forming zone leads to develop an isothermal tube drawing model. The drawn tube is then conveyed straight by rollers to further process of cutting, finishing, polishing and packaging at the end of the spinline. This manufacturing process is explained and illustrated in[3]. For the Vello process, we refer to[1-5].

The geometry of the drawn tube is illustrated in Fig. 1. The length of hot forming zone is denoted by L. Inner and outer radii, inside pressure, feeding and take up speeds, axis of symmetry of the tube during the production process are shown in Figure. All other parameters involved along with their numerical values are illustrated in Table 1.

Figure 1

Diagram of a glass tube during production process.

Table 1

Summary of parametric values appearing in the isothermal model (5).

Parameter	Symbol	Approximate value	Units
Feeding speed	v_0	1	mm/s
Drawing speed	v_L	12	mm/s
Length of the hot-forming zone	L	1	m
Input viscosity	μ_0	5×10^5	Pa s
Inside pressure	p_s	420	Pa
Density	ρ	2500	kg/m^3
Mean radius of the glass tube	R_0	30	mm
Initial area of the tube	A_0	1885	mm

The shaping parameters such as the wall thickness and cross-sectional area (or diameter) are the main characterizations of the drawn tube. In either of the manufacturing processes, the required shape of the tube can be maintained by the stream of the air gently blown through the mandrel. Insufficient quantity of the air blown through the mandrel can disturb the desired shape of the glass tube. As a result, this insufficient balance collapse the walls of the glass tube. Moreover, the geometry of the glass tube can be controlled by the parameters involved herein such as the glass temperature, the composition of the raw material, the pressure of the blowing air in the mandrel and the rate of draw. We have added a Fig. 1 to better understand the geometry of the drawn glass tube.

A variety of mathematical models both for the isothermal and non-isothermal tube drawing processes, with different levels of descriptions and needs, are available in the literature e.g., see[1-4] and[6-11]. In[1-3,7] numerical solutions of both the isothermal and non-isothermal models have been found and used to optimally control the geometry of the glass tube. In recent years, a variety of mathematical models representing physical and real world problems have been investigated for numerical solutions and stability analysis (for example see[12-24]). Thus, the concept of finding numerically accurate and exact solutions of real world problems has attracted attention from all over the world.

In this paper, we have considered an isothermal glass tube drawing model consisting of three coupled nonlinear partial differential equations of first and second order. The steady-state numerical solution of this model is required in order to investigate its stability. Recently, we have analyzed the stability of an isothermal tube drawing model and performed a complete mathematical analysis in[3] by incorporating a steady state numerical solution. To better understand the physical model and to examine the accuracy of the obtained numerical solution, our objective here is to develop an analytical solution of the steady-state isothermal model. To start with, it is not possible to determine an analytical solution of this type of model with the given initial and boundary conditions. When we solve the second order ordinary differential equation appearing in the steady-state model, the difficulty lies in determining the constants of integrations. To overcome this concern, we utilize a numerical technique called shooting method to convert the boundary value problem into a set of initial value problems. An integrating technique is implemented to develop an analytical solution of the converted initial value problems. Similar approach has been suggested in[10,22] to construct an analytical solution of a steady-state melt-spinning model. The semi-analytical approach developed in this paper has an advantage over other numerical methods that it may be applied directly to both linear and nonlinear systems of equations without the need for discretization, linearization, or perturbation.

The organization of the paper is as follows: In “Governing equations” section, we give a brief description of the mathematical model of an isothermal tube drawing process. A strategy converting a boundary value problem into a set of initial value problems is described in “Numerical approach” section. Analytical solution is developed in “Analytical solution” section and compared with the numerical solution in “Accuracy comparison” section. Conclusion is given in “Conclusions” section. We have included some leftover parts of the modeling of tube drawing process in the “Appendix”. In this section, we have included one example concerning the application of the analytical solution obtained through the implemented semi-analytic technique.

Governing equations

In the literature, different types of models for the drawing processes with different level of demands and descriptions are available. A considerable amount of work has been carried out by different researchers[1-4,6,8,9] and[11,25-28] to model the tube drawing process. In this section, we briefly explain the mathematical model for an isothermal tube drawing process.

To model the tube drawing process, we consider an incompressible Newtonian flow of a molten glass between two free surfaces and where and respectively denote the inner and outer radii of the glass tube, and assume that temperature remains constant throughout the forming zone. Glass tube during production process is illustrated in Fig. 1. In the draw-down zone, the surface tension force and the inertial force acting upon the molten glass are insignificant and hence can be neglected. This kind of flow is governed by the equations

The Eqs. (1a)–(1c) are taken from[4] and are known as the standard equations showing the axi-symmetric stokes flow. The first equation is the continuity equation and the last two equations are the momentum equations respectively in r and z directions. We denote the derivatives by subscripts r and z where z denotes the distance along the axis of the glass tube and r measures the distance perpendicular to it. The velocity of the molten glass is defined to be where u and v are the components of velocity along z and r direction respectively. The pressure, density and the acceleration due to gravity are denoted by p, and g respectively.

Now, at the free surfaces and , it is necessary to specify the stress conditions and the kinematic conditions.

On the inner and outer surfaces of the glass tube, the stress conditions are given as:where and are the unit normals on the surfaces and of the tube respectively defined asand is the inside pressure applied on the surface of the glass tube and is the stress tensor given as

The kinematic conditions are described asand thus the stress conditions (1d)–(1e) can be expanded to give where denotes the viscosity of the molten glass which remains constant at the given temperature.

To take benefit of the small parameters present in the problem, it is now convenient to convert the Eq. (1) into dimensionless form. For the dimensional quantities, the appropriate scaling is defined as:where . W is the width of the glass tube that has a very small value than the typical length of hot-forming zone L, U is the typical drawing speed and is the typical melt glass viscosity which has a constant value for the isothermal case. After dropping the bar notation, the system of governing equations in dimensionless form is given as follows: where is known as Stokes number.

The stress conditions (1h)–(1k) read as

The kinematics conditions (1f)–(1g) become

The above model can be simplified by means of an asymptotic expansion, in which the inverse aspect ratio is used as scaling parameter. Assuming the glass flow as a thin layer flow and ignoring the large aspect ratio of the flow, one derives the simplified equations to model the isothermal tube drawing process. In this derivation, the surface tension and inertial forces acting upon the molten glass have also been ignored due to their insignificant contributions. For details, we refer to[1-6] and[9,27]. However, we have included the remaining considerations of the modeling process in “Appendix A” and reach at the following model equations representing the isothermal process: equipped with the initial conditionsand the boundary conditionswhere and are the cross-sectional area and average radius of the glass tube at the time of entering the hot-forming zone, respectively. Acceleration due to gravity and density of the molten glass are denoted by g and , respectively. Average radius R of the tube is defined asEquations in system (3) give us cross-sectional area A, velocity v, and average radius R of the glass tube. Width W of the tube can be determined by the equation

Since the temperature throughout the forming zone of the tubing process remains constant, the viscosity of the melt glass also remains constant.

Dimensionless form

We introduce the following dimensionless quantitiesinto the model (3) to get the dimensionless equations where we have removed the asterisk notation andare the dimensionless parameters.

The initial and boundary conditions are respectively transformed to where is the draw ratio.

The steady state form of the model (4) is written as subject to the conditions:

The system (5) is a system of coupled nonlinear ordinary differential equations. We are interested in developing an analytical solution of the nonlinear system (5). However, we have already proved the existence and uniqueness of the solution of steady state isothermal tube drawing model (5) in[3]. For details, we refer the reader to[3].

Numerical approach

There are different techniques available in the literature to solve boundary value problems analytically. The obtained general solutions contains arbitrary constant which can be determined by using boundary conditions. However, there exist some problems (e.g. problem (5) and see[10,22]) where it is not possible to find integration constants with the given conditions. It is interesting to mention that integration constants can be found if boundary value problem (BVP) is transformed to an initial value problem (IVP).

Shooting method is a technique that transforms a BVP of the forminto an IVPwhere s is a guess to be determined in such a way that it shoots .

We apply RK4 method to the second order differential Eq. (6) to determine approximation for . Let the first guess be taken asRK4 method with gives us an approximation for . If the absolute error is less than some pre-assigned tolerance, we stop, otherwise we refine our guess by consideringand compute by RK4 method. If the error agrees to pre-assigned tolerance, we stop. In the case of disagreement, we make further guesses using the secant formulawhereis the error of approximation. RK4 method is continued each time with a new guess computed from (7) until agrees with the value to pre-set tolerance.

Once we able to compute a most suitable guess , we are able to transform the BVP to IVP, and hence an analytical solution of the problem of type given in (6) can be found.

Analytical solution

With the initial conditions (5d), the Eq. (5a) is solved to get

Equation (5b) can now be put in the form along with the conditions

With the usual techniques, it is not an easy task to develop an analytical solution of the BVP (9). We use a non-conventional approach to build an analytical solution of the given nonlinear model. We convert the BVP (9) into two IVPs. For this purpose, we setand

Then the IVPs corresponding to BVP (9) are defined as along with the conditions where is approximated value determined by shooting the value of v at 1 i.e., . (The process of finding the approximation (say) by shooting method is explained in “Numerical approach” section).

The solution of first order non-linear ordinary differential Eq. (11b) is determined to givewhere is a constant of integration and is computed to give

Solution (12) can be further re-written as

Separating the variables and then integrating, we reach atwhere is another constant of integration.

Solving for v, we obtain

Equation (15) in view of condition (5d) gives us the value of , i.e.

Equation (16) along with and respectively given in (13) and (17) represents an analytical solution of Eq. (5b).

Next we put the steady state solutions for A and v in Eq. (5c) to getwhere .

Equation (18) is a well-known Riccati’s equation[29] whose solution is supposed to be of the formwhere one solution , guessed by hit and trial procedure, is given asand y is any nonzero unknown function of z.

Differentiating (19) with respect to z and then putting the values of u and in Eq. (18), we reach atwhich is a first order linear ordinary differential equation in y. After few steps of calculations, the solution of this differential equation is found to bewhere is a constant of integration, and A(z), B(z) are given asUsing (19) and the condition (5d), we getwhereandSubstituting in (20), we obtainTherefore,whereand henceis an analytical solution of the Eq. (5c) along with and given respectively in (13) and (17). Thus, the solutions given in (8), (16) and (22) together constitute an analytical solution of the steady-state model (5).

Accuracy comparison

The developed analytical solution is plotted and compared with the numerical solution in Figs. 2, 3, 4. The parametric values used in the model are given in Table 1. The solutions obtained by RK4 method are reliable and competitive with the analytical solutions. We observe that the errors between numerical and analytical solutions for the steady state variables A, v and R are negligibly small (see Table 2). We also determine W, the width of the glass tube, by using the relation . The steady state analytical and numerical solutions for W and the corresponding error are also shown in the Fig. 5. It is observed from Table 2 that error between the numerical and analytical solution for A, v, and R is maximum at , , and respectively. Moreover, the error for the solution W is maximum at .

Figure 2

Analytical and numerical solutions for the steady state variable A and the corresponding errors between the solutions with discretization step size .

Figure 3

Analytical and numerical solutions for the steady state variable v and the corresponding errors between the solutions with discretization step size .

Figure 4

Analytical and numerical solutions for the steady state variable R and the corresponding errors between the solutions with discretization step size .

Table 2

Summary of errors in solutions for using step size .

z	Error in A	Error in v	Error in R	Error in W
0.0	0.00	0.00	0.00	0.00
0.1	1.08×10^-6	1.70×10^-7	1.16×10^-2	1.60×10^-6
0.2	1.21×10^-6	5.40×10^-7	1.07×10^-2	1.40×10^-6
0.3	9.67×10^-7	1.37×10^-6	9.90×10^-3	1.20×10^-6
0.4	7.89×10^-7	2.82×10^-6	9.20×10^-3	1.10×10^-6
0.5	6.72×10^-7	5.01×10^-6	8.70×10^-3	1.00×10^-6
0.6	5.92×10^-7	8.12×10^-6	8.20×10^-3	9.00×10^-7
0.7	5.37×10^-7	1.23×10^-5	7.80×10^-3	9.00×10^-7
0.8	4.99×10^-7	1.77×10^-5	7.50×10^-3	8.00×10^-7
0.9	4.72×10^-7	2.45×10^-5	7.20×10^-3	8.00×10^-7
1.0	4.53×10^-7	3.29×10^-5	6.90×10^-3	7.00×10^-7

Figure 5

Analytical and numerical solutions for the steady state variable W and the corresponding errors between the solutions with discretization step size .

Conclusions

In this work, we have presented a numerical based approach to find analytical solution of a boundary value problem where one cannot find constants of integration with the given boundary conditions. The approach is to convert a boundary value problem into initial value problems and then to develop an analytical solution of the resulting problems with the usual methods. To explain the approach, we have developed an analytical solution of a steady-state model of the isothermal tube drawing process with the help of shooting method. The obtained analytical solution is almost in agreement with the numerical solution that justifies our approach.