A Numerical Study of the Flow with Heat Transfer of a Pseudoplastic Fluid between Parallel Plates

One dimensional flow with heat transfer of a pseudoplastic fluid between two infinite horizontal parallel plates is investigated. The thermophysical properties of the fluid are assumed to be constant and numerical solution using the finite element method, along with the corresponding exact solution for the fluid velocity and the fluid temperature is obtained. The effect of variation of the governing parameters is studied using figures and tables. It is found that the numerical solution agrees well with the corresponding exact solution and that the fluid velocity, together with the fluid temperature, increases with increasing values of the governing parameters.


Introduction
Over the past few decades there has been a growing recognition of the fact that many fluids of industrial significance do not obey the Newtonian postulate of linear relationship between the shear stress and shear rate.Therefore, these fluids are known as non-Newtonian fluids.Common examples of such fluids are slurries, pastes, gels, molten plastics and lubricants containing polymer additives.Various food stuffs such as honey and tomato sauce, the biological fluids like blood and synovial fluid naturally found in cavities of synovial joints also belong to the general class of the non-Newtonian fluids [1,2].The simplest model of non-Newtonian fluids is the power law model and a special case of the power law equation, with , is the pseudoplastic model.Commonly encountered non-Newtonian fluids like polymer solutions, paper pulps, detergents, oils and greases may be classified by the pseudoplastic fluid model [3].The pseudoplastic fluids represent shear thinning fluids, and as far as we know, this class of non-Newtonian fluids received less attention in the literature.

< 1 n
Heat transfer processing of non-Newtonian fluids is one of the key components of the flow phenomenon in various industrial sectors including chemicals, petro chemicals, food industry, polymers and pharmaceuticals [3][4][5][6].Therefore, in the present paper we consider one-dimensional steady laminar flow, with heat transfer, of a pseudo plastic fluid between two fixed infinite horizontal parallel plates.The fluid motion is generated by the presence of a constant external pressure gradient along x-axis and the fluid velocity components along yand z-directions are zero.Thus, the momentum equations break down to a second order ordinary differential equation and we need two boundary conditions for finding a solution to this equation.These boundary conditions are obtained from the no slip condition and the condition of symmetry at the center line of the flow channel.We integrate the resulting momentum equation once in combination with the condition of symmetry to obtain a quadratic equation in the unknown shear rate d d u y .We apply the Finite Element Method (FEM) to solve this quadratic equation numerically and compare the numerical solution with the corresponding exact solution of the equation.For details on the implementation of the FEM the interested reader may consult from the references [7][8][9][10].To find exact solution of the problem, we follow the approach used by Deiber and Santa Cruz [11] and find two values of the shear rate by solving the resulting quadratic equation.However, we retain only one of these two values as the other does not satisfy the condition of symmetry.Then a second integration of / du dy d d u y in combination with the no-slip condition provides an exact solution for the fluid velocity.Once the fluid velocity solution is known, we may use it to find solution for the fluid temperature from the resulting energy equation.
The following discussion begins with the basic equations, governing the flow of an incompressible fluid, in section 2. The problem is, respectively, formulated and solved in sections 3 and 4. Results obtained from the numerical and exact solutions for various values of the governing parameters are presented and discussed in section 5. Finally some conclusions are mentioned in section 6.

Basic Equations
The basic equations governing the flow of an incompressible fluid are given by the laws of the conservation of mass, the conservation of momentum and the conservation of energy, see [1,2] where u is the velocity field, f is the body force, T is the stress tensor, is the fluid temperature,   is the constant fluid density, is the thermal conductivity,  p C is the specific heat, L is the gradient of u and D Dt is the material derivative.The stress tensor T for pseudo plastic fluid is defined by the constitutive equations, see [2].
where p is the reaction stress due to the constraint of incompressibility.The extra stress tensor S is defined by the expression, see [11]    where o  is the zero shear viscosity, 1  is the relaxation time and 1  is the material constant.The first Rivlin-Ericksen tensor 1 is defined as 1 A = T  A L L and the contravariant convected derivative denoted by superimposed over S is defined as

Problem Formulation
We where   where P  denotes the constant pressure gradient d d p x .Similarly, using (7), we see that  9) and ( 10), we introduce the following dimensionless variables.
where 0 is the bulk mean temperature and U is a typical fluid velocity.After omitting asterisks, the Equations ( 9) and ( 10) can be written in their respective dimensionless form as follows where is the Brinkman number, is the Prandtle and is the Eckert number.Finally, non-dimensionalising the remaining boundary condition on u and the two boundary conditions on  we obtain , and .
Pr Ec

Solution of the Problem
In this section we apply the FEM to solve Equation ( 12) subject to and the resulting energy Equation (13) subject to the boundary conditions and for finding fluid velocity and temperature fields, respectively.We also find exact solution for the sake of comparison.
 

Fluid Velocity
Clearly, for the choice = 0  , the Equation ( 12) reduces to the corresponding equation for viscous flow.Therefore, in order to maintain non-newtonian characteristics of the velocity Equation ( 12), we assume that 0   .
The Galerkin's finite element method provides an elegant approach for constructing approximations of solutions of boundary value problems.It provides a general and precise procedure for constructing basis functions.The main idea is that the basis functions can be defined piecewise over sub-regions of the domain called finite elements and that over any sub-domain, the basis functions can be chosen to be very simple functions such a polynomials of low degree.Therefore, the Galerkin's formulation in the finite element fashion requires, see [7][8][9][10], that we choose a suitable trial or basis approximation V to the true solution V that is applied locally over a typical finite element in the complete y domain.The simplest such solution is linear trial function, .From finite element analysis, rather than formulating the problem in terms of arbitrary constants 1 & 2 , we prefer to recast the bove linear trial function in terms of values of the dependent functions at nodes i & j [7,8], N N 2 is a vector of the interpolation or shape functions with are called linear Lagrange interpolation polynomials, now are the nodal values of the dependent variable .
V  Now applying quasi-liberalization, see [12], to first of the governing system of differential Equation ( 12) and (13), we obtain We follow the standard Galerkin's approach and choose the weighting function .Now integrating over the entire region, we get where ' d /dy  .For our particular problem, after substituting the trial functions in Equation ( 16), we have the following discretised form ) where e stands for element and n represents the total number of elements in the discretised region.In matrix notation, the system of Equations ( 17) can be written as where k is a stiffness matrix and f is a force vector, =    i y and = . Applying the assembly procedure given in the references [7-10] and using Equation (18) for n elements, we get The solution of the system (19) provides the numerical solution of Equation (12) for fluid velocity profile.The final arrangement of the stiffness matrix (19) is noteworthy, nonzero entries come into sight clustered near the main diagonal of the matrix.Outside this band of nonzero terms, all entries are zero.Matrices of this type are said to be banded.Also, it is seen in ( 19) that the stiffness matrix is symmetric, though we shall not always have symmetry, nevertheless in most physical problems based on conservation laws this symmetry will arise quite naturally.This symmetry will always be possible to attain in self-adjoint boundary value problem.These two properties of stiffness matrices play a vital role in the stratagem of programming finite element calculations.
For finding the exact solution of the Equation ( 12) we solve this equation for d d u y and obtain  22) is exact solution for the fluid velocity.Clearly this expression for the fluid velocity can have a physical meaning if it is defined only in the set of real numbers.Therefore, the solution given by Equation ( 22) is physically valid provided the condition 2 4 1 P   is satisfied.This further means that, if the pressure gradient P is fixed, the parameter  should be strictly set in the range .Similarly, if the parameter fixed at a positive value the P should take a value that ensures validity of the condition 0 < < 1 2 P  .

Fluid Temperature
Now, using the fluid velocity solution obtained from the system (19) and applying the same Galerkin's FEM formulation to the second of the system of governing differential Equations ( 12) and ( 13), we get which in matrix notation, can be written as where 1 f is a force vector.By applying the assembly procedure given in [7][8][9][10], and using (23) for n elements, we get The solution of the system (25) gives us the temperature distributions of the flow by imposing the boundary conditions.
Next, for evaluation of the temperature distribution, we substitute the fluid velocity solution (22) in the Equation (13) to obtain provided the condition is satisfied.Now integrating Equation (26) twice and applying the boundary conditions   Where

Results and Discussion
parameters  and P is presented in Figure 1(b).A satisfactory agreement of numerical solution with the corresponding exact solution may be observed from these figures.In Table 1 we show variation of the As mentioned earlier, for a fixed P, the parameter  is strictly set in the range

E
. From this table we observe that size of the error u is sufficiently small to confirm agrement of the numerical solution with the corresponding exact solution.Furthermore, we observe from Figures 1(a   From these figures we observe that the agreement of the numerical solution with the corresponding exact solution is very good.This agreement is further highlighted by sufficiently small size of the absolute error shown in Table 2 for    0.1, 0.6, 0.9 E  is small enough to illustrate accuracy of the numerical solution.Moreover, we observe that the fluid temperature increases with increasing values of P for fixed  and  are fixed. Then in Figure 4(a) we present variation of the FEM numerical solution for = 0.25, = 0.5 P  and  1.0,3.0,5.0  

Conclusions
One dimensional flow of a constant property pseudo plastic fluid with heat transfer between two infinite horizontal parallel plates is considered.Numerical solutions for the dimensionless fluid velocity and the fluid temperature are obtained by applying the finite element method.The corresponding exact solutions are also derived for the sake of comparison.It is found that the numerical technique based on the finite element method produces highly accurate solution both for the fluid velocity and the fluid temperature that agrees nicely with the corresponding exact solution.Moreover, the fluid velocity along with the fluid temperature increases with increasing values of the pressure gradient P, the non-Newtonian parameter  and/or the Brinkman number  .
.The corresponding exact solution for the same values of the parameters ,   and P is presented in Figure 4(b) for the sake of comparison.
From these figures we again observe that the agreement of the numerical solution for the fluid temperature with the corresponding exact solution is very good.We present the absolute error E  for = 0.25,  and in Table 4 and observe that the size of is satisfactorily small to confirm agreement of the numerical solution with the corresponding exact solution.Moreover, as we see from Further, we report without presenting that the numerical solution for the fluid temperature and the corresponding absolute error shows similar behavior with varying  E   and fixed P and  .Therefore, we conclude that the fluid temperature increases with increasing values of the non-Newtonian parameter  , the Brinkman number  and the pressure gradient P, when any two of these parameters are held fixed.

 1 ,
which are provided by the fluid temperature on the plates.For non-dimensionalization of Equations (

2
retain only one of the roots (20) because the other does not satisfy the condition of symmetry and write at the following exact solution for the fluid temperature.
value of the parameter  should lie in the range 0 < 1   .Similarly if the non-Newtonian parameter is taken to be =
for the fluid velocity with same values of  and P.


velocity magnitude again appears to show an increasing trend with increasing values of the pressure gradient P. Thus we conclude that the fluid velocity magnitude increases with increasing values of the non-Newtonian parameter and the pressure gradient P, when one of these is fixed.

Figure 3 (
a) illustrates variation of the numerical solution for the fluid temperature obtained by the finite element method with = exact solution is presented in Figure3(b)for same values of the parameters  0.1, 0.6, P   and P. From these figures we observe that the numerical solution for the fluid temperature exhibits good agreement with the corresponding exact solution.This agreement of the numerical solution with exact solution is further illustrated in

Figure 4 ,
again the fluid temperature shows an increasing trend with increasing A. M. Siddiqui, A. Zeb, Q. K. Ghori and A. M. Benharbit, "Homotopy Perturbation Method for Heat Transfer Flow of a Third Grade Fluid between Parallel Plates," Chaos Solitons and Fractals, Vol.36, No. 1, 2008, pp.182-192., when P and  are held fixed.

Table 3 .
This table presents variation of the absolute error