Geoscience Reference
In-Depth Information
thermal insulation at horizontal rigid surfaces,
∂T
∂z
=0 at
z
=
and a uniform mesh in the azimuthal direction according
to the Fourier series:
±
1,
φ
k
=
2
kπ
K
for
k
∈[
0,
K
]
.
and constant-temperature conditions at the vertical
sidewalls,
The time integration used is second-order accurate and
is based on a combination of Adams-Bashforth (AB) and
backward differentiation formula (BDF) schemes, cho-
sen for its good stability properties [
Vanel et al.
, 1986].
The resulting AB/BDF scheme is semi-implicit, and for
the transport equation of the velocity components in
equation (16.1),
T
=
±
1 t
r
=
±
1.
In the case of an open upper wall, planar free-surface
conditions are imposed, assuming the absence of verti-
cal deformation along this boundary, in agreement with
experimental observations for the control parameter val-
ues under consideration [
Harlander et al.
, 2011]:
∂T
4
f
l
+
f
l
−
1
2
δt
3
f
l
+1
−
∂z
=
∂V
r
∂z
=
∂V
φ
(f
l
)
(f
l
−
1
)
+2
N
−
N
=
V
z
=0 at
z
=1.
∂z
4
A
3
/
2
Ta
1
/
2
∇
l
+1
+
2
f
l
+1
+
F
l
+1
i
=
−∇
(16.4)
16.2.4. Solution Method
where
(f )
is a “global term” including the advection
terms
N(
V
)
and the Coriolis term,
F
i
corresponds to the
component of the forcing term
F
,
δt
is the time step,
and the superscript
l
refers to time level. The cross terms
in the diffusion part in the
(r
,
φ)
plane resulting from
the use of cylindrical coordinates in equation (16.1) are
treated within
N
A pseudospectral collocation Chebyshev method is
implemented in the meridional plane
(r
,
z)
, in association
with Galerkin-Fourier approximation in azimuth
φ
for
the three-dimensional flow regimes to solve the primitive
variable formulation described above. Each dependent
variable is expanded in the approximation space
P
NM
,
composed of Chebyshev polynomials of degrees less than
or equal to
N
and
M
in the
r
and
z
directions, respec-
tively, while Fourier series are introduced in the azimuthal
direction with
K
modes.
For each dependent variable
f
(
f
(f )
. The latter is discretized in time using
a second-order explicit AB scheme in order to main-
tain an overall second-order time accuracy with the BDF
scheme applied to the diffusion term, as shown in equa-
tion (16.4). An equivalent discretization applies for the
transport equation (16.3) of the temperature. For the ini-
tialstep,wehavetaken
f
−
1
=
f
0
. At each time step, the
problem then reduces to the solution of Helmholtz and
Poisson equations. We recall that the same time scheme is
also implemented for the integration of the axisymmetric
system, even though steady state solutions are sought.
An efficient projection scheme is introduced to solve
the coupling between the velocity and the pressure in
equation (16.1). This algorithm ensures a divergence-free
velocity field at each time step, maintains the order of
accuracy of the time scheme for each dependent variable
and does not require the use of staggered grids [
Hugues
and Randriamampianina
, 1998;
Raspo et al.
, 2002]. At
each time step, a preliminary Poisson equation for the
pressure, directly derived from the Navier-Stokes equa-
tions, is first solved before integrating the governing sys-
tem described above. It allows for a variation in time of
the normal pressure gradient at boundaries [
Hugues and
Randriamampianina
, 1998], which plays an important role
for time-dependent flows, in particular in the presence
of an open free surface. A complete diagonalization of
operators yields simple matrix products for the solution
of successive Helmholtz and Poisson equations at each
time step [
Haldenwang et al.
, 1984]. The computations of
eigenvalues, eigenvectors, and inversion of corresponding
N
≡
V
r
,
V
φ
,
V
z
,
T
,
p
),
it reads
f
NMK
(r
,
φ
,
z
,
t)
K/
2
−
1
N
M
f
nmk
(t)T
n
(r)T
m
(z)
exp
(ikφ)
,
=
n
=0
m
=0
k
=
−
K/
2
where
T
n
and
T
m
are Chebyshev polynomials of degrees
n
and
m
.
This approximation is applied at collocation points,
where the differential equations are assumed to be satisfied
exactly [
Gottlieb and Orszag
, 1977;
Canuto et al.
, 1987].
We have chosen the Chebyshev-Gauss-Lobatto distribu-
tion defined by a high concentration of points toward the
boundaries, well suited to handle the thin dynamical and
thermal boundary layers expected to develop at high val-
ues of the Taylor Ta and Rayleigh Ra numbers, which
scale as Ta
−
1
/
4
(Ekman layer, along horizontal walls) or
Ta
−
1
/
6
(Stewartson layer, along vertical cylinders) and
Ra
−
1
/
4
, respectively:
r
i
= cos
iπ
N
for
i
∈[
0,
N
]
,
z
j
= cos
jπ
M
for
j
∈[
0,
M
]
,
