# Linear and Nonlinear Stability Analysis of Double Diffusive Convection in a Maxwell Fluid Saturated Porous Layer with Internal Heat Source

Authors: Moli Zhao; Qiangyong Zhang; Shaowei Wang
Journal: Journal of Applied Mathematics
 (2014)

---

## Abstract

The onset of double diffusive convection is investigated in a Maxwell fluid saturated porous layer with internal heat source. The modified Darcy law for the Maxwell fluid is used to model the momentum equation of the system, and the criterion for the onset of the convection is established through the linear and nonlinear stability analyses. The linear analysis is obtained using the normal mode technique, and the nonlinear analysis of the system is studied with the help of truncated representation of Fourier series. The effects of internal Rayleigh number, stress relaxation parameter, normalized porosity, Lewis number, Vadasz number and solute Rayleigh number on the stationary, and oscillatory and weak nonlinear convection of the system are shown numerically and graphically. The effects of various parameters on transient heat and mass transfer are also discussed and presented analytically and graphically.


---

## Body

## 1. Introduction

Double diffusive convection in porous media without heat source has been intensively studied because of its application in different branches of science and engineering, such as underground disposal of nuclear wastes, groundwater pollution, contaminant transport in fluid-saturated soils, liquid gas storage, and food processing [1–4]. Since convective flow plays an important role in ore-forming systems, Zhao et al. have conducted extensive and systematic studies on double diffusive convection within the upper crust of the Earth during the last decade or so [5–9]. Their work has promoted the better understanding of ore-forming mechanisms, which are essential to explore new large ore deposits in the deep Earth [8, 9]. And the onset of convection due to internal heat source has become an interesting problem in various areas of geophysics and engineering under the situations of radioactive decay or a weak exothermic reaction within the porous material. Therefore, diffusive convection in porous media with internal heat source has attracted the attention of many authors like Nield and Bejan [10], Ingham and Pop [11], and Vafai [12] during the last several decades.The onset of thermodynamic instability in horizontal porous layer saturated with Newtonian fluid was first studied extensively on geological and engineering length scales [13–23]. Convection by internal heat sources has been studied in several papers including experimental and theoretical studies. The onset of convection in a fluid saturated horizontal layer of an anisotropic porous medium with internal heat source subjected to inclined temperature gradient has been considered [24]. The effect of internal heat source on the problem of triple diffusive convection has been analyzed by Straughan and Tracey [25]. Magyari et al. have carried out an analytical and numerical study about the effect of boundary-layer flows in a fluid-saturated porous medium with internal heat generation [26, 27]. Hill has investigated linear and nonlinear stability analyses of double diffusive convection in a fluid saturated porous layer with internal concentration source [28]. Bhadauria et al. have studied the linear and nonlinear thermal instability in an anisotropic saturated porous layer with internal heat source [29]. Bhadauria has made the stability analysis of convection in a binary fluid-saturated horizontal porous layer with internal heat source [30].Recently, viscoelastic fluid flow in porous media has attracted considerable attention, due to the large demands of such diverse fields as biorheology, geophysics, chemical industries, and petroleum industries. Wang and Tan have made the stability analysis of double diffusive convection in a Maxwell fluid saturated porous medium [31]. It is worthwhile to point out that the first viscoelastic rate type model, which is still used widely, is due to Maxwell. The onset of double diffusive convection in a viscoelastic saturated porous layer has been considered by many researchers (e.g., [31–36]). However, we are unaware of the double diffusive problem in any previous work which also introduces penetrative convection with internal heat source.In this paper, we focus on the linear and weakly nonlinear stability analyses in a viscoelastic fluid saturated porous layer with internal heat source using the Darcy-Maxwell model. The Dufour and Soret effects are ignored. The aim of the present paper is to study how onset criteria for stationary and oscillatory double diffusive convection are affected by the viscoelastic parameter and other parameters, as well as discussing their effects on heat and mass transfer.


## 2. Mathematical Model

Assuming that an infinite horizontal porous layer saturated with Maxwell fluid mixture with internal heat source, confined between the planes,z
=
0 and z
=
d. A Cartesian frame of reference is chosen with the origin in the lower boundary and the z-axis vertically upwards. For the purposes of the present paper, we consider that temperature difference Δ
T and concentration difference Δ
S are maintained between the lower and upper boundaries. In the Boussinesq approximation, the equation of state is
(1)
ρ
=
ρ
0
[
1
-
α
T
(
T
-
T
0
)
+
α
S
(
S
-
S
0
)
]
,
where T is the temperature, S is the concentration, T
0 and S
0 are temperature and concentration at the above plate, α
T and α
S are thermal and solutal expansion coefficients in the medium, ρ is the density of fluid, and ρ
0 is the density at some reference temperature T
0 and concentration S
0.Considering the vertically downward gravity forceg acting on it and neglecting the off-diagonal (Soret, Dufour, and cross-diffusion) contributions to the fluxes of the stratifying agencies, by using of the modified Darcy-Maxwell model [35], the governing system for double diffusive of Maxwell fluid in a porous layer with internal heat source can be represented by
(2)
∇
·
q
=
0
,
(
1
+
λ
1
∂
∂
t
)
(
ρ
0
ε
∂
q
∂
t
+
∇
p
-
ρ
g
)
+
μ
K
q
=
0
,
γ
∂
T
∂
t
+
q
·
∇
T
=
κ
T
∇
2
T
+
Q
(
T
-
T
0
)
,
ε
∂
S
∂
t
+
q
·
∇
S
=
κ
S
∇
2
S
.
The thermal boundary conditions are
(3)
T
=
T
0
+
Δ
T
at
z
=
0
,
T
=
T
0
at
z
=
d
,
S
=
S
0
+
Δ
S
at
z
=
0
,
S
=
S
0
at
z
=
d
,
where q
=
(
u
,
v
,
w
) is the Darcy velocity, p is the pressure, λ
1 is the relaxation time, Q is internal heat source, g is the gravitational acceleration, μ is the viscosity, while K and ε are the permeability and porosity of the medium, and γ
=
(
ρ
c
)
m
/
(
ρ
c
)
f is the ratio of heat capacities, κ
T and κ
S are effective thermal and solutal diffusivity of the medium, respectively.The basic state of the fluid is assumed to be quiescent, and the quantities of the basic state are given by(4)
q
b
=
0
,
p
=
p
b
(
z
)
,
T
=
T
b
(
z
)
,
S
=
S
b
(
z
)
,
which satisfy the following conditions:
(5)
d
p
b
d
z
=
-
ρ
b
g
,
d
2
S
b
d
z
2
=
0
,
κ
T
d
2
(
T
b
-
T
0
)
d
z
2
+
Q
(
T
b
-
T
0
)
=
0
.
Here the subscript b refers to the basic state. Then the steady state solutions are given by
(6)
T
b
=
T
0
+
Δ
T
sin
⁡
R
i
(
1
-
z
/
d
)
sin
⁡
R
i
,
S
b
=
S
0
+
Δ
S
(
1
-
z
d
)
,
p
b
=
p
0
-
ρ
0
g
z
×
(
1
-
α
T
Δ
T
sin
⁡
R
i
(
1
-
z
/
d
)
sin
⁡
R
i
+
α
S
Δ
S
(
1
-
z
d
)
)
,
where R
i
=
Q
d
2
/
κ
T is internal Rayleigh number. On the basic state, we superimpose perturbations in the following form:
(7)
q
=
q
b
+
q
′
,
p
=
p
b
+
p
′
,
T
=
T
b
+
T
′
,
S
=
S
b
+
S
′
,
ρ
=
ρ
b
+
ρ
′
,
where primes indicate perturbed quantities. Substituting (7) into (1)-(2) and using the basic state solutions, we obtain the following governing equations:
(8)
∇
·
q
′
=
0
,
(9)
(
1
+
λ
1
∂
∂
t
)
(
ρ
0
ε
∂
q
′
∂
t
+
∇
p
′
+
ρ
0
g
(
α
T
T
'
-
α
S
S
'
)
)
+
μ
K
q
′
=
0
,
(10)
γ
∂
T
′
∂
t
+
q
′
·
∇
T
'
+
w
′
∂
T
b
∂
z
=
κ
T
∇
2
T
'
+
Q
T
'
,
(11)
ε
∂
S
'
∂
t
+
q
′
·
∇
S
'
-
w
′
Δ
S
'
d
=
κ
S
∇
2
S
'
.
Eliminating the pressure term from the momentum equation by using curl-curl operator and introducing the following nondimensional transformation (after ignoring the primes for simplicity):
(12)
x
=
x
*
d
,
t
=
γ
d
2
κ
T
t
*
,
q
=
κ
T
d
q
*
,
p
=
μ
κ
T
K
p
*
,
T
=
(
Δ
T
)
T
*
,
S
=
(
Δ
S
)
S
*
,
we obtain the nondimensional governing equations (on dropping the asterisks for simplicity):
(13)
(
1
+
λ
∂
∂
t
)
(
η
Va
∂
∂
t
∇
2
w
-
Ra
∇
1
2
(
T
-
N
S
)
)
+
∇
2
w
=
0
,
(
∂
∂
t
-
∇
2
-
R
i
)
T
+
q
·
∇
T
+
w
f
(
z
)
=
0
,
(
η
∂
∂
t
-
1
Le
∇
2
)
S
+
q
·
∇
S
-
w
=
0
,
where λ
=
λ
1
κ
/
γ
d
2 the relaxation number, η
=
ε
/
γ the normalized porosity, Da
=
K
/
d
2 the Darcy number, Pr
=
μ
/
ρ
0
κ
T the Prandtl number, Va
=
ε
2
Pr
/
Da the Vadasz number, Ra
=
ρ
0
g
α
T
K
d
Δ
T
/
μ
κ
T the thermal Rayleigh number, N
=
α
S
Δ
S
/
α
T
Δ
T the buoyancy ratio, Le
=
κ
T
/
κ
S the Lewis number, and f
(
z
)
=
∂
T
b
/
∂
z, where T
b in nondimensionalized form is given by
(14)
T
b
=
sin
⁡
R
i
(
1
-
z
)
sin
⁡
R
i
.
The boundaries of the system considered here are impermeable isothermal and isosolutal. Hence, the boundary conditions for the perturbation variables are given by
(15)
w
=
T
=
S
=
0
at
z
=
0,1
.


## 3. Linear Stability Analysis

In this section, we discuss the linear stability analysis. According to the normal mode analysis [8, 34], convective motion is assumed to exhibit horizontal periodicity. Then the perturbed quantities can be assumed to be periodic waves of the form
(16)
(
W
T
S
)
=
(
W
(
z
)
Θ
(
z
)
Φ
(
z
)
)
exp
⁡
⁡
[
i
(
l
x
+
m
y
)
+
ω
t
]
,
where l and m are the wavenumbers in the horizontal plane and ω is the growth rate. Infinitesimal perturbations of the rest state may either damp or grow depending on the value of the parameter ω. Substituting (16) into (13) yields
(17)
(
1
+
λ
ω
)
(
η
ω
Va
(
D
2
-
a
2
)
W
+
a
2
Ra
(
Θ
-
N
Φ
)
)
+
(
D
2
-
a
2
)
W
=
0
,
W
f
(
z
)
+
(
-
D
2
+
a
2
+
ω
-
R
i
)
Θ
=
0
,
W
+
(
1
Le
(
D
2
-
a
2
)
-
η
ω
)
Φ
=
0
,
where D
=
d
/
d
z and a
2
=
l
2
+
m
2. Now, the boundary conditions become
(18)
W
=
Θ
=
Φ
=
0
at
z
=
0,1
.
To satisfy the boundary conditions (18), we assume the solutions of (17) in the following form:
(19)
(
W
(
z
)
Θ
(
z
)
Φ
(
z
)
)
=
(
W
0
Θ
0
Φ
0
)
sin
⁡
n
π
z
,
(
n
=
1,2
,
3
,
…
)
,
which is the most unstable mode when n
=
1, that is, the fundamental mode. Substituting the above equation into (17) yields
(20)
(
1
+
λ
ω
)
(
η
ω
Va
δ
2
W
0
-
a
2
Ra
(
Θ
0
-
N
Φ
0
)
)
+
δ
2
W
0
=
0
,
2
F
W
0
+
(
δ
2
+
ω
-
R
i
)
Θ
0
=
0
,
W
0
-
(
δ
2
Le
+
η
ω
)
Φ
0
=
0
,
where δ
2
=
π
2
+
a
2 is the total wave number and F
=
∫
0
1
f
(
z
)
sin
2
π
z
d
z. Rewrite the above equations in the matrix form as
(21)
(
δ
2
(
1
+
Λ
η
ω
Va
)
-
a
2
Ra
Λ
a
2
Ra
N
Λ
2
F
(
z
)
(
δ
2
+
ω
-
R
i
)
0
1
1
-
(
δ
2
Le
+
η
ω
)
)
(
W
0
Θ
0
Φ
0
)
=
(
0
0
0
)
,
where Λ
=
1
+
λ
ω. Considering the nontrivial solution of the above matrix equation, we require
(22)
Ra
=
δ
2
(
1
+
Λ
(
η
ω
/
Va
)
)
(
R
i
-
δ
2
-
ω
)
2
a
2
Λ
F
+
(
R
i
-
δ
2
-
ω
)
Ra
S
2
F
(
δ
2
+
Le
η
ω
)
,
where Ra
S
=
Ra
N
Le
=
ρ
0
g
α
S
K
d
Δ
S
/
μ
κ
S is the solutal Rayleigh number. For neutral solutions, we set ω
=
i
ω
i in (22) and rearrange the terms to get the Rayleigh number (dropping the subscript i of ω
i):
(23)
Ra
T
=
Δ
1
+
i
ω
Δ
2
,
where Δ
1 and Δ
2 are real part and imaginary part given by
(24)
Δ
1
=
δ
2
(
R
i
-
δ
2
-
λ
ω
2
+
(
η
/
Va
)
ω
2
+
λ
2
η
ω
4
/
Va
)
2
a
2
F
(
z
)
(
1
+
λ
2
ω
2
)
+
Ra
S
(
R
i
δ
2
-
δ
4
-
Le
η
ω
2
)
2
F
(
z
)
(
δ
4
+
Le
2
η
2
ω
2
)
,
Δ
2
=
(
+
λ
2
η
R
i
ω
2
-
λ
2
δ
2
η
ω
2
)
)
δ
2
(
R
i
η
-
η
δ
2
-
Va
(
λ
R
i
-
λ
δ
2
+
1
)
+
λ
2
η
R
i
ω
2
-
λ
2
δ
2
η
ω
2
)
)
×
(
2
a
2
F
(
z
)
(
1
+
λ
2
ω
2
)
Va
)
-
1
+
Ra
S
(
-
Le
η
R
i
+
Le
η
δ
2
-
δ
2
)
2
F
(
z
)
(
δ
4
+
Le
2
η
2
ω
2
)
.
Since Ra is a physical quantity, it must be real. Hence, from (24) it follows that either ω
=
0 (steady onset) or Δ
2
=
0 (ω
≠
0, oscillatory onset).


### 3.1. Stationary Convection

For the validity of principle of exchange of stabilities (i.e., steady case), we haveω
=
0 at the margin of stability. Then, the Rayleigh number at which marginally stable steady mode becomes
(25)
Ra
T
s
t
=
-
1
2
F
(
z
)
(
δ
2
(
δ
2
-
R
i
)
a
2
+
(
δ
2
-
R
i
)
Ra
S
δ
2
)
.
In the absence of the internal heat source; that is, Q
=
0, we have R
i
=
0 and F
(
z
)
=
-
1
/
2. So we obtain
(26)
Ra
T
s
t
=
(
π
2
+
a
2
)
2
a
2
+
Ra
S
,
which coincides with the known results [8, 10]. Furthermore, when a
S
=
0, the stationary Rayleigh number reduces to the classical result:
(27)
Ra
T
s
t
=
(
π
2
+
a
2
)
2
a
2
.
In addition, (26) gives the critical value Ra
T
,
c
s
t
=
4
π
2 for a
c
s
t
=
π.


### 3.2. Oscillatory Convection

For oscillatory onset,Δ
2
=
0 and ω
≠
0, and then
(28)
Ra
T
Osc
=
δ
2
(
R
i
-
δ
2
-
λ
ω
2
+
(
η
/
Va
)
ω
2
+
λ
2
η
ω
4
/
Va
)
2
a
2
F
(
z
)
(
1
+
λ
2
ω
2
)
+
Ra
S
(
R
i
δ
2
-
δ
4
-
Le
η
ω
2
)
2
F
(
z
)
(
δ
4
+
Le
2
η
2
ω
2
)
,
where the frequency ω in (28) is given by
(29)
b
2
(
ω
2
)
2
+
b
1
ω
2
+
b
0
=
0
.
Here
(30)
b
2
=
(
λ
2
η
R
i
-
λ
2
δ
2
η
)
Le
2
η
2
δ
2
,
b
1
=
(
R
i
η
-
η
δ
2
-
Va
(
λ
R
i
-
λ
δ
2
+
1
)
)
×
Le
2
η
2
δ
2
+
(
λ
2
η
R
i
-
λ
2
δ
2
η
)
δ
6
+
a
2
λ
2
Va
Ra
S
(
-
Le
η
R
i
+
Le
η
δ
2
-
δ
2
)
,
b
0
=
(
R
i
η
-
η
δ
2
-
Va
(
λ
R
i
-
λ
δ
2
+
1
)
)
δ
6
+
a
2
Va
Ra
S
(
-
Le
η
R
i
+
Le
η
δ
2
-
δ
2
)
.
Now, we try to find the positive solutions of (29). If there is none, then no oscillatory convection instability occurs. If there are two positive roots, then the minimum of (28) with ω
2 obtained by (29) gives the oscillatory neutral Rayleigh number. If we find only one positive solution, then substituting ω
2 into (27) yields the oscillatory Rayleigh number. After that, the effects of different parameters on the onset of oscillatory convection can be discussed from the Rayleigh number.


## 3.1. Stationary Convection

For the validity of principle of exchange of stabilities (i.e., steady case), we haveω
=
0 at the margin of stability. Then, the Rayleigh number at which marginally stable steady mode becomes
(25)
Ra
T
s
t
=
-
1
2
F
(
z
)
(
δ
2
(
δ
2
-
R
i
)
a
2
+
(
δ
2
-
R
i
)
Ra
S
δ
2
)
.
In the absence of the internal heat source; that is, Q
=
0, we have R
i
=
0 and F
(
z
)
=
-
1
/
2. So we obtain
(26)
Ra
T
s
t
=
(
π
2
+
a
2
)
2
a
2
+
Ra
S
,
which coincides with the known results [8, 10]. Furthermore, when a
S
=
0, the stationary Rayleigh number reduces to the classical result:
(27)
Ra
T
s
t
=
(
π
2
+
a
2
)
2
a
2
.
In addition, (26) gives the critical value Ra
T
,
c
s
t
=
4
π
2 for a
c
s
t
=
π.


## 3.2. Oscillatory Convection

For oscillatory onset,Δ
2
=
0 and ω
≠
0, and then
(28)
Ra
T
Osc
=
δ
2
(
R
i
-
δ
2
-
λ
ω
2
+
(
η
/
Va
)
ω
2
+
λ
2
η
ω
4
/
Va
)
2
a
2
F
(
z
)
(
1
+
λ
2
ω
2
)
+
Ra
S
(
R
i
δ
2
-
δ
4
-
Le
η
ω
2
)
2
F
(
z
)
(
δ
4
+
Le
2
η
2
ω
2
)
,
where the frequency ω in (28) is given by
(29)
b
2
(
ω
2
)
2
+
b
1
ω
2
+
b
0
=
0
.
Here
(30)
b
2
=
(
λ
2
η
R
i
-
λ
2
δ
2
η
)
Le
2
η
2
δ
2
,
b
1
=
(
R
i
η
-
η
δ
2
-
Va
(
λ
R
i
-
λ
δ
2
+
1
)
)
×
Le
2
η
2
δ
2
+
(
λ
2
η
R
i
-
λ
2
δ
2
η
)
δ
6
+
a
2
λ
2
Va
Ra
S
(
-
Le
η
R
i
+
Le
η
δ
2
-
δ
2
)
,
b
0
=
(
R
i
η
-
η
δ
2
-
Va
(
λ
R
i
-
λ
δ
2
+
1
)
)
δ
6
+
a
2
Va
Ra
S
(
-
Le
η
R
i
+
Le
η
δ
2
-
δ
2
)
.
Now, we try to find the positive solutions of (29). If there is none, then no oscillatory convection instability occurs. If there are two positive roots, then the minimum of (28) with ω
2 obtained by (29) gives the oscillatory neutral Rayleigh number. If we find only one positive solution, then substituting ω
2 into (27) yields the oscillatory Rayleigh number. After that, the effects of different parameters on the onset of oscillatory convection can be discussed from the Rayleigh number.


## 4. Nonlinear Analysis

In this section, we study the nonlinear stability analysis using minimal truncated Fourier series. For simplicity, we confine ourselves to two-dimensional rolls, so that all the physical quantities are independent ofy. Defining the stream function ψ such that u
=
∂
ψ
/
∂
z, w
=
-
∂
ψ
/
∂
x and substituting them into (9)–(11) and then eliminating the pressure term, nondimensionalizing the governing equations yields
(31)
(
1
+
λ
∂
∂
t
)
(
η
Va
∂
∂
t
(
∇
2
ψ
)
+
Ra
∂
T
∂
x
-
Ra
S
Le
∂
S
∂
x
)
+
∇
2
ψ
=
0
,
(
∂
∂
t
-
∇
2
-
R
i
)
T
=
∂
ψ
∂
x
f
(
z
)
+
∂
(
ψ
,
T
)
∂
(
x
,
z
)
,
η
∂
S
∂
t
-
∂
(
ψ
,
S
)
∂
(
x
,
z
)
=
1
Le
∇
2
S
-
∂
ψ
∂
x
.
We assume that the basic circulation remains undistorted but the temperature and concentration fields are distorted by the addition of a second harmonic with no x-dependence, and thus
(32)
ψ
=
A
1
(
t
)
sin
⁡
(
a
x
)
sin
⁡
(
π
z
)
,
(33)
T
=
B
1
(
t
)
cos
⁡
⁡
⁡
(
a
x
)
sin
⁡
⁡
(
π
z
)
+
B
2
(
t
)
sin
⁡
(
2
π
z
)
,
(34)
S
=
C
1
(
t
)
cos
⁡
⁡
(
a
x
)
sin
⁡
(
π
z
)
+
C
2
(
t
)
sin
⁡
(
2
π
z
)
,
where A
1, B
1, B
2, C
1, and C
2 are time dependent amplitudes and to be determined from the dynamics of the system. Substituting (32)–(34) into the coupled nonlinear system of partial differential equation (31) and equating coefficients for expression terms, we obtain the following nonlinear system:
(35)
d
X
d
t
=
D
,
where X
=
(
A
1
,
A
2
,
B
1
,
B
2
,
C
1
,
C
2
)
T, D
=
(
D
1
,
D
2
,
D
3
,
D
4
,
D
5
,
D
6
)
T with
(36)
D
1
=
A
2
,
D
2
=
-
Va
λ
η
δ
2
(
η
δ
2
Va
A
2
+
δ
2
A
1
+
a
Ra
B
1
-
a
Ra
S
Le
C
1
+
a
Ra
λ
D
3
-
a
λ
Ra
S
Le
D
5
(
η
δ
2
V
a
A
2
+
δ
2
A
1
+
a
R
a
B
1
-
a
R
a
S
L
e
C
1
)
,
D
3
=
2
a
F
(
z
)
A
1
-
π
a
A
1
B
2
+
(
R
i
-
δ
2
)
B
1
,
D
4
=
π
a
2
A
1
B
1
+
(
R
i
-
4
π
2
)
B
2
,
D
5
=
-
π
a
η
A
1
C
2
-
δ
2
η
Le
C
1
-
a
η
A
1
,
D
6
=
π
a
2
η
A
1
C
1
-
4
π
2
Le
η
C
2
.
The above ordinary differential equations are autonomous, which is difficult to obtain the solutions by use of analytical technique, so we solve the above system numerically. We note that the above system is uniformly bounded in time and possesses many properties to the full problem. Equation (34) must be dissipative. Therefore, the divergence of the flow is always negative. So we have
(37)
∂
A
˙
1
∂
A
1
+
∂
A
˙
2
∂
A
2
+
∂
B
˙
1
∂
B
1
+
∂
B
˙
2
∂
B
2
+
∂
C
˙
1
∂
C
1
+
∂
C
˙
2
∂
C
2
=
-
(
1
λ
-
2
R
i
+
(
4
π
2
+
δ
2
)
(
1
+
1
Le
η
)
)
.
The system is bounded and dissipative. Consequently, the trajectories are attracted to a set of measure zero in the phase space. In particular, they may be attracted to a fixed point, a limit cycle or, perhaps, a strange attractor. From (37), it can also be concluded that if a set of initial points in phase space occupy a region V
(
0
) at time t
=
0, then after some time t we have
(38)
V
(
t
)
=
V
(
0
)
exp
⁡
{
-
[
1
λ
-
2
R
i
+
(
4
π
2
+
δ
2
)
(
1
+
1
Le
η
)
]
t
}
,
which indicates that the volume decreases exponentially with time. Moreover, we can conclude that an increase in the value of relaxation parameter λ, Lewis number Le, and normalized porosity η is to enhance the dissipation. Furthermore, we observe that (37) has an important symmetry, for it is invariant under the transformation:
(39)
(
A
1
,
A
2
,
B
1
,
B
2
,
C
1
,
C
2
)
⟶
(
-
A
1
,
-
A
2
,
-
B
1
,
B
2
,
-
C
1
,
C
2
)
.


### 4.1. Steady Finite Amplitude Motions

The simplified model represented by (37) has the great advantage that steady finite amplitude solutions can be obtained at once and their stability can be investigated analytically. From qualitative prediction, we look into the possibility of an analytical solution. In the case of steady motions, setting the left hand side of (37) equal to zero yields(40a)
δ
2
A
1
+
a
Ra
B
1
-
a
Ra
s
Le
C
1
=
0
,
(40b)
2
a
F
(
z
)
A
1
-
π
a
A
1
B
2
+
(
R
i
-
δ
2
)
B
1
=
0
,
(40c)
π
a
2
A
1
B
1
+
(
R
i
-
4
π
2
)
B
2
=
0
,
(40d)
-
π
a
η
A
1
C
2
-
δ
2
η
Le
C
1
-
a
η
A
1
=
0
,
(40e)
π
a
2
η
A
1
C
1
-
4
π
2
Le
η
C
2
=
0
.On solving for the amplitudes in terms of A
1, we obtain
(41)
B
1
=
2
a
F
(
z
)
(
8
π
2
-
2
R
i
)
A
1
π
2
a
2
A
1
2
-
(
R
i
-
δ
2
)
(
8
π
2
-
2
R
i
)
,
C
1
=
-
8
a
Le
A
1
8
δ
2
+
a
2
Le
2
A
1
2
,
B
2
=
2
π
a
2
A
1
2
F
(
z
)
π
2
a
2
A
1
2
-
(
R
i
-
δ
2
)
(
8
π
2
-
2
R
i
)
,
C
2
=
-
a
2
Le
2
A
1
2
π
(
8
δ
2
+
a
2
Le
2
A
1
2
)
.
Substituting B
1 and C
1 into (39) yields
(42)
a
2
x
2
+
a
1
x
+
a
0
=
0
,
where x
=
A
1
2
/
8, and
(43)
a
2
=
δ
2
M
2
M
5
,
a
1
=
δ
4
M
2
-
δ
2
M
3
M
5
+
M
1
M
5
+
M
2
M
4
,
a
0
=
M
1
δ
2
-
M
3
M
4
-
M
3
δ
4
,
M
1
=
2
a
2
Ra
F
(
z
)
(
π
2
-
R
i
4
)
,
M
2
=
π
2
a
2
,
M
3
=
(
R
i
-
δ
2
)
(
π
2
-
R
i
4
)
,
M
4
=
a
2
Ra
S
,
M
5
=
a
2
Le
2
.
The required root of (42) is
(44)
A
1
2
8
=
-
a
1
+
a
1
2
-
4
a
2
a
0
2
a
2
.
If the discriminant equals zero, we obtain the expression for finite amplitude Rayleigh number Ra
T
F, which characterizes the onset of finite amplitude steady motions. The finite amplitude Rayleigh number can be expressed in the following form:
(45)
Ra
T
F
=
-
e
1
+
e
1
2
-
4
e
2
e
0
2
e
2
,
where
(46)
e
2
=
a
8
Le
4
F
2
(
z
)
(
π
2
-
R
i
4
)
2
-
4
δ
4
a
6
F
2
(
z
)
Le
2
(
π
2
-
R
i
4
)
2
,
e
1
=
2
a
4
Le
2
F
(
z
)
(
(
π
2
-
R
i
4
)
+
π
2
a
4
R
a
S
)
δ
4
π
2
a
2
-
δ
2
a
2
Le
2
(
R
i
-
δ
2
)
×
(
π
2
-
R
i
4
)
+
π
2
a
4
Ra
S
)
(
π
2
-
R
i
4
)
+
4
δ
2
a
4
F
(
z
)
Le
2
(
π
2
-
R
i
4
)
×
a
2
Ra
S
(
R
i
-
δ
2
)
(
π
2
-
R
i
4
)
+
δ
4
(
R
i
-
δ
2
)
(
π
2
-
R
i
4
)
,
e
0
=
(
(
π
2
-
R
i
4
)
+
π
2
a
4
R
a
S
(
δ
4
π
2
a
2
-
δ
2
a
2
L
e
2
(
R
i
-
δ
2
)
)
δ
4
π
2
a
2
-
δ
2
a
2
Le
2
(
R
i
-
δ
2
)
×
(
π
2
-
R
i
4
)
+
π
2
a
4
Ra
S
(
δ
4
π
2
a
2
-
δ
2
a
2
L
e
2
(
R
i
-
δ
2
)
)
2
.


### 4.2. Heat Transfer and Mass Transport

In the study of convection problems, the determination of heat transfer and mass transport play a very important role. LetH and J denote the rate of heat and mass transport per unit for the fluid phase. For
(47)
H
=
-
κ
T
〈
∂
T
total
∂
z
〉
|
z
=
0
,
J
=
-
κ
S
〈
∂
S
total
∂
z
〉
|
z
=
0
,
where the angular brackets correspond to a horizontal average and
(48)
T
total
=
T
0
-
Δ
T
z
d
+
T
(
x
,
z
,
t
)
,
S
total
=
S
0
-
Δ
S
z
d
+
S
(
x
,
z
,
t
)
.
Substituting (32)–(34) into (48) and using the resultant equation (48) yield
(49)
H
=
κ
T
Δ
T
d
(
1
-
2
π
B
2
)
,
J
=
κ
S
Δ
S
d
(
1
-
2
π
C
2
)
.
The Nusselt and Sherwood numbers are defined by
(50)
Nu
=
H
d
κ
T
Δ
T
=
1
-
2
π
B
2
,
Sh
=
J
d
κ
S
Δ
S
=
1
-
2
π
C
2
.
Substituting B
2 and C
2 of (41) into (51) gives
(51)
Nu
=
1
+
4
π
2
a
2
F
(
z
)
A
1
2
2
(
4
π
2
-
R
i
)
(
R
i
-
δ
2
)
-
π
2
a
2
A
1
2
,
Sh
=
1
+
2
a
2
Le
2
A
1
2
8
δ
2
+
a
2
Le
2
A
1
2
.


## 4.1. Steady Finite Amplitude Motions

The simplified model represented by (37) has the great advantage that steady finite amplitude solutions can be obtained at once and their stability can be investigated analytically. From qualitative prediction, we look into the possibility of an analytical solution. In the case of steady motions, setting the left hand side of (37) equal to zero yields(40a)
δ
2
A
1
+
a
Ra
B
1
-
a
Ra
s
Le
C
1
=
0
,
(40b)
2
a
F
(
z
)
A
1
-
π
a
A
1
B
2
+
(
R
i
-
δ
2
)
B
1
=
0
,
(40c)
π
a
2
A
1
B
1
+
(
R
i
-
4
π
2
)
B
2
=
0
,
(40d)
-
π
a
η
A
1
C
2
-
δ
2
η
Le
C
1
-
a
η
A
1
=
0
,
(40e)
π
a
2
η
A
1
C
1
-
4
π
2
Le
η
C
2
=
0
.On solving for the amplitudes in terms of A
1, we obtain
(41)
B
1
=
2
a
F
(
z
)
(
8
π
2
-
2
R
i
)
A
1
π
2
a
2
A
1
2
-
(
R
i
-
δ
2
)
(
8
π
2
-
2
R
i
)
,
C
1
=
-
8
a
Le
A
1
8
δ
2
+
a
2
Le
2
A
1
2
,
B
2
=
2
π
a
2
A
1
2
F
(
z
)
π
2
a
2
A
1
2
-
(
R
i
-
δ
2
)
(
8
π
2
-
2
R
i
)
,
C
2
=
-
a
2
Le
2
A
1
2
π
(
8
δ
2
+
a
2
Le
2
A
1
2
)
.
Substituting B
1 and C
1 into (39) yields
(42)
a
2
x
2
+
a
1
x
+
a
0
=
0
,
where x
=
A
1
2
/
8, and
(43)
a
2
=
δ
2
M
2
M
5
,
a
1
=
δ
4
M
2
-
δ
2
M
3
M
5
+
M
1
M
5
+
M
2
M
4
,
a
0
=
M
1
δ
2
-
M
3
M
4
-
M
3
δ
4
,
M
1
=
2
a
2
Ra
F
(
z
)
(
π
2
-
R
i
4
)
,
M
2
=
π
2
a
2
,
M
3
=
(
R
i
-
δ
2
)
(
π
2
-
R
i
4
)
,
M
4
=
a
2
Ra
S
,
M
5
=
a
2
Le
2
.
The required root of (42) is
(44)
A
1
2
8
=
-
a
1
+
a
1
2
-
4
a
2
a
0
2
a
2
.
If the discriminant equals zero, we obtain the expression for finite amplitude Rayleigh number Ra
T
F, which characterizes the onset of finite amplitude steady motions. The finite amplitude Rayleigh number can be expressed in the following form:
(45)
Ra
T
F
=
-
e
1
+
e
1
2
-
4
e
2
e
0
2
e
2
,
where
(46)
e
2
=
a
8
Le
4
F
2
(
z
)
(
π
2
-
R
i
4
)
2
-
4
δ
4
a
6
F
2
(
z
)
Le
2
(
π
2
-
R
i
4
)
2
,
e
1
=
2
a
4
Le
2
F
(
z
)
(
(
π
2
-
R
i
4
)
+
π
2
a
4
R
a
S
)
δ
4
π
2
a
2
-
δ
2
a
2
Le
2
(
R
i
-
δ
2
)
×
(
π
2
-
R
i
4
)
+
π
2
a
4
Ra
S
)
(
π
2
-
R
i
4
)
+
4
δ
2
a
4
F
(
z
)
Le
2
(
π
2
-
R
i
4
)
×
a
2
Ra
S
(
R
i
-
δ
2
)
(
π
2
-
R
i
4
)
+
δ
4
(
R
i
-
δ
2
)
(
π
2
-
R
i
4
)
,
e
0
=
(
(
π
2
-
R
i
4
)
+
π
2
a
4
R
a
S
(
δ
4
π
2
a
2
-
δ
2
a
2
L
e
2
(
R
i
-
δ
2
)
)
δ
4
π
2
a
2
-
δ
2
a
2
Le
2
(
R
i
-
δ
2
)
×
(
π
2
-
R
i
4
)
+
π
2
a
4
Ra
S
(
δ
4
π
2
a
2
-
δ
2
a
2
L
e
2
(
R
i
-
δ
2
)
)
2
.


## 4.2. Heat Transfer and Mass Transport

In the study of convection problems, the determination of heat transfer and mass transport play a very important role. LetH and J denote the rate of heat and mass transport per unit for the fluid phase. For
(47)
H
=
-
κ
T
〈
∂
T
total
∂
z
〉
|
z
=
0
,
J
=
-
κ
S
〈
∂
S
total
∂
z
〉
|
z
=
0
,
where the angular brackets correspond to a horizontal average and
(48)
T
total
=
T
0
-
Δ
T
z
d
+
T
(
x
,
z
,
t
)
,
S
total
=
S
0
-
Δ
S
z
d
+
S
(
x
,
z
,
t
)
.
Substituting (32)–(34) into (48) and using the resultant equation (48) yield
(49)
H
=
κ
T
Δ
T
d
(
1
-
2
π
B
2
)
,
J
=
κ
S
Δ
S
d
(
1
-
2
π
C
2
)
.
The Nusselt and Sherwood numbers are defined by
(50)
Nu
=
H
d
κ
T
Δ
T
=
1
-
2
π
B
2
,
Sh
=
J
d
κ
S
Δ
S
=
1
-
2
π
C
2
.
Substituting B
2 and C
2 of (41) into (51) gives
(51)
Nu
=
1
+
4
π
2
a
2
F
(
z
)
A
1
2
2
(
4
π
2
-
R
i
)
(
R
i
-
δ
2
)
-
π
2
a
2
A
1
2
,
Sh
=
1
+
2
a
2
Le
2
A
1
2
8
δ
2
+
a
2
Le
2
A
1
2
.


## 5. Results and Discussion

### 5.1. Linear Stability Analysis

The linear stability analysis of double diffusive convection in a binary Maxwell fluid saturated porous layer with internal heat source has been studied analytically. In this section, we discuss the effects of the parameters in the governing equations on the onset of the double diffusive convection numerically and graphically. Figure1 shows the neutral curves for stationary and oscillatory mode for fixed values of Va
=
12, Ra
S
=
500, λ
=
0.01, R
i
=
3, Le
=
7, and η
=
0.8 with variation in one of the parameters. From Figure 1, it can be found that the increase in the value of normalized porosity parameter η, stress relaxation time λ, and Lewis number Le decreases the oscillatory critical Rayleigh number. This indicates that those parameters advance the onset of double diffusive convection. Figure 1(d) depicts the stationary Rayleigh number increases with an increasing solutal Rayleigh number Ra
S, which indicates that the effect of the solutal Rayleigh number is to enhance the stability of the system. On the other hand, the oscillatory Rayleigh number decreases with the increasing solutal Rayleigh number, which means that the solutal Rayleigh number has a stabilizing effect on the oscillatory onset. In Figure 1(f), we observe that the stationary Rayleigh number decreases with the increase in internal Rayleigh number R
i, which indicates that the internal Rayleigh number destabilizes the system. Additionally, the increasing internal Rayleigh number decreases the minimum of oscillatory Rayleigh number, which means that the internal Rayleigh number has a destabilizing effect on the oscillatory onset. Figure 1(e) shows the effect of Vadasz number on the critical Rayleigh number of stationary mode for fixed values of other parameters. When a is small, the increasing Va decreases the critical Rayleigh number, indicating that the effect of the increasing Va is to stabilize the system. However, for large value of the wavenumber, the facts are just the opposite.Variation ofRa
T with respect to wave number a for neutral stability curves for different values of parameters.
(a)
(b)
(c)
(d)
(e)
(f)


### 5.2. Nonlinear Analysis

The weak nonlinear analysis provides the quantification of heat and mass transport. The effects of various parameters on the rate of heat and mass transfer are shown in Figure2. The values of parameters are fixed at R
i
=
3, Ra
S
=
100, η
=
0.8, λ
=
0.01, and Le
=
3 with variation in one of the parameters. Figures 2(a) and 2(b) show that an increase in the value of Lewis number Le increases the value of Nu and Sh. Thus, the effect of Lewis number is to increase the rate of heat and mass transfer. Figures 2(c) and 2(d) present that both the rate of heat and mass transfer decrease with the increasing solutal Rayleigh number. Figures 2(e) and 2(f) show that an increase in the value of the internal Rayleigh number R
i increases both the rate of heat and mass transfer.Variation in Nusselt numberNu and Sherwood number Sh with respect to Rayleigh number R
a
T for steady nonlinear stability at different values of parameters.
(a)
(b)
(c)
(d)
(e)
(f)Using the Runge-Kutta method with suitable initial conditions, we solve the autonomous system numerically given by (37). For fixed parameters R
i
=
3, Ra
S
=
100, η
=
0.8, λ
=
0.1, and Le
=
3 with variation in one of the parameters, the results of Nusselt number Nu and Sherwood number Sh with respect to time t are graphically shown in Figures 3 and 4, respectively. Although the maximum amplification of Nu and Sh occurs at the beginning time, it generates spatial oscillations of increasing frequency. Eventually the oscillations reach a steady state when time is large enough. Figures 3 and 4 show that the parameters have insignificant effect on heat transfer; however, the effects of the parameters on mass transfer are smaller. Notice by inspection of Figures 3 and 4 that Nu and Sh are considerably reduced as the value of relaxation time λ increases, and so are η and R
S. However, the internal Rayleigh number Ra
i increases the rate of heat and mass transfer.Variation in Nusselt numberNu with respect to time t for unsteady nonlinear stability at different values of parameters.
(a)
(b)
(c)
(d)Variation in Nusselt numberSh with respect to time t for unsteady nonlinear stability at different values of parameters.
(a)
(b)
(c)
(d)


## 5.1. Linear Stability Analysis

The linear stability analysis of double diffusive convection in a binary Maxwell fluid saturated porous layer with internal heat source has been studied analytically. In this section, we discuss the effects of the parameters in the governing equations on the onset of the double diffusive convection numerically and graphically. Figure1 shows the neutral curves for stationary and oscillatory mode for fixed values of Va
=
12, Ra
S
=
500, λ
=
0.01, R
i
=
3, Le
=
7, and η
=
0.8 with variation in one of the parameters. From Figure 1, it can be found that the increase in the value of normalized porosity parameter η, stress relaxation time λ, and Lewis number Le decreases the oscillatory critical Rayleigh number. This indicates that those parameters advance the onset of double diffusive convection. Figure 1(d) depicts the stationary Rayleigh number increases with an increasing solutal Rayleigh number Ra
S, which indicates that the effect of the solutal Rayleigh number is to enhance the stability of the system. On the other hand, the oscillatory Rayleigh number decreases with the increasing solutal Rayleigh number, which means that the solutal Rayleigh number has a stabilizing effect on the oscillatory onset. In Figure 1(f), we observe that the stationary Rayleigh number decreases with the increase in internal Rayleigh number R
i, which indicates that the internal Rayleigh number destabilizes the system. Additionally, the increasing internal Rayleigh number decreases the minimum of oscillatory Rayleigh number, which means that the internal Rayleigh number has a destabilizing effect on the oscillatory onset. Figure 1(e) shows the effect of Vadasz number on the critical Rayleigh number of stationary mode for fixed values of other parameters. When a is small, the increasing Va decreases the critical Rayleigh number, indicating that the effect of the increasing Va is to stabilize the system. However, for large value of the wavenumber, the facts are just the opposite.Variation ofRa
T with respect to wave number a for neutral stability curves for different values of parameters.
(a)
(b)
(c)
(d)
(e)
(f)


## 5.2. Nonlinear Analysis

The weak nonlinear analysis provides the quantification of heat and mass transport. The effects of various parameters on the rate of heat and mass transfer are shown in Figure2. The values of parameters are fixed at R
i
=
3, Ra
S
=
100, η
=
0.8, λ
=
0.01, and Le
=
3 with variation in one of the parameters. Figures 2(a) and 2(b) show that an increase in the value of Lewis number Le increases the value of Nu and Sh. Thus, the effect of Lewis number is to increase the rate of heat and mass transfer. Figures 2(c) and 2(d) present that both the rate of heat and mass transfer decrease with the increasing solutal Rayleigh number. Figures 2(e) and 2(f) show that an increase in the value of the internal Rayleigh number R
i increases both the rate of heat and mass transfer.Variation in Nusselt numberNu and Sherwood number Sh with respect to Rayleigh number R
a
T for steady nonlinear stability at different values of parameters.
(a)
(b)
(c)
(d)
(e)
(f)Using the Runge-Kutta method with suitable initial conditions, we solve the autonomous system numerically given by (37). For fixed parameters R
i
=
3, Ra
S
=
100, η
=
0.8, λ
=
0.1, and Le
=
3 with variation in one of the parameters, the results of Nusselt number Nu and Sherwood number Sh with respect to time t are graphically shown in Figures 3 and 4, respectively. Although the maximum amplification of Nu and Sh occurs at the beginning time, it generates spatial oscillations of increasing frequency. Eventually the oscillations reach a steady state when time is large enough. Figures 3 and 4 show that the parameters have insignificant effect on heat transfer; however, the effects of the parameters on mass transfer are smaller. Notice by inspection of Figures 3 and 4 that Nu and Sh are considerably reduced as the value of relaxation time λ increases, and so are η and R
S. However, the internal Rayleigh number Ra
i increases the rate of heat and mass transfer.Variation in Nusselt numberNu with respect to time t for unsteady nonlinear stability at different values of parameters.
(a)
(b)
(c)
(d)Variation in Nusselt numberSh with respect to time t for unsteady nonlinear stability at different values of parameters.
(a)
(b)
(c)
(d)


## 6. Conclusion

Linear and nonlinear analysis of double diffusive convection in a Maxwell fluid saturated porous layer with internal heat source, which is heated and salted from below, is investigated analytically and numerically. The linear analysis is analyzed using the normal mode technique. On the other hand, the nonlinear analysis of the system is established through a truncated form of the Fourier series. The effects of physical parameters in governing equations, such as relaxation time, Lewis number, normalized porosity parameter, Vadasz number, solutal Rayleigh number, and internal Rayleigh number, on stationary, oscillatory convection, and heat and mass transfer are shown graphically and the following conclusions are drawn: Vadasz numberVa, relaxation parameter λ, and Lewis number Le advance the onset of oscillatory convection. The internal Rayleigh number R
i has a destabilizing effect for the governing system, both stationary and oscillatory mode, and the solutal Rayleigh number Ra
S stabilizes the system. The rate heat and mass transfer increase with the increasing internal Rayleigh number, but Nu and Sh are reduced as the solutal Rayleigh number, stress relaxation time, and normalized porosity parameter increase.



---

*Source: 489279-2014-04-09.xml*