# Nonlinear Magnetoconvection in a Sparsely Packed Porous Medium

Authors: A. Benerji Babu; Ragoju Ravi; S. G. Tagare
Journal: International Journal of Geophysics
 (2011)

---

## Abstract

Linear and weakly nonlinear properties of magnetoconvection in a sparsely
packed porous medium are investigated. We have obtained the values of Takens-Bogdanov bifurcation points and codimension two bifurcation points by plotting
graphs of neutral curves corresponding to stationary and oscillatory convection
for different values of physical parameters relevant to magnetoconvection in a
sparsely packed porous medium near a supercritical pitchfork bifurcation. We
have derived a nonlinear two-dimensional Ginzburg-Landau equation with real
coefficients by using Newell-Whitehead (1969) method. The effect of the parameter
values on the stability mode is investigated and shown the occurrence of
secondary instabilities namely, Eckhaus and Zigzag instabilities. We have studied
Nessult number contribution at the onset of stationary convection. We have also
derived two nonlinear one-dimensional coupled Ginzburg-Landau-type equations
with complex coefficients near the onset of oscillatory convection at a supercritical
Hopf bifurcation and discussed the stability regions of standing and travelling
waves.


---

## Body

## 1. Introduction

Magnetoconvection in a porous medium uniformly heated from below is of considerable interest in geophysical fluid dynamics, as this phenomena may occur within the mushy layer of Earth's outer core. Earth's outer core consists of molten iron and lighter alloying element, sulphur in its molten form. This lighter alloying element present in the liquid phase is released as the new iron freezes due to supercooling onto the solid inner core. Hence we get mushy layer near the inner core boundary where the problem becomes convective instability in a porous medium [1]. The effect of geomagnetic field on the magnetoconvection instability is of interest in geophysics, particular in the study of Earth's interior where the molten liquid Iron is electrically conducting, which can become convectively unstable as a result of differential diffusion.Magnetoconvection in an electrically conducting fluid in a nonporous medium has been studied extensively [2–8]. However, magnetoconvection in a porous medium has not received any attention inspite of its application in geophysical fluid dynamics problems. Palm et al., [9] investigated Rayleigh-Benard convection problem in a porous medium. Brand and Steinberg [10, 11] investigated convecting instabilities in binary liquid in a porous medium; However, Plam et al. [9] and Brand et al. have made use of Darcy's law (-ν∇2V¯ is replaced by KV¯ where K is the permeability of a porous medium. for nonporous medium K is infinity). They have also not considered usual convective nonlinearity. It is well known that Darcy's law breaks down in situations where in other effects like viscous shear and inertia come into play. In fact Darcy's law is applicable to densely packed porous medium. An alternative to Darcy's equation is Brinkman equation and is of the form(1)∇′ρ′-ρ′g¯=-μKV¯′+μe∇′2V¯′,
where μ is the fluid viscosity and μe is the effective fluid viscosity. Brinkman model is valid for a sparsely packed porous medium wherein there is more window fluid to flow so that the distortion of velocity give rise to the usual shear force. Lapwood [13] was the first to suggest the inclusion of convective term (V¯′·∇′)V¯′ in the momentum equation and study the Rayleigh-Benard convection in a sparsely packed porous medium. Recently, Tagare and Benerji [14] have investigated the problem of nonlinear convection in a sparsely packed porous medium due to thermal and compositional buoyancy.In this paper we investigate the problem of magnetoconvection in a sparsely packed porous medium. The multiplicity of control parameters makes this system an interesting one for the study of hydrodynamic stability, bifurcation and turbulence [15]. Rudraiah [16] and  Rudraiah and Vortmeyer [17] have studied both linear and steady nonlinear magnetoconvection in a sparsely packed porous medium using Brinkman model but they have taken effective viscosity μe same as fluid viscosity μ. However, experiments show that the ratio of effective viscosity μe to fluid viscosity μ takes the value ranging from 0.5 to 10.9 [18]. In Section 2, we write basic dimensionless equations in Boussinesq approximation for magnetoconvection in a sparsely packed medium by using a momentum equation with effective viscosity different from fluid viscosity. In Section 3, we study linear stability analysis. In Section 4.1, by using multiple-scale analysis of Newell and Whitehead [19], we derive two-dimensional nonlinear Ginzburg-Landau equation in complex amplitude A(X,Y,T) with real coefficients near the super critical pitchfork bifurcation. In Section 4.2, we show the occurrence of secondary instabilities such as Eckhaus instability and Zigzag instability. We have also considered the effect of Nusselt number on heat transport by magnetoconvection in a sparsely packed porous medium. In Section 5, we derive two nonlinear one-dimensional coupled Ginzburg-Landau type equations with complex coefficients near the onset of oscillatory convection at a supercritical Hopf bifurcation. Following Matthews and Rucklidge [20], we have dropped slow space dependence in X and obtained two coupled ordinary differential equations in A1R and A1L and discussed the stability regions of travelling and standing waves. By obtaining a one-dimensional Ginzburg-Landau equation in complex amplitude A(X,Y,T) with complex coefficients near a supercritical Hopf bifurcation, we have shown the condition for occurrence of Benjamin-Feir-type instability [21] for travelling and standing waves. In Section 6, we write conclusions of the paper.


## 2. Basic Equations

We consider an electrically and thermally conducting fluid saturating an infinite horizontal layer of a sparsely packed isotropic porous medium of depthd with a uniform magnetic field H0 in the vertical z-direction. This layer is heated from below, the upper and lower bounding surfaces of the layer are assumed to be stress-free. Physical properties of the fluid are assumed to be constant, except for the density in the buoyancy term, so that the Boussinesq approximation is valid. The temperature difference across the stress-free boundaries is ΔT′ and the flow in the sparsely packed porous medium is governed by the Darcy-Lapwood-Brinkman model. The relevant basic equations are(2)∇′⋅V¯′=0,∇′⋅H¯′=0,(3)ρ0′[1ϕ∂V¯′∂t′+1ϕ2(V¯′⋅∇′)V¯′]-μm4π[H0′∂H¯′∂z′+(H¯′⋅∇′)H¯′]=-∇(p′+μm8π|H′|2+μmH04π2Hz′)+ρ′g¯-μKV¯′+μe∇′2V¯′,(4)M∂T′∂t′+(V¯′⋅∇′)T′=κ∇′2T′,(5)ϕ∂H¯′∂t′=∇′×(V¯′×H0′êz)+∇′×(V¯′×H¯′)+η∇′2H¯′.
The fluid density ρ′ is described by(6)ρ′=ρ0′[1-α(T′-Tb′)],
where α=-ρ0′-1(∂ρ′/∂T′) is thermal expansion coefficient and ρ0′ is mean fluid density. Here p′ is pressure, V¯′ is mean fluid velocity, T′ is temperature, H¯′ is magnetic field, ϕ is porosity, g is acceleration due to gravity, K is permeability of porous medium, μe is coefficient of effective fluid viscosity, κ is thermal diffusivity, μm is magnetic permeability, and η is magnetic diffusivity. Equation (3) is known as Darcy-Lapwood-Brinkman equation and is valid for 0.8<ϕ<1. Givler and Altobelli [18] shown that the range of Λ=μe/μ varies from 0.5 to 10.9. M is dimensionless heat capacity and is defined as the ratio of the effective heat capacity of the porous medium to the heat capacity (ρ′Cp)f of the fluid. In a nonporous medium, ϕ=M=Λ=1 and K→∞ and (3) reduces to Navier-Stokes equation. In this paper, for sparsely packed porous medium, we consider M=0.9,ϕ=0.9. The conduction state is characterized by(7)V¯s′=0,Ts′=T0′-(ΔT′d)z′,
and we take the temperature perturbation as θ′=T′-Ts′. We use the scaling(8)x=x′d,y=y′d,z=z′d,t=t′Md2/κ,u=u′κ/Md,v=v′κ/Md,w=w′κ/Md,θ=θ′ΔT′,P=P′ρ0′M-2κ2d-2,H¯=H¯′κH0/η.
Here Md2/κ is thermal diffusion time in a porous medium. Using (6) and (8), we can write basic dimensionless equations for magnetoconvection in a porous medium as(9)∇⋅V¯=0,∇⋅H¯=0,(10)1M2ϕPr1[∂V¯∂t+1ϕ(V¯⋅∇)V¯]-QPr2Pr1(H¯⋅∇)H¯-Q∂H¯∂z=-∇(PMPr1+Q2Pr2Pr1|H¯|2+QHz)-1MDaV¯+ΛM∇2V¯+Rθêz,(11)∂θ∂t+1M(V¯⋅∇)θ=wM+∇2θ,(12)ϕPr2Pr1∂H¯∂t-M∇2H¯=∇×(V¯×êz)+Pr2Pr1∇×(V¯×H¯).
The dimensionless parameters required for the description of the motion are Rayleigh number R=gαΔTd3/κν, thermal Prandtl number Pr1=ν/κ, magnetic Prandtl number Pr2=ν/η, Chandrasekhar number Q=μmH02d2/4πρ0νη, and Darcy number Da=κ/d2. The Curl of (10) gives(13)(1M2ϕPr1∂∂t+1MDa-ΛM∇2)ω¯-Q∂J¯∂z-R(∇×(θêz))=QPr2Pr1[∇×(H¯⋅∇)H¯]-1M2ϕ2Pr1[∇×(V¯⋅∇)V¯],
where vorticity ω¯=∇×V¯, current J¯=∇×H¯ and (14)∇×[(V¯⋅∇)V¯]=[(V¯⋅∇)ω¯-(ω¯⋅∇)V¯],∇×[(H¯⋅∇)H¯]=[(H¯⋅∇)J¯-(J¯⋅∇)H¯].
The Curl of (13) in turn gives, after use of (9),(15)(1M2ϕPr1∂∂t+1MDa-ΛM∇2)∇2V¯-R[(∇2θ)êz-∇(∂θ∂z)]-Q∂∂z(∇2H¯)=1M2ϕ2Pr1∇×[(V¯⋅∇)ω¯-(ω¯⋅∇)V¯]-QPr2Pr1∇×[(H¯⋅∇)J¯-(J¯⋅∇)H¯].
Now taking the scalar product of (12), (13), and (15) with êz, we get,(16)(ϕPr2Pr1∂∂t-M∇2)Hz-∂w∂z=Pr2Pr1êz⋅[∇×(V¯×H¯)].(17)(1M2ϕPr1∂∂t+1MDa-ΛM∇2)∇2ω¯z+Q∂Jz∂z-R∂θ∂x=QPr2Pr1êz⋅∇×(H¯⋅∇)H¯-1M2ϕ2Pr1êz⋅[∇×(V¯⋅∇)V¯],(18)(1M2ϕPr1∂∂t+1MDa-ΛM∇2)∇2w-R∇h2θ-Q∂∂z(∇2Hz)=1M2ϕ2Pr1êz⋅[∇×[(V¯⋅∇)ω¯-(ω¯⋅∇)V¯]]-QPr2Pr1êz⋅[∇×[(H¯⋅∇)J¯-(J¯⋅∇)H¯]].
Geophysically acceptable velocities of propagating instabilities corresponding to geometric scalar variations occur only Pr2/Pr1>1 (where instabilities develop in ohmic diffusion timescale d2/η), Pr2/Pr1=2 and 5, when the turbulent is present in the Earth's outer core. In the case of Pr2/Pr1≪1 the instabilities are extremely slow depending on the thermal diffusion timescale d2/κ. Using (11), (18), and (16) can be reduced to a form(19)Lw=N,
where(20)L=(DϕDPr1-Q∂2∂z2)D∇2-RM∇h2Dϕ,(21)N=QD∇2Pr2Pr1∂∂z[(H¯⋅∇)w-(V¯⋅∇)Hz]+DDϕêz⋅[[1M2ϕ2Pr1∇×[(V¯⋅∇)ω¯-(ω¯⋅∇)V¯]]-[QPr2Pr1∇×[(H¯⋅∇)J¯-(J¯⋅∇)H¯]]1M2ϕ2Pr1]-RM∇h2Dϕ(V¯⋅∇)θ,
here(22)D=(∂∂t-∇2),Dϕ=(ϕPr2Pr1∂∂t-M∇2),DPr1=(1M2ϕPr1∂∂t+1MDa-ΛM∇2),∇h2=∂2∂x2,∇2=(∂2∂x2+∂2∂z2).Boundary Conditions
We assume that fluid is contained betweenz=0 and z=1, where z=0 corresponds to boundary of solid iron core with Earth's mushy layer and z=1 corresponds to boundary of Earth's mushy layer with Earth's outer liquid core. For perfectly conducting boundary with temperature, we have
(23)θ=0,Hz=0on⁡z=0,z=1∀x,y.
Also the normal component of the velocity would vanish on z=0, z=1, that is,
(24)w=0on⁡z=0,z=1∀x,y.
However, there are two more conditions to be imposed on velocity depending on the nature of the surface. In this paper we consider free-free boundary conditions, that is, on surfaces the tangential stresses vanish, which is equivalent to
(25)Pxz=μ(∂u∂z+∂w∂x)=0,Pyz=μ(∂v∂z+∂w∂y)=0,
where μ=γρ0 is dynamic viscosity. Since w vanishes for x,y on z=0, z=1, it follows that ∂u/∂z=∂v/∂z=0 on a free surface z=0,z=1. Hence from equation of continuity we have ∂2w/∂z2=0 on z=0,z=1 for all x, y. In this paper we have considered only the idealized stress-free conditions on the surface and vanishing of temperature fluctuations. Thus w=D2w=D4w=0 at z=0,1. w and its even derivatives vanish at z=0 and z=1.


## 3. Linear Stability Analysis

We perform a linear stability analysis of the problem by substituting(26)w=W(z)eiqx+pt,
into linearized version of (19) is ℒw=0, and obtaining an equation(27)[ΛMpM2ϕPr1(D2-q2-p)(MD2-Mq2-pϕPr2Pr1)(D2-q2)×{ΛM(D2-q2)-1MDa-pM2ϕPr1}]W=[-Rq2M(MD2-Mq2-pϕPr2Pr1)+Q(D2-q2)(D2-q2-p)D2-Rq2M]W.
We consider stress-free boundary conditions, then W=D2W=0 on z=0, z=1 for all x,y. Thus we can assume W=sinπz.SubstitutingW=sinπz and p=iω into (27), we get(28)R=Mq2[A1+iω(A2ω2+A3)],(29)A1=K[δϕ2Pr2MDaPr1{ΛMPr2Pr1-ω2ϕΛMPr2Pr1ΛMPr2Pr1Mδ6(1Da+Λδ2)+Mδ4×(Qπ2-ω2Pr2M2Pr12-ω2MϕPr1-ω2ϕΛMPr2Pr1)}+ω2Pr2Pr1{Pr2Pr1(δ4ϕ2ΛM-ω2MϕPr12)×δ4(Λϕ+1MPr1)+Qπ2ϕ+δϕ2Pr2MDaPr1}],(30)A2=K[δ2(Pr2Pr1)2(ϕ2ΛM+ϕM2Pr1)+ϕ2Pr22MDaPr12],(31)A3=K[δ4(MΛ+1ϕPr1)+Mδ2Da+Qπ2(M-ϕPr2Pr1)],
where 𝒦=δ2(M4δ4+ω2ϕ2Pr22/Pr12)-1 and δ2=(π2+q2), from relation equation (30), A2>0.


### 3.1. Stationary Convection (ω=0)

Substitutingω=0 in (28), we get(32)Rs=δs2qs2[δs2(1Da+δs2Λ)+Qπ2],
here Rs is the value of the Rayleigh number for stationary convection. The minimum value of Rs is obtained for qs=qsc. where(33)2Λ(qπ)6+(3Λ+1π2DDa)(qπ)4=Λ+Qπ2+1π2DDa.
The wave number is identical to that for the single component fluid, while the threshold for the onset of stationary convection at pitchfork bifurcation is given by (34) with qs=qsc,(34)Rsc=δsc2qsc2[δsc2(1Da+δsc2Λ)+Qπ2],
where δsc2=π2+qsc2. Thus the magnetic field inhibits the onset of stationary convection.


### 3.2. Oscillatory Convection (ω2>0)

For the oscillatory convection(ω≠0) and from (28), R will be complex. But the physical meaning of R requires it to be real. The condition that R is real implies that imaginary part of (28) is zero, that is,(35)A2ω2+A3=0,
where A2 and A3 are given by (30) and (31). For oscillatory convection ω2=-A3/A2>0 since A2>0, for oscillatory convection A3<0. For A3=0, (35) implies that ω=0 is a double zero corresponding to Takens-Bogdanov bifurcation point. For oscillatory convection, we have(36)ω2=-[δo4(MΛ+1/ϕPr1)+Mδo2/Da+Qπ2(M-ϕPr2/Pr1)][δo2(Pr2/Pr1)2(ϕ2Λ/M+ϕ/M2Pr1)+ϕ2Pr22/MDaPr12],where δo2=π2+qo2. A necessary condition for ω2>0 is(37)Pr2Pr1>1ϕ.
However, this is not sufficient condition and one must have in addition(38)Q>δo4(MΛ+1/ϕPr1)+Mδo2/Daπ2(ϕPr2/Pr1-M).
At Takens-Bogdanov bifurcation point Ro(qo)=Rs(qs)=Rc(qc), qo=qs=qc, and ω2=0 is a double zero at Q=Qc(qc) where(39)Q=δc4(MΛ+1/ϕPr1)+Mδc2/Daπ2(ϕPr2/Pr1-M),q=qc.
The Takens-Bogdanov bifurcation point occurs when the neutral curves for Hopf and pitchfork bifurcation meet and only a single wave number is present, namely, qo=qs=qc. If qc>qsc then for all q<qc the first instability to set in is an oscillatory convection. The asymptotic values of qc and qsc for large Chandrasekhar number (Q→∞) are(40)qc⟶[(ϕPr2-MPr1)Qπ2(MΛPr1+1/ϕ)]1/4,qsc⟶(Qπ4M2Λ)1/6.
From the monotonic dependence of qc and qsc on Q, we may conclude that for Pr2>Pr1, there exists a Q(Pr1,Pr2,M,Λ,ϕ,Da) such that for Q<Q(Pr1,Pr2,M,Λ,ϕ,Da) the onset of first instability will be stationary convection at pitchfork bifurcation while for Q>Q(Pr1,Pr2,M,Λ,ϕ,Da) it will be oscillatory convection at Hopf bifurcation. Q(Pr1,Pr2,M,Λ,ϕ,Da) and for Q=Q(Pr1,Pr2,M,Λ,ϕ,Da), we have(41)Rct=Roc(qoc)=Rsc(qsc)butqoc≠qsc,
above condition (41) gives codimension-two bifurcation point. However, there is no simple formula to give Q(Pr1,Pr2,M,Λ,ϕ,Da) at the codimension-two bifurcation point by assuming R as an independent variable, such kind of interesting result is not available in Chandrasekhar [2]. In Figures 1 and 2, each solid line stands for stationary convection (pitchfork bifurcation) and dotted line stands for oscillatory convection (Hopf bifurcation). In Figures 1 and 2, we have showed the effect of several physical parameters, like Q, Pr1, Pr2, Λ, M, ϕ, and Da on the onset of both stationary convection and oscillatory convection when a physical parameter increases for the remaining fixed parameters, the onset of instabilities increases, that is, the onset of stationary convection and oscillatory convection inhibit when a parameter increases with the remaining fixed parameters.Numerically calculated marginal stability curves are plotted in(R,q)-plane for Pr1=1, Pr2=1.5, Da=1500, Λ=0.85, ϕ=0.9, and M=0.9,  (a)  Q=104, (b)  Q=106, (c)  Q=108, and  (d)  Q=1010, then the onset of stationary convection and the onset of oscillatory convection increase (stationary convection represented by solid lines and oscillatory convection represented by dotted lines).
(a)(b)(c)(d)Neutral curves for the stationary bifurcation (solid lines) and for the Hopf bifurcation (dashed lines) near the codimension two point forQ=2000, Da=1500, Λ=0.85, ϕ=0.9, Pr1=1, and M=0.9, (a)  Pr2=1.27, (b)  Pr2=1.3, (c)  Pr2=1.35  x-axis wave number, y-Rayleigh numbers Rs,Ro.
(a)(b)(c)


## 3.1. Stationary Convection (ω=0)

Substitutingω=0 in (28), we get(32)Rs=δs2qs2[δs2(1Da+δs2Λ)+Qπ2],
here Rs is the value of the Rayleigh number for stationary convection. The minimum value of Rs is obtained for qs=qsc. where(33)2Λ(qπ)6+(3Λ+1π2DDa)(qπ)4=Λ+Qπ2+1π2DDa.
The wave number is identical to that for the single component fluid, while the threshold for the onset of stationary convection at pitchfork bifurcation is given by (34) with qs=qsc,(34)Rsc=δsc2qsc2[δsc2(1Da+δsc2Λ)+Qπ2],
where δsc2=π2+qsc2. Thus the magnetic field inhibits the onset of stationary convection.


## 3.2. Oscillatory Convection (ω2>0)

For the oscillatory convection(ω≠0) and from (28), R will be complex. But the physical meaning of R requires it to be real. The condition that R is real implies that imaginary part of (28) is zero, that is,(35)A2ω2+A3=0,
where A2 and A3 are given by (30) and (31). For oscillatory convection ω2=-A3/A2>0 since A2>0, for oscillatory convection A3<0. For A3=0, (35) implies that ω=0 is a double zero corresponding to Takens-Bogdanov bifurcation point. For oscillatory convection, we have(36)ω2=-[δo4(MΛ+1/ϕPr1)+Mδo2/Da+Qπ2(M-ϕPr2/Pr1)][δo2(Pr2/Pr1)2(ϕ2Λ/M+ϕ/M2Pr1)+ϕ2Pr22/MDaPr12],where δo2=π2+qo2. A necessary condition for ω2>0 is(37)Pr2Pr1>1ϕ.
However, this is not sufficient condition and one must have in addition(38)Q>δo4(MΛ+1/ϕPr1)+Mδo2/Daπ2(ϕPr2/Pr1-M).
At Takens-Bogdanov bifurcation point Ro(qo)=Rs(qs)=Rc(qc), qo=qs=qc, and ω2=0 is a double zero at Q=Qc(qc) where(39)Q=δc4(MΛ+1/ϕPr1)+Mδc2/Daπ2(ϕPr2/Pr1-M),q=qc.
The Takens-Bogdanov bifurcation point occurs when the neutral curves for Hopf and pitchfork bifurcation meet and only a single wave number is present, namely, qo=qs=qc. If qc>qsc then for all q<qc the first instability to set in is an oscillatory convection. The asymptotic values of qc and qsc for large Chandrasekhar number (Q→∞) are(40)qc⟶[(ϕPr2-MPr1)Qπ2(MΛPr1+1/ϕ)]1/4,qsc⟶(Qπ4M2Λ)1/6.
From the monotonic dependence of qc and qsc on Q, we may conclude that for Pr2>Pr1, there exists a Q(Pr1,Pr2,M,Λ,ϕ,Da) such that for Q<Q(Pr1,Pr2,M,Λ,ϕ,Da) the onset of first instability will be stationary convection at pitchfork bifurcation while for Q>Q(Pr1,Pr2,M,Λ,ϕ,Da) it will be oscillatory convection at Hopf bifurcation. Q(Pr1,Pr2,M,Λ,ϕ,Da) and for Q=Q(Pr1,Pr2,M,Λ,ϕ,Da), we have(41)Rct=Roc(qoc)=Rsc(qsc)butqoc≠qsc,
above condition (41) gives codimension-two bifurcation point. However, there is no simple formula to give Q(Pr1,Pr2,M,Λ,ϕ,Da) at the codimension-two bifurcation point by assuming R as an independent variable, such kind of interesting result is not available in Chandrasekhar [2]. In Figures 1 and 2, each solid line stands for stationary convection (pitchfork bifurcation) and dotted line stands for oscillatory convection (Hopf bifurcation). In Figures 1 and 2, we have showed the effect of several physical parameters, like Q, Pr1, Pr2, Λ, M, ϕ, and Da on the onset of both stationary convection and oscillatory convection when a physical parameter increases for the remaining fixed parameters, the onset of instabilities increases, that is, the onset of stationary convection and oscillatory convection inhibit when a parameter increases with the remaining fixed parameters.Numerically calculated marginal stability curves are plotted in(R,q)-plane for Pr1=1, Pr2=1.5, Da=1500, Λ=0.85, ϕ=0.9, and M=0.9,  (a)  Q=104, (b)  Q=106, (c)  Q=108, and  (d)  Q=1010, then the onset of stationary convection and the onset of oscillatory convection increase (stationary convection represented by solid lines and oscillatory convection represented by dotted lines).
(a)(b)(c)(d)Neutral curves for the stationary bifurcation (solid lines) and for the Hopf bifurcation (dashed lines) near the codimension two point forQ=2000, Da=1500, Λ=0.85, ϕ=0.9, Pr1=1, and M=0.9, (a)  Pr2=1.27, (b)  Pr2=1.3, (c)  Pr2=1.35  x-axis wave number, y-Rayleigh numbers Rs,Ro.
(a)(b)(c)


## 4. Onset of Stationary Convection at Supercritical Pitchfork Bifurcation

### 4.1. Derivation of Two-Dimensional Nonlinear Ginzburg-Landau Equation Using Newell-Whitehead [19] Method

In this section the evolution of a general pattern is developed by means of a multiple scale analysis used by Newell and Whitehead [19]. A small amplitude convection cell is imposed on the basic flow. If this amplitude is of the size O(ϵ) then the interaction of the cell with itself forces a second harmonic and mean state correction of size O(ϵ2) and then in turn drives an O(ϵ3) correction to the fundamental component of the imposed roll. A solvability criteria for this correction yields the one-dimensional nonlinear Ginzburg-Landau equation of the complex valued amplitude A(X,Y,T) of the imposed disturbance with real coefficients. To simplify the problem we assume the formulation of cylindrical rolls with axis parallel to y-axis, so that y-dependence disappears from (19). The z-dependence is contained entirely in the sine and cosine functions, which ensures that stress-free boundary conditions are satisfied. We use the expansion parameter ϵ as(42)ϵ2=R-RscRsc.
For the values of R close to threshold value Rsc that is, ϵ≪1, the structure of the slow length scales will be insensitive to ϵ, but a slow modulation in space and time is possible by making use of the band of the unstable solutions and linear growth rate is likely to saturate due to nonlinear effects. This behavior can be analyzed by writing solutions of (9)–(12) in power series ϵ as(43)f=ϵf0+ϵ2f1+ϵ3f2+⋯,
where f=f(u,v,w,θ,Hx,Hy,Hz) with the first approximation is given by the eigenvector of the linearized problem:(44)u0=iπqsc[A(X,Y,T)eiqscxcos⁡πz-c.c.],v0=0,w0=A(X,Y,T)eiqscxsin⁡πz+c.c.,θ0=1Mδsc2[A(X,Y,T)eiqscxsin⁡πz+c.c.],Hx0=-iπ2Mqscδsc2[A(X,Y,T)eiqscxsin⁡πz-c.c.],Hy0=0,Hz0=πMδsc2[A(X,Y,T)eiqscxcos⁡πz+c.c.],
where δsc2=π2+qsc2, here C.C. stands for complex conjugate, eiqscxsinπz is the critical mode for the linear problem at R=Rsc and q=qsc. The complex amplitude A(X,Y,T) depends on the slow variables X, Y, Z, and T to be scaled by introducing multiple scales(45)X=ϵx,Y=ϵ1/2y,Z=z,T=ϵ2t,
and these formally separate the fast and slow dependent variables in f. It should be noted that the difference in scaling in the two directions reflects the inherent symmetry breaking of instability which was chosen here with wave vector in x-direction. The differential operators can be expressed as(46)∂∂x⟶∂∂x+ϵ∂∂X,∂∂y⟶ϵ1/2∂∂Y,∂∂z⟶∂∂Z,∂∂t⟶ϵ2∂∂T
with the assumption (46), the operators (20) and (21) are transformed into a set of linear in homogeneous equations. The solvability conditions for the latter yields the amplitude equation using (44) in the linear operator (20) can be written as(47)L=L0+ϵL1+ϵ2L2+⋯,
where(48)L0=-Λ∇8+1Da∇6+Q∇4∂2∂z2+Rsc∇2∂2∂x2,(49)L1=(2∂2∂x∂X+∂2∂Y2)×[3Da∇4-4Λ∇6+2Q∇2∂2∂z2+Rsc∇2+Rsc∂2∂x2],(50)L2=∂∂T[(Λ+1MϕPr1+ϕΛMPr2Pr1)∇6-(1Da+ϕPr2Pr11MDa)∇4-Q∇2∂2∂z2-RscMϕPr2Pr1∂2∂x2]+∂2∂X2×[3Da∇4-4Λ∇6+2Q∇2∂2∂z2+Rsc∇2+Rsc∂2∂x2]+(4∂4∂x2∂X2+12∂2∂Y2)×[-6Λ∇4+3Da∇2+Q∂2∂z2+Rsc].
Similarly nonlinear term 𝒩 is given by(51)N=ϵ2N0+ϵ3N1+⋯,
substituting (47), (51), and (43) into (19), we get by equating the coefficients of ϵ,ϵ2,ϵ3;(52)L0w0=0,(53)L0w1+L1w0=N0,(54)L0w2+L1w1+L2w0=N1.
Equation (48) gives the critical Rayleigh number for the onset of stationary convection(55)Rsc=δsc2qsc2[δsc4Λ+1Daδsc2+Qπ2].
In (53), 𝒩0=0, ℒ1w0=0 and hence w1=0. From equation of continuity we find that u1=0. The relevant equations for θ1 and Hz1 are(56)(∂∂t-∇2)θ1=w1M-1M[u0∂θ0∂x+w0∂θ0∂z],
form (56) and (44), we get(57)θ1=-12πM2δsc2|A|2sin⁡2πz.
Equation (12) gives relevant equation for Hz1 as(58)(ϕPr2Pr1-M∇2)Hz1=∂w1∂z+Pr2Pr1∂∂x[w0Hx0-u0Hz0].
From (58) and (44), we get(59)Hz1=Pr2Pr1⋅π24δsc2qsc2[A2e2iqscx+c.c.].
Similarly we have Hx1=0,Hy1=0. The solvability criterion of (54) gives the amplitude equation which can be written as(60)λ0∂A∂T-λ1(∂∂X-i2qsc∂2∂Y2)2-λ2A+λ3|A|2A=0,
where(61)λ0=(Λ+1MϕPr1+ϕPr2Pr1ΛM)δsc6+(1Da+ϕPr2Pr11MDa)δsc4+Qπ2δsc2-RscMPr2Pr1qsc2ϕ,λ1=4qsc2[6Λδsc4+3Daδsc2+Qπ2-Rsc],λ2=Rscqsc2δsc2,λ3=QPr22Pr12π4M2qsc2(qsc2-π2)+Rsc2M2qsc2.
Equation (60) is two-dimensional, nonlinear time-dependent Ginzburg-Landau equation describing the effect of magneticfield in a sparsely packed porous medium near the onset of stationary convection at supercritical pitchfork bifurcation. Here λ0 is always positive for Pr2/Pr1<1/ϕ and for any Q but if Pr2/Pr1>1/ϕ then λ0 is positive only if Q<Qc. Thus for supercritical pitchfork bifurcation λ0 is always positive. For Pr2/Pr1>1/ϕ, λ0 decreases as Q increases and becomes zero at Q=Qc. λ1 and λ2 are always positive. λ3 is positive only if(62)Q<Rscqsc42π4(π2-qsc2)⋅Pr12Pr22.
The pitchfork bifurcation is supercritical if λ3>0 and subcritical if λ3<0. At λ3=0, we get tricritical bifurcation point [22] (see Figure 3). Dropping the time dependence from (60), we get(63)d2AdX2+λ2λ1(1-λ3λ2|A|2)A=0,
since λ1>0, the solution of (63) is given by(64)A(X)=A0tanh(XΛ1),
where(65)A0=(λ2λ3)1/2,Λ1=(2λ1λ2)1/2.Figure 3
Above figure is plotted forDa=1500, Λ=0.85, ϕ=0.9, M=0.9, and Pr2=1.λ3 is the nonlinear coefficient of Ginzburg-Landau equation at the onset of stationary convection. The pitchfork bifurcation is supercritical if λ3>0, subcritical if λ3<0 and λ3=0 on the curve.


### 4.2. Long Wavelength Instabilities (Secondary Instabilities)

The secondary Instabilities arising in nonequilibrium systems do not exhibit strict symmetries but may show spatially slow deformations of the cellular structures. Further, there are secondary instabilities like Eckhaus and Zigzag instabilities, such phenomena are studied using evolution equations for amplitudes which are slowly varying in time as well as in space. These envelope equations can be derived by the method of Newell and Whitehead [19]. The two-dimensional Ginzburg-Landau equation (60), can be written in fast variables x,y,t, and A(X,Y,T)=A(x,y,t)/ϵ, as(66)λ0∂A∂t-λ1(∂∂x-i2qsc∂2∂y2)2A-ϵ2λ2A+λ3|A|2A=0.
In order to study the properties of a structure with a given phase winding number δk, we substitute(67)A(x,y,t)=A1(x,y,t)eiδkx,
into (66) and we obtain(68)λ0∂A1∂t=(ϵ2λ2-λ1(δk)2)A1+2iλ1δk(∂∂x-i2qsc∂2∂y2)A1+λ1(∂∂x-i2qsc∂2∂y2)2A1-λ3|A1|2A1=0.
The steady-state uniform solution of (68) is(69)A1=A1o=[(ϵ2λ2-λ1(δk)2)λ3-1]1/2.
Let ũ(x,y,t)+iṽ(x,y,t) be an infinitesimal perturbation from a uniform steady-state solution A1o given by (69). Now substituting (70)A1=A1o=[(ϵ2λ2-λ1(δk)2)λ3-1]1/2+ũ+iṽ,
into (68) and equating real and imaginary parts, we obtain (71)λ0∂ũ∂t=[-14qsc2∂4∂y4-2(ϵ2λ2-λ1(δk))2+λ1(∂2∂x2+δkqsc∂2∂y2-14qsc2∂4∂y4)]ũ-(2λ1δk-λ1qsc∂2∂y2)∂ṽ∂x,λ0∂ṽ∂t=(2λ1δk-λ1qsc∂2∂y2)∂ũ∂x+λ1(∂2∂x2+δkqsc∂2∂y2-14qsc2∂4∂y4)ṽ.
We analyze (71) by using normal modes of the form(72)ũ=UeStcos⁡(qxx)cos⁡(qyy),ṽ=VeStsin⁡(qxx)cos⁡(qyy).
Putting (72) in (71) we get, (73)[λ0S+2(ϵ2λ2-λ1(δk)2)+χ1]U+λ1qxχ2V=0,λ1qxχ2U+(λ0S+χ1)V=0.
Here χ1=λ1[qx2+(qy2δk)/qsc+qy4/4qsc2], χ2=(2δk+qy2/qsc). On solving (73) we get,(74)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0χ1]+[2(ϵ2λ2-λ1(δk)2)+χ1]ψ1-qx2λ1χ2=0,
whose roots (S±) are real. Here (S±) is defined as(75)(S±)=-1λ02{(2λ0(ϵ2λ2-λ1(δk)2)2+λ12qx2χ22)1/2(2λ0(ϵ2λ2-λ1(δk)2)+λ0χ1)±(2λ0(ϵ2λ2-λ1(δk)2)2+λ12qx2χ22)1/2}.
Solution S(-) is clearly negative, thus the corresponding mode is stable and if S(+) is positive then rolls can be unstable. Symmetry considerations help us to restrict the study of S(+) to a domain qx≥0,qy≥0.


#### 4.2.1. Longitudinal Perturbations and Eckhaus Instability

Insertingqy=0 into (75), we get(76)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0λ1qx2]+λ1qx2[2(ϵ2λ2-3λ1(δk)2)+qx2]=0,
since the roots are real and their sum always negative, the pattern is stable as long as both roots are negative, that is, their product is positive. The cell pattern becomes unstable when the product is negative, that is, when  (77)qx2≤2(3λ1δk2-ϵ2λ2),
for this requires |δk|≥(ϵ2λ2/3λ1), this condition defines the domain of Eckhaus instability. The above condition implies that the most unstable wave vector tends to zero, when |δk|→(ϵ2λ2/3λ1).


#### 4.2.2. Transverse Perturbations and Zigzag Instability

Let us considerqx=0 into (75), we get(78)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0χ1y]+[2(ϵ2λ2-λ1(δk)2)+χ1y]χ1y=0,
where χ1y=λ1(qy2δk/qsc+qy4/4qsc2). The two eigenmodes are uncoupled and we have S(-),(79)S(-)=-2(ϵ2λ2-λ1(δk)2)-λ1qscδkqy2-λ14qsc2qy2<0,
for one of them. The other is amplified when  (80)S(+)=-λ1qy2(δk+qy24qsc)>0.
This implies that δk<0, the above condition defines the domain of the Zigzag Instability. When δk→0 from below the wave vector qy of the instability also tends to zero while the growth rate varies as qy2. We have studied the effect of magnetic field on long wavelength instabilities. We have observed that Eckhaus instability and Zigzag instability regions increases when Q increases (see Figure 4).Figure 4
Numerically computed secondary instability regions of Eckhaus instability (E), Zigzag instability (Z), and stable regions (S) are plotted in(λ2/λ1,δqs)-plane for Q=2000, Da=1500, Λ=0.85, ϕ=0.9, M=0.9, Pr1=1, and Pr2=2. As |δqs| increases then the secondary instability regions increases.


### 4.3. Heat Transport by Convection

The maximum of steady amplitudeA is denoted by |Amax| which is given as(81)|Amax⁡|=(ϵ2λ2λ3)1/2.
Equation (81) is obtained from (64) with tanh(X/Λ1)=1. We use |Amax| to calculate Nusselt number Nu. To discuss the heat transfer near the neutral region, we express it through the Nusselt number is defined as Nu=Hd/κΔT, which is the ratio of the heat transported across any layer to the heat which would be transported by conduction alone. Here H is the rate of heat transfer per unit area and is defined as(82)H=-〈∂Ttotal∂z′〉z′=0.
In (82), angular brackets correspond to a horizontal average. The Nusselt number Nu can be calculated in terms of amplitude A and is given as(83)Nu=1+ϵ2δsc2|Amax⁡|2.
From (83), we get conduction for R≤Rsc and convection for R>Rsc. Since the amplitude equation is valid for λ3>0, which is possible for R>Rsc (supercritical pitchfork bifurcation), thus we get Nu>1 for R>Rsc. We get convection for Nu>1 and conduction for Nu≤1. In stationary convection Nu increases implies that heat conducted by steady mode increases. In the problem of double diffusive convection in porous medium with magnetic field, Nu depends on Pr1,Pr2,Λ,M,ϕ,Da, and Q. We have computed Nu for different values of Q, for some fixed values of other parameters and observed that Nu increases as Q decreases (see Figures 5(a) and 5(b)). This implies that magnetic field inhibits the heat transport. The parameters Pr1, Pr2, Λ, M, ϕ, and Da show the same result as Q shows on Nu.Graph (a) is plotted forQ=1000 and graph (b) is plotted for Q=3000 for the fixed values of Da=1500, Λ=0.85, ϕ=0.9, Pr2=1, Pr1=2, and M=0.9. in (Nu,R/Rsc)-plane. In graphs (a) and (b), as R/Rsc increases then Nu increases.
(a)(b)


## 4.1. Derivation of Two-Dimensional Nonlinear Ginzburg-Landau Equation Using Newell-Whitehead [19] Method

In this section the evolution of a general pattern is developed by means of a multiple scale analysis used by Newell and Whitehead [19]. A small amplitude convection cell is imposed on the basic flow. If this amplitude is of the size O(ϵ) then the interaction of the cell with itself forces a second harmonic and mean state correction of size O(ϵ2) and then in turn drives an O(ϵ3) correction to the fundamental component of the imposed roll. A solvability criteria for this correction yields the one-dimensional nonlinear Ginzburg-Landau equation of the complex valued amplitude A(X,Y,T) of the imposed disturbance with real coefficients. To simplify the problem we assume the formulation of cylindrical rolls with axis parallel to y-axis, so that y-dependence disappears from (19). The z-dependence is contained entirely in the sine and cosine functions, which ensures that stress-free boundary conditions are satisfied. We use the expansion parameter ϵ as(42)ϵ2=R-RscRsc.
For the values of R close to threshold value Rsc that is, ϵ≪1, the structure of the slow length scales will be insensitive to ϵ, but a slow modulation in space and time is possible by making use of the band of the unstable solutions and linear growth rate is likely to saturate due to nonlinear effects. This behavior can be analyzed by writing solutions of (9)–(12) in power series ϵ as(43)f=ϵf0+ϵ2f1+ϵ3f2+⋯,
where f=f(u,v,w,θ,Hx,Hy,Hz) with the first approximation is given by the eigenvector of the linearized problem:(44)u0=iπqsc[A(X,Y,T)eiqscxcos⁡πz-c.c.],v0=0,w0=A(X,Y,T)eiqscxsin⁡πz+c.c.,θ0=1Mδsc2[A(X,Y,T)eiqscxsin⁡πz+c.c.],Hx0=-iπ2Mqscδsc2[A(X,Y,T)eiqscxsin⁡πz-c.c.],Hy0=0,Hz0=πMδsc2[A(X,Y,T)eiqscxcos⁡πz+c.c.],
where δsc2=π2+qsc2, here C.C. stands for complex conjugate, eiqscxsinπz is the critical mode for the linear problem at R=Rsc and q=qsc. The complex amplitude A(X,Y,T) depends on the slow variables X, Y, Z, and T to be scaled by introducing multiple scales(45)X=ϵx,Y=ϵ1/2y,Z=z,T=ϵ2t,
and these formally separate the fast and slow dependent variables in f. It should be noted that the difference in scaling in the two directions reflects the inherent symmetry breaking of instability which was chosen here with wave vector in x-direction. The differential operators can be expressed as(46)∂∂x⟶∂∂x+ϵ∂∂X,∂∂y⟶ϵ1/2∂∂Y,∂∂z⟶∂∂Z,∂∂t⟶ϵ2∂∂T
with the assumption (46), the operators (20) and (21) are transformed into a set of linear in homogeneous equations. The solvability conditions for the latter yields the amplitude equation using (44) in the linear operator (20) can be written as(47)L=L0+ϵL1+ϵ2L2+⋯,
where(48)L0=-Λ∇8+1Da∇6+Q∇4∂2∂z2+Rsc∇2∂2∂x2,(49)L1=(2∂2∂x∂X+∂2∂Y2)×[3Da∇4-4Λ∇6+2Q∇2∂2∂z2+Rsc∇2+Rsc∂2∂x2],(50)L2=∂∂T[(Λ+1MϕPr1+ϕΛMPr2Pr1)∇6-(1Da+ϕPr2Pr11MDa)∇4-Q∇2∂2∂z2-RscMϕPr2Pr1∂2∂x2]+∂2∂X2×[3Da∇4-4Λ∇6+2Q∇2∂2∂z2+Rsc∇2+Rsc∂2∂x2]+(4∂4∂x2∂X2+12∂2∂Y2)×[-6Λ∇4+3Da∇2+Q∂2∂z2+Rsc].
Similarly nonlinear term 𝒩 is given by(51)N=ϵ2N0+ϵ3N1+⋯,
substituting (47), (51), and (43) into (19), we get by equating the coefficients of ϵ,ϵ2,ϵ3;(52)L0w0=0,(53)L0w1+L1w0=N0,(54)L0w2+L1w1+L2w0=N1.
Equation (48) gives the critical Rayleigh number for the onset of stationary convection(55)Rsc=δsc2qsc2[δsc4Λ+1Daδsc2+Qπ2].
In (53), 𝒩0=0, ℒ1w0=0 and hence w1=0. From equation of continuity we find that u1=0. The relevant equations for θ1 and Hz1 are(56)(∂∂t-∇2)θ1=w1M-1M[u0∂θ0∂x+w0∂θ0∂z],
form (56) and (44), we get(57)θ1=-12πM2δsc2|A|2sin⁡2πz.
Equation (12) gives relevant equation for Hz1 as(58)(ϕPr2Pr1-M∇2)Hz1=∂w1∂z+Pr2Pr1∂∂x[w0Hx0-u0Hz0].
From (58) and (44), we get(59)Hz1=Pr2Pr1⋅π24δsc2qsc2[A2e2iqscx+c.c.].
Similarly we have Hx1=0,Hy1=0. The solvability criterion of (54) gives the amplitude equation which can be written as(60)λ0∂A∂T-λ1(∂∂X-i2qsc∂2∂Y2)2-λ2A+λ3|A|2A=0,
where(61)λ0=(Λ+1MϕPr1+ϕPr2Pr1ΛM)δsc6+(1Da+ϕPr2Pr11MDa)δsc4+Qπ2δsc2-RscMPr2Pr1qsc2ϕ,λ1=4qsc2[6Λδsc4+3Daδsc2+Qπ2-Rsc],λ2=Rscqsc2δsc2,λ3=QPr22Pr12π4M2qsc2(qsc2-π2)+Rsc2M2qsc2.
Equation (60) is two-dimensional, nonlinear time-dependent Ginzburg-Landau equation describing the effect of magneticfield in a sparsely packed porous medium near the onset of stationary convection at supercritical pitchfork bifurcation. Here λ0 is always positive for Pr2/Pr1<1/ϕ and for any Q but if Pr2/Pr1>1/ϕ then λ0 is positive only if Q<Qc. Thus for supercritical pitchfork bifurcation λ0 is always positive. For Pr2/Pr1>1/ϕ, λ0 decreases as Q increases and becomes zero at Q=Qc. λ1 and λ2 are always positive. λ3 is positive only if(62)Q<Rscqsc42π4(π2-qsc2)⋅Pr12Pr22.
The pitchfork bifurcation is supercritical if λ3>0 and subcritical if λ3<0. At λ3=0, we get tricritical bifurcation point [22] (see Figure 3). Dropping the time dependence from (60), we get(63)d2AdX2+λ2λ1(1-λ3λ2|A|2)A=0,
since λ1>0, the solution of (63) is given by(64)A(X)=A0tanh(XΛ1),
where(65)A0=(λ2λ3)1/2,Λ1=(2λ1λ2)1/2.Figure 3
Above figure is plotted forDa=1500, Λ=0.85, ϕ=0.9, M=0.9, and Pr2=1.λ3 is the nonlinear coefficient of Ginzburg-Landau equation at the onset of stationary convection. The pitchfork bifurcation is supercritical if λ3>0, subcritical if λ3<0 and λ3=0 on the curve.


## 4.2. Long Wavelength Instabilities (Secondary Instabilities)

The secondary Instabilities arising in nonequilibrium systems do not exhibit strict symmetries but may show spatially slow deformations of the cellular structures. Further, there are secondary instabilities like Eckhaus and Zigzag instabilities, such phenomena are studied using evolution equations for amplitudes which are slowly varying in time as well as in space. These envelope equations can be derived by the method of Newell and Whitehead [19]. The two-dimensional Ginzburg-Landau equation (60), can be written in fast variables x,y,t, and A(X,Y,T)=A(x,y,t)/ϵ, as(66)λ0∂A∂t-λ1(∂∂x-i2qsc∂2∂y2)2A-ϵ2λ2A+λ3|A|2A=0.
In order to study the properties of a structure with a given phase winding number δk, we substitute(67)A(x,y,t)=A1(x,y,t)eiδkx,
into (66) and we obtain(68)λ0∂A1∂t=(ϵ2λ2-λ1(δk)2)A1+2iλ1δk(∂∂x-i2qsc∂2∂y2)A1+λ1(∂∂x-i2qsc∂2∂y2)2A1-λ3|A1|2A1=0.
The steady-state uniform solution of (68) is(69)A1=A1o=[(ϵ2λ2-λ1(δk)2)λ3-1]1/2.
Let ũ(x,y,t)+iṽ(x,y,t) be an infinitesimal perturbation from a uniform steady-state solution A1o given by (69). Now substituting (70)A1=A1o=[(ϵ2λ2-λ1(δk)2)λ3-1]1/2+ũ+iṽ,
into (68) and equating real and imaginary parts, we obtain (71)λ0∂ũ∂t=[-14qsc2∂4∂y4-2(ϵ2λ2-λ1(δk))2+λ1(∂2∂x2+δkqsc∂2∂y2-14qsc2∂4∂y4)]ũ-(2λ1δk-λ1qsc∂2∂y2)∂ṽ∂x,λ0∂ṽ∂t=(2λ1δk-λ1qsc∂2∂y2)∂ũ∂x+λ1(∂2∂x2+δkqsc∂2∂y2-14qsc2∂4∂y4)ṽ.
We analyze (71) by using normal modes of the form(72)ũ=UeStcos⁡(qxx)cos⁡(qyy),ṽ=VeStsin⁡(qxx)cos⁡(qyy).
Putting (72) in (71) we get, (73)[λ0S+2(ϵ2λ2-λ1(δk)2)+χ1]U+λ1qxχ2V=0,λ1qxχ2U+(λ0S+χ1)V=0.
Here χ1=λ1[qx2+(qy2δk)/qsc+qy4/4qsc2], χ2=(2δk+qy2/qsc). On solving (73) we get,(74)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0χ1]+[2(ϵ2λ2-λ1(δk)2)+χ1]ψ1-qx2λ1χ2=0,
whose roots (S±) are real. Here (S±) is defined as(75)(S±)=-1λ02{(2λ0(ϵ2λ2-λ1(δk)2)2+λ12qx2χ22)1/2(2λ0(ϵ2λ2-λ1(δk)2)+λ0χ1)±(2λ0(ϵ2λ2-λ1(δk)2)2+λ12qx2χ22)1/2}.
Solution S(-) is clearly negative, thus the corresponding mode is stable and if S(+) is positive then rolls can be unstable. Symmetry considerations help us to restrict the study of S(+) to a domain qx≥0,qy≥0.


### 4.2.1. Longitudinal Perturbations and Eckhaus Instability

Insertingqy=0 into (75), we get(76)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0λ1qx2]+λ1qx2[2(ϵ2λ2-3λ1(δk)2)+qx2]=0,
since the roots are real and their sum always negative, the pattern is stable as long as both roots are negative, that is, their product is positive. The cell pattern becomes unstable when the product is negative, that is, when  (77)qx2≤2(3λ1δk2-ϵ2λ2),
for this requires |δk|≥(ϵ2λ2/3λ1), this condition defines the domain of Eckhaus instability. The above condition implies that the most unstable wave vector tends to zero, when |δk|→(ϵ2λ2/3λ1).


### 4.2.2. Transverse Perturbations and Zigzag Instability

Let us considerqx=0 into (75), we get(78)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0χ1y]+[2(ϵ2λ2-λ1(δk)2)+χ1y]χ1y=0,
where χ1y=λ1(qy2δk/qsc+qy4/4qsc2). The two eigenmodes are uncoupled and we have S(-),(79)S(-)=-2(ϵ2λ2-λ1(δk)2)-λ1qscδkqy2-λ14qsc2qy2<0,
for one of them. The other is amplified when  (80)S(+)=-λ1qy2(δk+qy24qsc)>0.
This implies that δk<0, the above condition defines the domain of the Zigzag Instability. When δk→0 from below the wave vector qy of the instability also tends to zero while the growth rate varies as qy2. We have studied the effect of magnetic field on long wavelength instabilities. We have observed that Eckhaus instability and Zigzag instability regions increases when Q increases (see Figure 4).Figure 4
Numerically computed secondary instability regions of Eckhaus instability (E), Zigzag instability (Z), and stable regions (S) are plotted in(λ2/λ1,δqs)-plane for Q=2000, Da=1500, Λ=0.85, ϕ=0.9, M=0.9, Pr1=1, and Pr2=2. As |δqs| increases then the secondary instability regions increases.


## 4.2.1. Longitudinal Perturbations and Eckhaus Instability

Insertingqy=0 into (75), we get(76)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0λ1qx2]+λ1qx2[2(ϵ2λ2-3λ1(δk)2)+qx2]=0,
since the roots are real and their sum always negative, the pattern is stable as long as both roots are negative, that is, their product is positive. The cell pattern becomes unstable when the product is negative, that is, when  (77)qx2≤2(3λ1δk2-ϵ2λ2),
for this requires |δk|≥(ϵ2λ2/3λ1), this condition defines the domain of Eckhaus instability. The above condition implies that the most unstable wave vector tends to zero, when |δk|→(ϵ2λ2/3λ1).


## 4.2.2. Transverse Perturbations and Zigzag Instability

Let us considerqx=0 into (75), we get(78)λ02S2+2S[2λ0(ϵ2λ2-λ1(δk)2)+λ0χ1y]+[2(ϵ2λ2-λ1(δk)2)+χ1y]χ1y=0,
where χ1y=λ1(qy2δk/qsc+qy4/4qsc2). The two eigenmodes are uncoupled and we have S(-),(79)S(-)=-2(ϵ2λ2-λ1(δk)2)-λ1qscδkqy2-λ14qsc2qy2<0,
for one of them. The other is amplified when  (80)S(+)=-λ1qy2(δk+qy24qsc)>0.
This implies that δk<0, the above condition defines the domain of the Zigzag Instability. When δk→0 from below the wave vector qy of the instability also tends to zero while the growth rate varies as qy2. We have studied the effect of magnetic field on long wavelength instabilities. We have observed that Eckhaus instability and Zigzag instability regions increases when Q increases (see Figure 4).Figure 4
Numerically computed secondary instability regions of Eckhaus instability (E), Zigzag instability (Z), and stable regions (S) are plotted in(λ2/λ1,δqs)-plane for Q=2000, Da=1500, Λ=0.85, ϕ=0.9, M=0.9, Pr1=1, and Pr2=2. As |δqs| increases then the secondary instability regions increases.


## 4.3. Heat Transport by Convection

The maximum of steady amplitudeA is denoted by |Amax| which is given as(81)|Amax⁡|=(ϵ2λ2λ3)1/2.
Equation (81) is obtained from (64) with tanh(X/Λ1)=1. We use |Amax| to calculate Nusselt number Nu. To discuss the heat transfer near the neutral region, we express it through the Nusselt number is defined as Nu=Hd/κΔT, which is the ratio of the heat transported across any layer to the heat which would be transported by conduction alone. Here H is the rate of heat transfer per unit area and is defined as(82)H=-〈∂Ttotal∂z′〉z′=0.
In (82), angular brackets correspond to a horizontal average. The Nusselt number Nu can be calculated in terms of amplitude A and is given as(83)Nu=1+ϵ2δsc2|Amax⁡|2.
From (83), we get conduction for R≤Rsc and convection for R>Rsc. Since the amplitude equation is valid for λ3>0, which is possible for R>Rsc (supercritical pitchfork bifurcation), thus we get Nu>1 for R>Rsc. We get convection for Nu>1 and conduction for Nu≤1. In stationary convection Nu increases implies that heat conducted by steady mode increases. In the problem of double diffusive convection in porous medium with magnetic field, Nu depends on Pr1,Pr2,Λ,M,ϕ,Da, and Q. We have computed Nu for different values of Q, for some fixed values of other parameters and observed that Nu increases as Q decreases (see Figures 5(a) and 5(b)). This implies that magnetic field inhibits the heat transport. The parameters Pr1, Pr2, Λ, M, ϕ, and Da show the same result as Q shows on Nu.Graph (a) is plotted forQ=1000 and graph (b) is plotted for Q=3000 for the fixed values of Da=1500, Λ=0.85, ϕ=0.9, Pr2=1, Pr1=2, and M=0.9. in (Nu,R/Rsc)-plane. In graphs (a) and (b), as R/Rsc increases then Nu increases.
(a)(b)


## 5. Oscillatory Convection at the Supercritical Hopf Bifurcation

The existence of a threshold (critical value of Rayleigh number for the onset of oscillatory convectionR=Roc) and a cellular structure (critical wave number q=qoc) are main characteristics of the oscillatory convection. In this section we treat the region near the onset of oscillatory convection. Here the axis of cylindrical rolls is taken as y-axis, so that y-dependence disappears from equation ℒw=𝒩. The z-dependence contained entirely in sine and cosine functions which ensure that the free-free boundary conditions are satisfied. The purpose of this section is to derive coupled one-dimensional nonlinear time-dependent Ginzburg-Landau type equations near the onset of oscillatory convection at supercritical Hopf bifurcation. We introduce ϵ as(84)ϵ2=Ro-RocRoc≪1.
We assume that(85)w0=[A1Lei(qocx+ωoct)+A1Rei(qocx-ωoct)+c.c.]sin⁡πz
is a solution to linearized equation ℒw=0, which satisfies free-free boundary conditions. Here A1L denotes the amplitude of left travelling wave of the roll and A1R denotes the amplitude of right travelling wave of the roll, which depends on slow space and time variables [23](86)X=ϵx,τ=ϵt,T=ϵ2t,
and assume that A1L=A1L(X,τ,T), A1R=A1R(X,τ,T). The differential operators can be expressed as(87)∂∂x⟶∂∂x+ϵ∂∂X,∂∂t⟶∂∂t+ϵ∂∂τ+ϵ2∂∂T.
The solution of basic equations can be sought as power series in ϵ,(88)f=ϵf0+ϵ2f1+ϵ3f2+⋯,
where f=f(u,v,w,θ,Hx,Hy,Hz) with the first approximation given by eigenvector of the linearized problem:(89)u0=iπqoc[A1Lei(qocx+ωoct)+A1Rei(qocx-ωoct)-c.c.]cos⁡πz,v0=0,Hy0=0,θ0=1M[1e1A1Lei(qocx+ωoct)+1e1*A1Rei(qocx-ωoct)+c.c.]×sin⁡πz,Hx0=-iπ2qoc[1e2A1Lei(qocx+ωoct)+1e2*A1Rei(qocx-ωoct)-c.c.]×sin⁡πz,Hz0=π[1e2A1Lei(qocx+ωoct)+1e2*A1Rei(qocx-ωoct)+c.c.]×cos⁡πz.
where δoc2=(π2+qoc2), e1=(δoc2+iωoc), and e2=(Mδoc2+iωocϕPr2/Pr1), here e1* and e2* are complex conjugate of e1 and e2.We expand the linear operatorℒ and nonlinear term 𝒩 as the following power series(90)L=L0+ϵL1+ϵ2L2+⋯,N=ϵ2N0+ϵ3N1+⋯,
substituting (87) and (88) into ℒw=𝒩, for each order of ϵ, we get(91)L0w0=0,(92)L0w1+L1w0=N0,(93)L0w2+L1w1+L2w0=N1.
Here(94)L0=(DϕDPr1-Q∂2∂z2)D∇2-RM∂2∂x2Dϕ,L1=∂∂τF1+2∂2∂xXF2,L2=∂F1∂T+∂4∂x2X2[Q∂2∂z2(M∇2-(Dϕ+MD))DPr1+ΛD∇2-ΛMDDϕ+ΛMDϕ∇2+Q∂2∂z2+R]+∂2F2∂X2+2∂2∂xX∂∂τ×[-ϕRPr2MPr1ϕRPr2MPr1ϕRPr2MPr1DϕDPr1-(M+ϕPr2Pr1)DPr1∇2-(ϕPr2Pr1ΛM+1MϕPr1)D∇2+1M2ϕPr1DDϕ+ϕPr2Pr1DDPr1-(ΛM+1M2ϕPr1)Dϕ∇2-Q∂2∂z2-ϕRPr2MPr1]+∂2∂τ2[ϕPr2Pr1DPr1∇2+Pr2M2Pr12D∇2+1M2ϕPr1Dϕ∇2-RM∇h2Dϕ],
where(95)F1=(DϕDPr1+ϕPr2Pr1DDPr1+1M2ϕPr1DDϕ)∇2-Q∇2∂2∂z2-ϕRPr2MPr1∇h2,F2=(DDϕ-Dϕ∇2-MD∇2)DPr1-ΛMDDϕ∇2+Q∇2∂2∂z2-QD∂2∂z2-RMDϕ+R∇h2.
Equation (91) is linear problem. We get critical Rayleigh number for the onset of oscillatory convection by using the zeroth-order solution w0 in (91). At O(ϵ2),𝒩0=0 and ℒ1w0=0 gives(96)∂A1L∂τ-vg∂A1L∂X=0,∂A1R∂τ-vg∂A1R∂X=0,
where vg=(∂ω/∂q)q=qoc is the group velocity and is real. Hence from (92), we get w1=0. From equation of continuity we find that u1=0. Substituting the zeroth-order and first-order approximation into (56) and (58) we get,(97)θ1=-πM2[(|A1L|2+|A1R|2)t1+2e1e4J1+2e1*e4*J1*]sin⁡2πz,v1=0,Hy1=0,Hx1=iπ2MqocPr2Pr1[(1e2-1e2*)(|A1L|2-|A1R|2)]sin⁡2πz,Hz1=2π2Pr2Pr1[1e2e5A1L2e2i(qocx+ωoct)+1e2*e5*A1R2e2i(qocx-ωoct)+14Mqoc2(1e2+1e2*)A1LA1Re2iqocx+c.c.],
where t1=(1/4π2)(1/e1+1/e1*), J1=A1LA1R*e2iωoct, e4=(4π2+2iωoc), and e5=(4Mqoc2+2iϕωocPr2/Pr1) and e4*, e5* and J1* are complex conjugate of e4,e5 and J1, respectively.Equation (93) is solvable when ℒ0w0=0, one requires that its right-hand side be orthogonal to w0, which is ensured that if the coefficients of sinπz in 𝒩1-ℒ2w0 are equal to zero. This implies that(98)Λ0∂A1L∂T+Λ1(∂∂τ-vg∂∂X)A2L-Λ2∂2A1L∂X2-Λ3A1L+Λ4|A1L|2A1L+Λ5|A1R|2A1L=0,Λ0∂A1R∂T+Λ1(∂∂τ-vg∂∂X)A2R-Λ2∂2A1R∂X2-Λ3A1R+Λ4|A1R|2A1R+Λ5|A1L|2A1R=0,
where(99)Λ0=(1M2ϕPr1e1e2+e2e3+ϕPr2Pr1e1e3+Qπ2)δoc2-Rocqoc2ϕPr2MPr1,Λ1=δoc2[e3ϕPr2Pr1+e1Pr2M2Pr12+e2M2ϕPr1],Λ2=4qoc2[e2e3+Me3δoc2+Λe1δoc2+ΛMe1e2+Me1e3+ΛMe2δoc2+Qπ2-R],Λ3=RMqoc2e2,Λ4=-Qπ4MPr22Pr12e1(1e2-1e2*)-2Qπ4Pr22Pr12e1e5(1δoc2+π2-3qoc2e2*)+RM3π2qoc2t1e2,Λ5=Qπ4MPr22Pr12e1(1e2-1e2*)-Qπ4Mqoc2Pr22Pr12e1(π2-qoc2)(1e2+1e2*)+RM3π2qoc2e2(t1+2e1e4).
Here e3=(iω/M2ϕPr+(Λ/M)δ2+1/MDa). It should be noted that A1L and A1R are of order ϵ and A2L and A2R are of order ϵ2. If ωoc=0 in Λ0,Λ2,Λ3, and Λ4 then these expressions match with the coefficients λ0, λ1, λ2, and λ3 of Ginzburg-Landau equation at the onset of stationary convection. From (96), we get A1L(ξ′,T) and A1R(η′,T), where ξ′=vgτ+X,η′=vgτ-X. Equations (98) can be written as(100)2vgΛ1∂A2L∂η′=-Λ0∂A1L∂T+Λ2∂A1L∂X2+λ3A1L-(Λ4|A1L|2+Λ5|A1R|2)A1L,(101)2vgΛ1∂A2R∂η′=-Λ0∂A1R∂T+Λ2∂A1R∂X2+λ3A1R-(Λ4|A1R|2+Λ5|A1L|2)A1R.Letξ′ϵ[0,l1], η′ϵ[0,l2] where l1,l2 are periods of A1L,A1R, respectively. Expansion (88) remains asymptotic for times t=O(ϵ-2) only if an appropriate solvability condition holds. This condition obtained integrating (100) over η′ and (101) over ξ′, we get(102)Λ0∂A1L∂T=Λ2∂A1L∂X2+λ3A1L-(Λ4|A1L|2+Λ5|A1R|2)A1L,(103)Λ0∂A1R∂T=Λ2∂A1R∂X2+λ3A1R-(Λ4|A1R|2+Λ5|A1L|2)A1R.


### 5.1. Travelling Wave and Standing Wave Convection

To study the stability regions of travelling waves and standing waves, Coullet et al. [24]. we proceed as follows.On dropping slow variableX from (102) and (103), we get a pair of first ODE's(104)dA1LdT=Λ3Λ0A1L-Λ4Λ0A1L|A1L|2-Λ5Λ0A1L|A1R|2,(105)dA1RdT=Λ3Λ0A1R-Λ4Λ0A1R|A1R|2-Λ5Λ0A1R|A1L|2.
Put(106)β′=Λ3Λ0,γ′=-Λ4Λ0,δ′=-Λ5Λ0.
Then (104) and (105) take the following form(107)dA1LdT=β′A1L+γ′A1L|A1L|2+δ′A1L|A1R|2,(108)dA1RdT=β′A1R+γ′A1R|A1R|2+δ′A1R|A1L|2.
Consider A1L=aLeiϕL and A1R=aLeiϕR (we can write a complex number in the amplitude and phase form), where aL=|A1L|, ϕL=arg(A1L)=tan-1(Im(A1L)/Re(A1L)) and aR=|A1R|, ϕR=arg(A1R)=tan-1(Im(A1R)/Re(A1R)), here aL, aR, ϕL, and ϕR are functions of time T since A1L and A1R are functions of T. Thus aL and aR are positive functions. Substituting the definitions of A1L and A1R and β′=β1+iβ2, γ′=γ1+iγ2, δ′=δ1+iδ2 into (107) and (108) we get,(109)daLdT=β1aL+γ1aL|aL|2+δ1aL|aR|2,(110)dϕLdT=β2+γ2|aL|2+δ2|aR|2,(111)daRdT=β1aR+γ1aR|aR|2+δ1aR|aL|2,(112)dϕRdT=β2+γ2|aR|2+δ2|aL|2.
Equations (109) and (111) not contain phase term, so we take these two equations for the future discussions. We have (109) and (111) as(113)daLdT=β1aL+γ1aL3+δ1aR2,daRdT=β1aR+γ1aR3+δ1aL2,
since aL and aR are positive functions. Put(114)daLdT=F1(aL,aR),daRdT=F2(aL,aR).
Now we discuss the stability of equilibrium points of (114). We get four equilibrium points like (aL,aR)=(0,0) (conduction state), (aL,aR)=(aL,0) (aL= amplitude of left travelling waves, here we get F2=0, and we get one condition from F1=0 i.e., aL2=-β1/γ1(=|A1L|2)), (aL,aR)=(0,aR) (aR= amplitude of right travelling waves, here we get F1=0, and we get one condition from F2=0 i.e., aR2=-β1/γ1(=|A1R|2)), and for aL≠0 and aR≠0 we get (aL,aR)=(-β1/(γ1+δ1),-β1/(γ1+δ1)) (this gives condition for standing waves. At standing waves we have AL=AR, so aL=aR). For the pair of (104) and (105), we do not get aL≠aR≠0 (modulated waves). Now the Jacobian of F1 and F2 is given by (115)(∂F1∂aL∂F1∂aR∂F2∂aL∂F2∂aR).If real parts of all eigenvalues of the Jacobian are negative at an equilibrium point, then that point is a stable equilibrium (Lyapunov's theorem or principle of linearized stability). Some valuable conditions for travelling waves and standing waves are travelling waves are stable ifβ1>0,γ1<0 and δ1<γ1<0. Standing waves are stable if β1>0,γ1<0 and (i) if δ1>0, then -γ1>δ1>0 and (ii) if δ1<0, then -γ1>-δ1>0.The stability regions of travelling waves and standing waves are summarized in Figure6. Here E is total amplitude and defined as E=aL2+aR2. We do not distinguish between left travelling waves and right travelling waves. For rest state (steady state) E=0, for travelling waves E=-β1/γ1, for standing waves E=-2β1/(γ1+ς1). Travelling waves are supercritical if γ1<0 and standing waves are supercritical if γ1+ς1<0. Figure 6(a) is drawn for stable travelling wave conditions and Figure 6(b) is drawn for stable standing wave conditions in (β1,E)-plane. The symbols (-,-) and (+,+) in Figures 6(a) and 6(b) indicate that both roots of Jacobian are negative and at least one root is positive between two roots. In Figures 6(a) and 6(b), travelling wave solution and standing wave solution bifurcate simultaneously from the steady-state solution (β1≥0 at this bifurcation point).(a, b, c) are typical diagrams showing the stability of equilibrium solutions SS (steady state), SW (standing waves), and TW (travelling waves). On solid lines equilibrium solutions are stable and on dotted lines they are unstable.
(a)(b)(c)In these Figures6(a) and 6(b), steady-state solution is stable for β1<0 and unstable β1>0. These figures show that for β1>0 both travelling waves and standing waves are supercritical. When travelling waves and standing waves bifurcate supercritically then at most one solution among travelling waves and standing waves will be stable. Thus, for β1>0 (Figure 6(a)) travelling waves are stable and (Figure 6(b)) standing waves are stable. In more detail we reproduce results of the stability analysis of equilibrium solutions in Figure 6(c), which is plotted in (γ1,ς1)-plane. From this figure we can observe that travelling waves are subcritical for γ1>0 and standing waves are subcritical for γ1+ς1>0. In Figure 7, We study the stability regions of travelling waves and standing waves at the onset of Hopf bifurcation. The stability regions of standing waves and travelling waves increases when Pr2/Pr1 increases for fixed parameters. For a fixed Pr1 if we get initially travelling waves at the onset of oscillatory convection then they are replaced by standing waves as Q increases.Figures (a–d), are plotted forPr2/Pr1=6,12,18, and 24, respectively, and for fixed parameters Da=1500, Λ=0.85, ϕ=0.9, and M=0.9. Stability regions of steady state (SS), travelling waves (TW), and standing waves (SW) are plotted (Q,Pr1)-plane.
(a)(b)(c)(d)


### 5.2. Long Wavelength Instabilities for the Onset of Travelling Wave Convection (Benjamin-Feir Instability)

For right travelling waveAR(X,T)=A(X,T) and AL(X,T)=0, for left travelling wave AR(X,T)=0 and AL(X,T)=A(X,T). Thus for travelling waves we get a single amplitude equation from (102) and (103), given as(116)Λ0∂A∂T-Λ2∂2A∂X2-Λ3A+Λ4|A|2A=0.
For standing waves A1L(X,T)=A1R(X,T)=A(X,T) and we get a single amplitude equation from (102) and (103), given as(117)Λ0∂A∂T-Λ2∂2∂X2A-Λ3A+(Λ4+Λ5)|A|2A=0.
Equation (117) possesses a family of planar wave solutions and solutions containing phase singular points, which describes weakly nonlinear wave phenomena [25]. We study the Benjamin-Feir instability of travelling waves from complex Ginzburg- Landau equation (116) can be written as(118)∂A∂T=ξ∂2A∂X2+βA+γ|A|2A,
where ξ=ξ1+iξ2,β=β1+iβ2,γ=γ1+iγ2. The phase winding solutions are obtained by substituting A=Ãoei(δqoX-δωT) into (118), and equating real and imaginary parts we get(119)|Ão|2=ξ1δqo2-β1γ1-1,δω=ξ2δqo2-β2+γ2(β1-ξ1δqo2)γ1-1.
Here Ão is constant and δqo=qX-qoc. We consider a modulated solution in the form: A(X,T)=Ã(X,T)ei(δqoX-δωT). Substituting the modulated into (118) which gives(120)∂Ã∂T=(γ1+iγ2)[(β1-δqo2ξ1)γ1+|Ã|2]Ã+(γ1+iγ2)(∂2∂X2+2iδqo∂∂X)Ã.
It is possible to conduct a general investigation of the linear stability of A(X,T), but this is very difficult task, and therefore our primary concern here is to treat the stability of the uniformly oscillating solution Ão. Inserting Ã=Ão+ũ+iṽ into (120) and equating real and imaginary parts we get(121)∂ũ∂T=-2(β1-δqo2ξ1)ũ+ξ1(∂2ũ∂X2-2δqo∂ṽ∂X)-ξ2(2δqo∂ũ∂X+∂ṽ∂X2),(122)∂ṽ∂T=-2γ2(β1-δqo2ξ1)γ1ũ+ξ1(2δqo∂ũ∂X+∂2ṽ∂X2)+ξ2(∂2ũ∂X2-2δqo∂ṽ∂X).
Consider (ũ,ṽ)=(U,V)eSTcosqXX and S in the growth rate of disturbances. Using solutions of ũ,ṽ, and δqo=0 into (121) and (122) we get,(123)(S+2β1+ξ1qX2)U-qX2ξ2V=0,(124)(S+qX2ξ1)V+(2β1γ2γ1-1+qX2ξ2)U=0,
solving (123) and (124), we get(125)S2+2S(β1+ξ1qX2)+qX2ξ1(2β1+ξ1qX2)+qX2ξ2(2β1γ2γ1-1+qX2ξ2).
There will be an instability only when a root of (125) is possible, that is,(126)2β1(ξ1+γ2ξ2γ1-1)+qX2(ξ12+ξ22)<0,β1>0 when travelling waves or standing waves are stable. The instability of waves against long wavelength longitudinal modes is often called the Benjamin-Feir instability. Thus we get Benjamin-Feir instability for travelling waves when ξ1+γ2ξ2/γ1<0. Similarly by considering (118) instead of (117) and proceeding in the same way we get Benjamin-Feir instability for standing waves when ξ1+(γ2+δ2)ξ2/(γ1+δ1)<0.


## 5.1. Travelling Wave and Standing Wave Convection

To study the stability regions of travelling waves and standing waves, Coullet et al. [24]. we proceed as follows.On dropping slow variableX from (102) and (103), we get a pair of first ODE's(104)dA1LdT=Λ3Λ0A1L-Λ4Λ0A1L|A1L|2-Λ5Λ0A1L|A1R|2,(105)dA1RdT=Λ3Λ0A1R-Λ4Λ0A1R|A1R|2-Λ5Λ0A1R|A1L|2.
Put(106)β′=Λ3Λ0,γ′=-Λ4Λ0,δ′=-Λ5Λ0.
Then (104) and (105) take the following form(107)dA1LdT=β′A1L+γ′A1L|A1L|2+δ′A1L|A1R|2,(108)dA1RdT=β′A1R+γ′A1R|A1R|2+δ′A1R|A1L|2.
Consider A1L=aLeiϕL and A1R=aLeiϕR (we can write a complex number in the amplitude and phase form), where aL=|A1L|, ϕL=arg(A1L)=tan-1(Im(A1L)/Re(A1L)) and aR=|A1R|, ϕR=arg(A1R)=tan-1(Im(A1R)/Re(A1R)), here aL, aR, ϕL, and ϕR are functions of time T since A1L and A1R are functions of T. Thus aL and aR are positive functions. Substituting the definitions of A1L and A1R and β′=β1+iβ2, γ′=γ1+iγ2, δ′=δ1+iδ2 into (107) and (108) we get,(109)daLdT=β1aL+γ1aL|aL|2+δ1aL|aR|2,(110)dϕLdT=β2+γ2|aL|2+δ2|aR|2,(111)daRdT=β1aR+γ1aR|aR|2+δ1aR|aL|2,(112)dϕRdT=β2+γ2|aR|2+δ2|aL|2.
Equations (109) and (111) not contain phase term, so we take these two equations for the future discussions. We have (109) and (111) as(113)daLdT=β1aL+γ1aL3+δ1aR2,daRdT=β1aR+γ1aR3+δ1aL2,
since aL and aR are positive functions. Put(114)daLdT=F1(aL,aR),daRdT=F2(aL,aR).
Now we discuss the stability of equilibrium points of (114). We get four equilibrium points like (aL,aR)=(0,0) (conduction state), (aL,aR)=(aL,0) (aL= amplitude of left travelling waves, here we get F2=0, and we get one condition from F1=0 i.e., aL2=-β1/γ1(=|A1L|2)), (aL,aR)=(0,aR) (aR= amplitude of right travelling waves, here we get F1=0, and we get one condition from F2=0 i.e., aR2=-β1/γ1(=|A1R|2)), and for aL≠0 and aR≠0 we get (aL,aR)=(-β1/(γ1+δ1),-β1/(γ1+δ1)) (this gives condition for standing waves. At standing waves we have AL=AR, so aL=aR). For the pair of (104) and (105), we do not get aL≠aR≠0 (modulated waves). Now the Jacobian of F1 and F2 is given by (115)(∂F1∂aL∂F1∂aR∂F2∂aL∂F2∂aR).If real parts of all eigenvalues of the Jacobian are negative at an equilibrium point, then that point is a stable equilibrium (Lyapunov's theorem or principle of linearized stability). Some valuable conditions for travelling waves and standing waves are travelling waves are stable ifβ1>0,γ1<0 and δ1<γ1<0. Standing waves are stable if β1>0,γ1<0 and (i) if δ1>0, then -γ1>δ1>0 and (ii) if δ1<0, then -γ1>-δ1>0.The stability regions of travelling waves and standing waves are summarized in Figure6. Here E is total amplitude and defined as E=aL2+aR2. We do not distinguish between left travelling waves and right travelling waves. For rest state (steady state) E=0, for travelling waves E=-β1/γ1, for standing waves E=-2β1/(γ1+ς1). Travelling waves are supercritical if γ1<0 and standing waves are supercritical if γ1+ς1<0. Figure 6(a) is drawn for stable travelling wave conditions and Figure 6(b) is drawn for stable standing wave conditions in (β1,E)-plane. The symbols (-,-) and (+,+) in Figures 6(a) and 6(b) indicate that both roots of Jacobian are negative and at least one root is positive between two roots. In Figures 6(a) and 6(b), travelling wave solution and standing wave solution bifurcate simultaneously from the steady-state solution (β1≥0 at this bifurcation point).(a, b, c) are typical diagrams showing the stability of equilibrium solutions SS (steady state), SW (standing waves), and TW (travelling waves). On solid lines equilibrium solutions are stable and on dotted lines they are unstable.
(a)(b)(c)In these Figures6(a) and 6(b), steady-state solution is stable for β1<0 and unstable β1>0. These figures show that for β1>0 both travelling waves and standing waves are supercritical. When travelling waves and standing waves bifurcate supercritically then at most one solution among travelling waves and standing waves will be stable. Thus, for β1>0 (Figure 6(a)) travelling waves are stable and (Figure 6(b)) standing waves are stable. In more detail we reproduce results of the stability analysis of equilibrium solutions in Figure 6(c), which is plotted in (γ1,ς1)-plane. From this figure we can observe that travelling waves are subcritical for γ1>0 and standing waves are subcritical for γ1+ς1>0. In Figure 7, We study the stability regions of travelling waves and standing waves at the onset of Hopf bifurcation. The stability regions of standing waves and travelling waves increases when Pr2/Pr1 increases for fixed parameters. For a fixed Pr1 if we get initially travelling waves at the onset of oscillatory convection then they are replaced by standing waves as Q increases.Figures (a–d), are plotted forPr2/Pr1=6,12,18, and 24, respectively, and for fixed parameters Da=1500, Λ=0.85, ϕ=0.9, and M=0.9. Stability regions of steady state (SS), travelling waves (TW), and standing waves (SW) are plotted (Q,Pr1)-plane.
(a)(b)(c)(d)


## 5.2. Long Wavelength Instabilities for the Onset of Travelling Wave Convection (Benjamin-Feir Instability)

For right travelling waveAR(X,T)=A(X,T) and AL(X,T)=0, for left travelling wave AR(X,T)=0 and AL(X,T)=A(X,T). Thus for travelling waves we get a single amplitude equation from (102) and (103), given as(116)Λ0∂A∂T-Λ2∂2A∂X2-Λ3A+Λ4|A|2A=0.
For standing waves A1L(X,T)=A1R(X,T)=A(X,T) and we get a single amplitude equation from (102) and (103), given as(117)Λ0∂A∂T-Λ2∂2∂X2A-Λ3A+(Λ4+Λ5)|A|2A=0.
Equation (117) possesses a family of planar wave solutions and solutions containing phase singular points, which describes weakly nonlinear wave phenomena [25]. We study the Benjamin-Feir instability of travelling waves from complex Ginzburg- Landau equation (116) can be written as(118)∂A∂T=ξ∂2A∂X2+βA+γ|A|2A,
where ξ=ξ1+iξ2,β=β1+iβ2,γ=γ1+iγ2. The phase winding solutions are obtained by substituting A=Ãoei(δqoX-δωT) into (118), and equating real and imaginary parts we get(119)|Ão|2=ξ1δqo2-β1γ1-1,δω=ξ2δqo2-β2+γ2(β1-ξ1δqo2)γ1-1.
Here Ão is constant and δqo=qX-qoc. We consider a modulated solution in the form: A(X,T)=Ã(X,T)ei(δqoX-δωT). Substituting the modulated into (118) which gives(120)∂Ã∂T=(γ1+iγ2)[(β1-δqo2ξ1)γ1+|Ã|2]Ã+(γ1+iγ2)(∂2∂X2+2iδqo∂∂X)Ã.
It is possible to conduct a general investigation of the linear stability of A(X,T), but this is very difficult task, and therefore our primary concern here is to treat the stability of the uniformly oscillating solution Ão. Inserting Ã=Ão+ũ+iṽ into (120) and equating real and imaginary parts we get(121)∂ũ∂T=-2(β1-δqo2ξ1)ũ+ξ1(∂2ũ∂X2-2δqo∂ṽ∂X)-ξ2(2δqo∂ũ∂X+∂ṽ∂X2),(122)∂ṽ∂T=-2γ2(β1-δqo2ξ1)γ1ũ+ξ1(2δqo∂ũ∂X+∂2ṽ∂X2)+ξ2(∂2ũ∂X2-2δqo∂ṽ∂X).
Consider (ũ,ṽ)=(U,V)eSTcosqXX and S in the growth rate of disturbances. Using solutions of ũ,ṽ, and δqo=0 into (121) and (122) we get,(123)(S+2β1+ξ1qX2)U-qX2ξ2V=0,(124)(S+qX2ξ1)V+(2β1γ2γ1-1+qX2ξ2)U=0,
solving (123) and (124), we get(125)S2+2S(β1+ξ1qX2)+qX2ξ1(2β1+ξ1qX2)+qX2ξ2(2β1γ2γ1-1+qX2ξ2).
There will be an instability only when a root of (125) is possible, that is,(126)2β1(ξ1+γ2ξ2γ1-1)+qX2(ξ12+ξ22)<0,β1>0 when travelling waves or standing waves are stable. The instability of waves against long wavelength longitudinal modes is often called the Benjamin-Feir instability. Thus we get Benjamin-Feir instability for travelling waves when ξ1+γ2ξ2/γ1<0. Similarly by considering (118) instead of (117) and proceeding in the same way we get Benjamin-Feir instability for standing waves when ξ1+(γ2+δ2)ξ2/(γ1+δ1)<0.


## 6. Conclusions

In this paper we have considered both linear and weakly nonlinear analysis of magnetoconvection in a sparsely packed porous medium in Earth's outer core by using free-free (stress-free) boundary conditions. Even though free-free boundary conditions cannot be achieved in laboratory, one can use it in geophysical fluid dynamic applications to Earth's outer core since they allow simple trigonometric eigenfunctions. Our goal is to identify the region of parameter values, for which rolls emerge at the onset of convection.Following Chandrasekhar [2], we have described the stationary convection and oscillatory convection as curves Rs(q) and Ro(q,Pr2) versus wave numbers. The critical wave numbers for stationary convection and oscillatory convection are qsc=qoc=π/2. For the problem of magnetoconvection in a sparsely packed porous medium, we get Takens-Bogdanov bifurcation point and codimension-two bifurcation point. In the case of linear theory both marginal and overstable motions are discussed. In Figures 1 and 2, is shown that the effect of Chandrasekhar number and porous parameter is to make the system more stable. By drawing stability boundaries in the Rayleigh number plane it is shown that the effect of magnetic field and porous parameter is to decrease the region of stabilities. In the nonlinear equation (60), λ0=0 gives the Takens-Bogdanov bifurcation point at qs=qsc and when λ0=0, (60) is not valid. The pitchfork bifurcation is supercritical if λ3>0 subcritical if λ3<0. and we get tricritical point if λ3=0. We have obtained from (60), long wave length instabilities, namely, Eckhaus and Zigzag instabilities. From (60) which is valid only for λ3>0, we have calculated Nusselt number Nu and studied heat transport by convection. We have also derived two one-dimensional nonlinear coupled Ginzburg-Landau type equations, namely, (98) at the onset of oscillatory convection at supercritical Hopf bifurcation. We have computed stability regions of SW and TW at both Hopf bifurcation. The conditions for SW and TW are AL=AR and AL=0 or AR=0, respectively. TW exist if |AL|2=-β1/γ1>0 and they are supercritical if γ1<0. SW exist if |AL|2=|AR|2=-β1/(γ1+δ1)>0 and SW are supercritical if γ1+δ1<0. When both SW and TW are supercritical then at most one equilibrium solution is stable. At Takens-Bogdanov bifurcation point we get both TW and SW. By deriving one-dimensional Ginzburg-Landau equations with complex coefficients, namely, (116) and (117), we have shown the existence of Benjamin-Feir-type of instability for both TW and SW. Near the Takens-Bogdanov bifurcation point the conducting state becomes unstable against both stationary and oscillatory mode, that is, the real parts of two eigenvalues pass through zero simultaneously. This violates the assumption made for deriving amplitude equations (60) and  (98). Instead a new equation, which is second order in time, has to be used near the Takens-Bogdanov bifurcation point.



---

*Source: 207123-2011-12-07.xml*