Stability Characterization of Three Porous Layers Model in the Presence of Transverse Magnetic Field

The current study concerns, the effect of a horizontal magnetic field on the stability of three horizontal finite layers of immiscible fluids in porous media. The problem examines few representatives of porous media, in which the porous media are assumed to be uniform, homogeneous and isotropic. The dispersion relations are derived using suitable boundary and surface conditions in the form of two simultaneous Mathieu equations of damping terms having complex coefficients. The stability conditions of the perturbed system of linear evolution equations are investigated both analytically and numerically and stability diagrams are obtained. The stability diagrams are discussed in detail in terms of various parameters governing the flow on the stability behavior of the system such as the streaming velocity, permeability of the porous medium and the magnetic properties. In the special case of uniform velocity, the fluid motion has been displayed in terms of streamlines concept, in which the streamlines contours are plotted. In the uniform velocity motion, a fourth order polynomial equation with complex coefficients is obtained. According to the complexity of the mathematical treatments, when the periodicity of the velocity is taken into account, the method of multiple scales is applied to obtain stability solution for the considered system. It is found that a stability effect is found for increasing, the magnetic permeability ratio, the magnetic field, and the permeability parameter while the opposite influence is observed for increasing the upper layer velocity.


Introduction
The flow instability of a plane interface between two superposed fluids of different densities through porous media is of considerable interest for petroleum engineers and in geophysical fluid dynamicists.On the other hand, the phenomenon of resonance is a fundamental one in mathematics and physics.This phenomenon is obtained from nonlinear interactions among a few (two or three) wavetrains that was initially pointed out by (Phillips, 1960) and (Longuet-Higgins, 1962), followed by experimental verifications (Mcgoldrick et al., 1966) and (Longuent-Higgins & Smith, 1966).They indicated that, under certain conditions, conspicuous energy transfer takes place from primary waves to a tertiary wave, produced through the third-order interaction.A series of studies for hydrodynamics stability have been initiated by many authors, for example, the weakly nonlinear stability is employed to analyze the interfacial phenomenon of two magnetic fluids in porous media.has investigated in paper (El-Dib & Ghaly, 2003).The method of multiple scale expansion is employed in order to obtain a dispersion relation for the first-order problem and nonlinear Ginzburg-Landau equation, for the higherorder problem, describing the behavior of the system in a nonlinear approach.In paper (Zakaria et al., 2009), the instability properties of streaming superposed conducting fluids through porous media under the influence of uniform magnetic field have been investigated, where the system is composed of a middle fluid sheet of finite thickness embedded between two semi-infinite fluids.
A good account of hydrodynamic stability problems has also been given in papers (Drazin & Reid, 1981), (Joseph, 1976), (Sisoev et al., 2009), (Sadiq et al., 2010) and (Rosensweig, 1985).The authors in the paper (Funada & Joseph, 2001) have discussed the instability of viscous potential flow in a horizontal rectangular channel.Bhatia (1974) has studied the influence of viscosity on the stability of the plane interface separating two incompressible superposed fluids of uniform densities, when the whole system is immersed in a uniform horizontal magnetic field.He has developed the stability analysis for two fluids of equal kinematic viscosities and different uniform densities.Li et al., (2007) have examined the electrohydrodynamic stability of the interface between two superposed viscous fluids in a channel subjected to a normal electric field.The long wave linear stability analysis is performed within the generic OrrSommerfeld framework for both perfect and leaky dielectrics.The approach proposed in paper (Tseluiko & Blyth, 2008) is limited to study the gravitydriven flow of a liquid film down an inclined wall with periodic indentations in the presence of a normal electric field.Kumar and Singh (2006) have investigated the stability of a plane interface separating two viscoelastic (Rivlin-Ericksen) superposed fluids in the presence of suspended particles.They concluded that the system is stable for stable configuration and unstable for unstable configuration in the presence of suspended particles.Khan et al., (2007) have demonstrated the analytical solutions for the magnetohydrodynamic flow of an Oldroyd-B fluid through a porous medium.They obtained the expressions for the velocity field and the tangential stress by means of the Fourier sine transform.In his study of hydromagnetic parametric resonance instability of two superposed conducting fluids in porous medium El-Sayed (2007) illustrated Rayleigh-Taylor instability of a heavy fluid supported by a lighter one through porous medium, in the presence of a uniform, horizontal and oscillating magnetic field .The fluids are taken as viscous (obeying Darcy's law), uniform, incompressible, and infinitely conducting, where he amplitude of the oscillating part of the field is taken to be small compared with its steady part.
In view of these, the goal here is to develop a mathematical model for a streaming fluid sheet embedded between two bounded fluid layers in the presence of porous media.The fluids are subjected to a horizontal magnetic field.The results illustrated in this work may be of interest for the fluid dynamic aspects of an aquifer or an oil well and in petroleum reservoirs, spray coating process, plastics manufacture, metal powder production and in microchips fabrication.The rest of this paper is organized as follows.The next section lays down geometry and mathematical formulation of the problem and a sketch of the system under consideration.Also in this section the equations of motion and boundary conditions are derived.The third section is prepared to linear stability analysis.In the fourth, the dispersion relations, for a periodic velocity, have been derived to control the surface wave propagation and the streamlines contours are plotted an discussed.Also, in this section, the perturbation scheme using the multiple scales analysis, and numerical estimation for stability configuration have been discussed.The results and some important conclusions are outlined in last section of this work.

Geometry and Mathematical Formulation
The physical configuration is shown in Fig. 1.We consider a parallel flow consists of a model of a liquid sheet sandwiched between two bounded fluid layers, fully saturated, uniform, homogeneous and isotropic porous media with constant permeability.The co-ordinates x and y spans the horizontal and vertical directions.In order to relax the mathematical calculations, we consider that the fluids are assumed to be immiscible incompressible and have constant properties.The system is considered to be influenced by the gravity force which acts in the negative y-direction.The upper fluid is referred to as fluid 1, while the lower one is denoted by fluid 3, whereas the middle layer is distinguish by fluid 3. Fluid r (r = 1, 2, 3) is assumed to have density ρ r , a pressure function of the fluids p r and magnetic permeability μ r .
There are two parallel interfaces between the three fluids are assumed to be well defined and initially flat and forms the plane y = (−1) l+1 L, l = 1, 2 and the instantaneous perturbed interface height is y=h l (x, t) is along the y direction.Suppose that the layers are moving with velocity V 0r = V 0r cos t î, where V 0r and are constants.The unit vectors î and ĵ are in x− and y− directions.In addition, the system is initially assumed to be stressed by a uniform magnetic field of intensity H 0r along the direction tangential to the flat interface y = L.
Fluid motion is governed by a set of nonlinear partial differential equations expressing conservation of mass, momentum and the field energy.On the other hand, flow in a porous medium is described by Darcy's law that relates the movement of fluid to the pressure gradients acting on a parcel of fluid and we consider media that are initially uniform so that motion is of homogeneous fluids in a homogeneous medium.Hence, the basic equations governing the motion of an incompressible liquid through porous medium are comes from the combination of the momentum equation and Darcy's law (see (Sisoev et al., 2009), (Sadiq et al., 2010), (Rosensweig, 1985), (Funada & Joseph, 2001) and (Bhatia, 1974) ): For incompressible fluid flow, the continuity equation is based on the principle of conservation of mass, which reads Here, in these equations the symbol ∇ ≡ (∂x, ∂y) is the horizontal gradient operator, the vector g = (0, −g) is the acceleration due to gravity.The resistance is modelled by the ratio μ 0r = μ0r /Q r , where the permeability Q r describes the ability of the fluid to flow through the porous medium and μ0r is the fluid viscosity measures the resistance of fluid to shearing that is necessary for flow.In addition, it should be noticed that Q r and μ0r are assumed to be constant.Introducing the velocity potential ϕ r (x, y, t) of the perturbed motion such that the total fluid velocity is given by where, u r = −∂ x ϕ r and v r = −∂ y ϕ r are the velocity components due to disturbances, thus ϕ r will satisfy Laplace's equation We will assume that there are no free currents at the surface of separation in the equilibrium state, and hence in a magnetoquasi-static system with negligible displacement current, Maxwell's equations will be reduced to Gauss and Ampére laws, which can be expressed as Here, H r refers to the magnetic field intensity vector and μ r H r is defined as the magnetic induction vector.Since, in magnetic fluids the magnetic energy greatly exceeds electric energy storage and where the propagation times of electromagnetic waves are complete in relatively short times compared to those of interest to us.Thus, the construction of a potential function χ r (x, y, t), can be representable as the gradient of the scalar potential such that So, we can introduce the magnetic potential χ r that satisfies zero curl for a constant permittivity and therefore this potential satisfies the Laplace's equation: The above equations are general governing equations of motion for flow of a homogeneous fluid in an isotropic medium with constant viscosity and permeability.These equations can be solved to investigate a magnetic fluid flow through a porous medium.To complete the formulation of the problem, we must supplement of the momentum and magnetic field equations with the corresponding boundary conditions.These constraints are information about the solutions at the upper and lower boundaries and at the interfaces between the fluids and thus the system adopted here meets the following boundary conditions (see (Joseph, 1976), (Sisoev et al., 2009), (Sadiq et al., 2010) and (Rosensweig, 1985) ).
On the interface y = (−1) l+1 L + h l (x, t), l = 1, 2 it is natural to impose the kinematic boundary conditions, the kinematic condition expresses the fact that the interface always comprises the same fluid particles, and therefore the function h l (x, t) whose graph defines the interface satisfies simultaneously In addition kinematic relation follows from the assumption that the normal component of the velocity vector in each of the phases of the system is continuous at the dividing surface: Here, n l is the exterior pointing normal unit vector to the interfaces which has the form , where F l (x, y, t) is the surface geometry defined by the locus of points satisfying the relation F l = y− (−1) l+1 L + h l (x, t) .The boundary conditions on the upper and lower plates, in which the plates are assumed to be rigid and kept constant, this implies that: Since we are dealing with the case of pure magnetization effects, the continuity of the normal and the tangential components of the magnetic displacement at the interface obeyed Maxwell's conditions, and thus we have The surface force that accounted by the stress tensor is the Maxwell tensor that describes the stress field induced in the material due to electrostatic forces whose expression is where, I is the identity tensor.Furthermore the dynamical boundary condition, where the normal stresses are balanced by the amount of the surface tension is where, it is assumed that the fluid interfaces have surface tension coefficient γ l(l+1) .In deed, sharp interfaces between the fluids may not exist.Rather, there is an ill-defined transition region in which the two fluids intermix.The width of this transition zone is usually small compared with the other characteristic length of the motion; hence, for the purpose of the mathematical analysis, we will assume that the fluids are separated by sharp interfaces.In addition, the above boundary conditions are prescribed at the interface y = (−1) l+1 L + h l (x, t).It is necessary to express all the physical quantities involved in terms of Taylor series expansion about y = (−1) l+1 L.

Stability Analysis and the Solution Method
In illustrating the problem in the light of linear perturbation, the second order as well as the higher-orders terms containing the elevation parameter h l (x, t) are neglected.So, to treat the stabilization of the problem under consideration, the amplitude of waves formed on the fluid sheet is assumed to be small and a finite disturbances are introduced into the equation of motion and continuity equation as well as the above boundary conditions.Hence, for a small departure from the equilibrium state, every physical perturbed quantity may be expressed as functions of both the horizontal and vertical co-ordinates as well as time: Here, this analysis based on a normal modes technique, where k is the wave number, which is assayed to be real and positive and c.c. represents complex conjugate of the preceding terms and S stands for ϕ and χ.These, expansions are introduced into the governing equations and the relevant boundary conditions.The linearized terms in these perturbed quantities are only maintained in view of the linear stability theory (Chandrasekhar, 1961) and (Moatimid, 2003).To perform a linear stability analysis of the present problem, the interfaces between the three fluids will be assumed to be perturbed about their equilibrium locations to cause displacements of the material particles of the fluid system.Consider the effect of small wave disturbances to the interfaces y = (−1) l+1 L, propagating in the positive x−direction.Assuming that the surface deflections are given by y where and ξ 1 and ξ 2 are arbitrary time-dependent functions which determine the behavior of the amplitude of the disturbances on the interfaces.It should be noted here the linear term, which comes from the nonlinear term that appear in the left hand side in Eq. ( 1) is neglected in this analysis.This is because of the averaging process through which Darcy's Eq.
(1) has been derived.For low Reynolds number flows O(1) and thus this term can be ignored.So, in the hydrodynamic description by inserting Eq. (3) into Eq.(1), we obtain the pressure in terms of the velocity potential such that Substituting Eqs. ( 13) and (18) into Eq.( 14), the balance at the dividing surfaces can have the relation Actually, the bulk solutions are written in accordance with the interface deflection given by ( 17) and in view of a standard Fourier decomposition, these solutions can be similarly expressed as To investigate and discuss the boundary-value problem cited above, we depend on the hydrodynamic stability analysis that given in the book (Chandrasekhar, 1961).It constitutes a homogeneous system of equations and boundary conditions for explaining the factors governing the surface wave's propagation.In view of the above boundary conditions, the solution of Laplace's equation yields the distribution of the velocity potential ϕ r and the magnetic potential χ r in the three layers.
On substituting Eq. ( 20) into Laplace's Eq. ( 4); the resulting solutions in the three fluid phases with the aid of the above kinematic boundary conditions can be obtained as where, Inserting (21) into Laplace's Eq. ( 7), the resulting solutions in view of the previous Maxwell's conditions will give χ 1 (x, y, t) = iH 01 sinh k(y + L 1 ) Γ 11 (μ 1 − 1)ξ 1 (t) + Γ 12 (μ 2 − 1)ξ 2 (t) exp(ikx) + c.c, (25) where, the coefficients Γ s that appear in these relations are clear from the context.These distribution are derived in the light of the linearized form of the appropriate boundary conditions and it is equivalent to those obtained before by Rosensweig (1985).On the other hand, the units in the above solutions are removed by using the dimensionless quantities as, the stream velocity and the velocity potential function are made dimensionless using √ Lg and L √ Lg, while the applied field and the magnetic potential are made dimensionless by Lgρ 2 /μ 2 and L Lgρ 2 /μ 2 , respectively.In addition the viscosity ρ 2 L 2 g, permeability of the porous medium L 2 Q and the time by L/g .And by using the symbols μ1 = μ 1 /μ 2 , μ2 = μ 3 /μ 2 the density ρ1 = ρ 1 /ρ 2 , ρ2 = ρ 3 /ρ 2 , the Weber number W l = γ l /L 2 gρ 2 , (l = 1, 2).

Derivation of the Dispersion Relations
Our goal in this section is to study effect of general surface deformations on the onset of a periodic velocity applied to the fluid sheet.Equations that determine the surface deflections are called the characteristic equations.The whole system will be reduced to a two coupled partial differential equations with periodic coefficients, in the elevation parameter ξ l as a dependent variable, by substituting Eqs.(22-27) into the normal stress tensor (19).The dependence on the potential velocity ϕ r , the magnetic potential function χ r and the fluid pressure function p r is replace by their equivalence in terms of the amplitude ξ l , finally after a straightforward calculations, one obtains the coupled equations where the symbol D t = d/dt refers to the derivative with respect to time t.The coefficients that appear in these equations are real and depend on the physical parameters of the problem.The mathematical formulas of these coefficients are lengthy and not included here.However, they are available upon request from the author.Eqs. ( 28) and ( 29) are two coupled Mathieu equations having damping terms and complex coefficients.By making use of these equations, the stability behavior of the fluid sheet is controlled.For a uniform stream, the periodicity of the stream will be absent.Therefore, wave propagation is excited by using the electro-capillarity technique.Hence, in the limiting case of tending to zero in the above system, the damped Mathieu equations then become a linear differential equations with constant coefficients.It can be satisfied by a growth rate solution, which may be written as where ξl is the constant of integration.Substituting this equation into the above system of Mathieu equations the dispersion equation of the perturbed motion is then where the coefficients α's are clear from the context.It should be noted that (31) represents a complex linear dispersion relation that is satisfied by values of ˆ and k.In the limiting case of nonporous media, the above dispersion relation will reduce to those obtained by Rosensweig (1985).It is clear that the surface waves propagating along the interfaces separating between the inviscid fluids will only be stable if all the roots of (31) are real.Otherwise, there are at least two roots (complex conjugate) and thereby the interfacial inviscid waves are unstable.

Streamlines Distribution
An important concept in the study of fluid dynamics concerns the idea of streamlines.Streamlines in the physical domain are a family of curves (sometimes called curvilinear) that are instantaneously tangent to the velocity vector of the flow, resulting in a rectangular computational region.These show the direction a massless fluid element will travel in at any point in time.In other meaning a streamline is a path traced out by a massless particle as it moves with the flow.In order to study the concept of these lines, let us define a stream function, ψ r of the time and space coordinates, where the following relationships between the velocity components (u r and v r ) and stream function have been used u r = ∂ y ψ r , v r = −∂ x ψ r , which automatically satisfies the continuity equation.Eq. ( 2) together with these equations give the relation between the stream function and the velocity potential as ψ r = − ∂ y ϕ r dx.In the following we plot the stream lines through the stream function to show its concept on the stability of the movement of the waves, in which the streamlines are effective tools to visualize a qualitative impression of the flow behavior during the motion.
In Figs. 2 and 3, the streamlines profile is plotted by fixing the value of all the physical parameters except for one parameter has varying values for comparison, snapshots of instantaneous streamlines of the stream function, are shown in these graphs.
Figure 2. Streamlines contours for a system having the parameters  The influence of the magnetic permeability ratio μ2 is presented throughout the parts of Fig. 2 for a system having the parameters given in the caption of Fig. 2. The inspection of Fig. 2(a), where μ2 =0.2 reveals that the flow consists of three cells (contours) consisting of two clockwise (left and right, negative values of streamlines) and the middle contours is anti clockwise ( positive values of streamlines) circulations.In part (b) of this graph, the value of the magnetic permeability ratio μ2 is increased to be 0.5.It is worthwhile to notice that the three streamlines contours are shifted to the left side, until reduce to two in the part 3, when μ2 is increased to the value 0.7.A conclusion that may be made from the comparison among the parts (a-c) of Fig. 2 is that the magnetic permeability ratio μ2 leads to crowd in the concentration of the streamlines in the movement of the fluids.
Fig. 3 illustrates streamlines under the same values considered in the above system of Fig. 2, but at H 01 =20, while the velocity of the upper layer V 01 has three different values for the sake of comparison.In Fig. 3(a), where V 01 =3, it is shown that the flow consists of three contours, the middle cell is clockwise has negative values of streamlines and left and right contours are anti clockwise have positive values of streamlines.In the part (b) of this graph the velocity V 01 is increased to the value 6, it is observed from this figure that the increasing of the velocity leads to crowd of the streamlines cells, which are shifted above and replaced by another three ones in different orientations.The value V 01 =9 is added to the part (c) of Fig. 3 , we noticed that the contours of the streamlines is contracted at the center, until it is divided into three contours to the top of this graph.A general conclusion of Figs. 2 and 3 is that the magnetic permeability ratio μ2 has a stabilizing influence on the movement on the waves, while the opposite effect is found for increasing the upper layer velocity.

Periodicity and Numerical Results
Due to the periodicity of the velocity, the stability picture has changed dramatically and hence we return to the general form of dispersion relations ( 28) and ( 29).The nature of the solution of these equations will govern the fluctuations of the amplitude of the interface deflection and will, therefore, determine the parametric excitation of the interfacial waves.
According to the complexity of the mathematical problem we use a perturbation technique, one of this technique is the method of multiple time scales (Nayfeh, 1979) which has been successfully used to treat similar these equations, since the solutions and the properties of Eqs. ( 28) and ( 29) are unknown.Applying the method of multiple scales, where the independent variable t can be extended to introduce alternative independent variables: t n = t ,n = 0, 1, 2, where the parameter represents a small dimensionless parameter characterizing the steepness ratio of the wave.Thus, we define t 0 as the variables appropriate to fast variations and t 1 , t 2 as the slow variables.The differential operators can now be expressed as the derivative expansions: where t 0 is the time of the lowest order.For the small dimensionless parameter , we can characterize the amplitude of the periodic force which is defined as V 0r = Ṽ0r .The analysis then follows the usual perturbation procedure and suppression of the secular terms except that is now more convenient to write the solution in a complex form.Now, let the dependent variables ξ l be expanded in the form ξ l (t, ) = ξ 0l (t 0 , t 1 ) + ξ 1l (t 0 , t 1 ) + ..., l = 1, 2.
Inserting ( 32) and ( 33) into ( 28) and ( 29) and equating coefficients of like powers of (because each of the ξ l are independent of ) yields simpler inhomogeneous equations, which can be solved successively with knowledge of the solutions of the previous orders.Uniform solutions are required to eliminate the secular terms.This elimination produces the solvability conditions corresponding to the terms containing the factor exp(i ˆ t 0 ), in which the solvability condition is divided into two cases.The first is valid in the non-resonant case in which the frequency is away from the frequency ˆ .Otherwise the resonance arises when the frequency approaches the frequency ˆ .Hence the solvability condition in the non-resonant case is ( This condition show that the motion is stable if In the resonance case when the frequency approaches the frequency ˆ , we introducing a detuning parameter λ (1) defined by and hence the solvability conditions are Figure 4. Illustrated in the plane (λ (1) − k) according to equations ( 42), for a system having Figure 5.The graph is constructed for λ (1) versus k, for the same system given in Fig. 4, with with the permeability parameter Q 2 = 2.5, 3.5 and 4.5 of the partitions (a), (b) and (c), respectively.
In order to examine the effect of the permeability parameter in the middle layer Q 2 on the stability criteria, numerical calculations are made in the parts of Fig. 5 The graph shown in the plane (λ (1) − k) are achieved for three values of the permeability parameter Q 2 =2.5, 3.5, and 4.5, corresponding to the partitions (a), (b) and (c) respectively, where the other quantities are held fixed.The inspection of the stability diagram of the parts of Fig. 5 reveals that the increase of the permeability parameter leads to increase in the width of the stability regions, while the unstable areas are decrease.The conclusion that may be drawn here is that the permeability parameter has a stabilizing influence on the stability behavior of the waves.In the parts (a), (b) and (c) of Fig. 6, we repeat the same diagrams as illustrated in Fig. 4, with a change in the value of the effect of the magnetic field H 01 = 10, 20, 30, while the other parameters are fixed.Applying the above stability constraints to separate the stable and the unstable regions, we notice that the magnetic field plays a stabilizing role in the stability criteria.This influence may be physically interpreted as suggesting that some of the kinetic energy of the waves has been transferred to the magnetic field.

Conclusions
Theoretical and numerical analysis of linear stability of a fluid sheet of finite thickness embedded between two bounded layers of fluids through porous media are carried out.The system is under the influence of a horizontal magnetic field with a periodic stream.Two linear dispersion equations of Mathieu type have been derived and involved parametric coefficients as well as parametric imaginary damping term.These equations are used to control the stability of the fluid sheet motion.The stability analyses have been investigated by using the multiple timescales method.Consequently, a mathematical simplification is desired to relieve this complication for the Mathieu equation.The transition curves separating the stable region from unstable regions are identified.The analysis recovers the key numerical findings and provides qualitative understanding.Numerical calculation of the stability of the system are made where the physical parameters are put in the dimensionless form.Stability diagrams are plotted and discussed for different sets of physical parameters.The stability examination yields the following results.In the case of uniform velocity, the streamlines contours are plotted an discussed, where the results show that the magnetic permeability ratio has a stabilizing influence on the movement on the waves, while the opposite influence is found for increasing the velocity of the upper layer.When the periodicity of the velocity is taken into accunt, the method of multiple scales is applied to obtain stability solution for the considered system.It is found that the velocity of the lower layer has a destabilizing effect whereas the permeability parameter and the magnetic

Figure 1 .
Figure 1.The physical model and coordinate system.

Figure 3 .
Figure 3.The same system as that considered in Fig. 2, while the velocity V 01 = 3, 6 and 9 of the parts (a), (b), and (c), respectively.