Fowler A. Mathematical Geoscience

Подождите немного. Документ загружается.

500 8 Mantle Convection

In the delamination layer, a balance of advection with conduction implies

ψ ∼1/δ, and because buoyancy (



θdψ)is advected round the corner, we have

. (8.181)

Finally, we must balance shear stress with buoyancy in the delamination layer (by

analogy with the boundary layer beneath a rigid lid). Thus p

∼ τ

, p

∼ Raεθ,

and since ∂/∂x ∼1/ε (since ν =δ/ε) and ∂/∂z ∼1/δ, we have in the delamination

layer

∼δ

Ra,p∼εδRa, (8.182)

and together with ψ ∼1/δ, ∂/∂z ∼ 1/δ

, the deﬁnition τ

≈ ηψ

( τ

) implies

∼ τ ∼ η/δ

∼ Λ/τ

n−1

, hence τ

∼ τ ∼ (Λ/δ

)

1/n

, and combining this with

(8.182) and (8.179), we ﬁnd

Ra =



2(n−1)



1/n

. (8.183)

From (8.178), (8.179), (8.180), (8.181) and (8.183), we ﬁnally ﬁnd

δ =



n−1

n+1



1/5

,Λ=Ra

3(n−1)

(2n+3)(n−1)

5(n+1)

. (8.184)

The second of these deﬁnes η

, and ensures that the sub-lithospheric viscosity is

O(η

).From(8.175), we have

Λ =

/κ)

n−1

exp(1/ε). (8.185)

Combining this with (8.184), we ﬁnd



/κ)

n−1

exp



∗





∗



(2n+3)(n−1)

(n+1)



αρ



3(n−1)



2n+3

(8.186)

As discussed above, we have T

≈T

, and more precisely, we can deﬁne



1 +

εδ



, (8.187)

where the unknown O(1) constant φ

is chosen in an analysis of the basal boundary

layer (cf. (8.200) below). In effect we can take T

as known. Taking all the values for

the constants which we have used earlier, we ﬁnd that the expected value of η

for

the Earth is 1.4 ×10

Pa s. Extraordinarily, this is exactly the sort of value which

is thought to be appropriate in the Earth’s asthenosphere. The theory smells right.

8.5 Variable Viscosity 501

Delamination Layer We would now progress through the separate regions: core,

basal boundary layer, plumes, etc. These are much the same as for the constant

viscosity case, but the delamination layer (and the slab above) is different, and we

begin with it.

In the delamination layer, we rescale the variables as follows, based on the dis-

cussion above:

z =1 −νs −δζ, T =1 +εθ, τ =δ

RaT

∗

,η=ε

n−1

n+1

ψ =

,p=εδRaP, τ

=δ

RaT

,τ

;

(8.188)

thus

∂

∂x

→

∂

∂x

−



∂

∂ζ

∂

∂z

→−

∂

∂ζ

, (8.189)

and at leading order we ﬁnd (assuming δ ε)

−s



∂P

∂ζ

=−

∂T

∂ζ

−

∂P

∂ζ

=θ,

=−2N

∂

∂ζ

=−N

∂

∂ζ

∂Ψ

∂ζ

∂θ

∂x

−

∂Ψ

∂x

∂θ

∂ζ

∂

∂ζ

∗

=|T

N =

∗n−1

−θ

(8.190)

We write T

=−S; then (8.190) can be shrunk to the three coupled equations

∂S

∂ζ



θ,

∂

∂ζ

=|S|

n−1

∂Ψ

∂ζ

∂θ

∂x

−

∂Ψ

∂x

∂θ

∂ζ

∂

∂ζ

(8.191)

We anticipate suitable boundary conditions for these equations as follows. As

ζ →−∞, the viscosity increases exponentially as we enter the lid, and we an-

ticipate negligible ﬂow and a conductive temperature gradient. Thus

Ψ →0,θ

→Γ (8.192)

502 8 Mantle Convection

as ζ →−∞. Below the lid lies the isothermal core, and hence we can suppose

S →S

∞

,θ→0asζ →+∞, (8.193)

where S

∞

must be chosen. Because p,τ

,τ

∼1/δ

in the core, and because (from

(8.180) and (8.181))

Ra ε

, (8.194)

we ﬁnd, using the deﬁnition of δ in (8.184), that the delamination layer scaled vari-

ables P and T

take values

P ∼νε

n+1

∼ε

n+1

(8.195)

in the core. Thus T

must reach O(ε

n+1

) 1 in the core, and an obvious choice for

∞

is thus S

∞

= 0. Despite its appeal, we shall ﬁnd that this is incorrect, because

of a further buoyant layer below the delamination layer which arises from the up-

welling plume. Before considering the plume structure, we solve for the temperature

in the lid.

Stagnant Lid Temperature Anticipating that ψ  1 (fuller discussion follows

later), then ∇

T =0 in the lid. We rescale

z =1 −νZ, (8.196)

so that

+ν

=0, (8.197)

and at leading order

T =T

+(1 −T

)

, (8.198)

where T

is the non-dimensional surface temperature. It follows from the

deﬁnition of (8.192) that

Γ =



∂T

∂Z



Z=s

1 −T

. (8.199)

Sub-lithospheric Flow In the plumes and the basal thermal boundary layer,

θ ∼ δ

/δ  1, so that viscosity is approximately constant, and the ﬂow is directly

comparable to the isoviscous case. However, the upper boundary for the core ﬂow

is now one of no slip, so that the ﬂow is not symmetric, and the buoyancy in each

plume may be different.

Put

θ =

φ, (8.200)

8.5 Variable Viscosity 503

Fig. 8.4 The style of

stagnant lid convection for a

strongly variable viscosity

ﬂuid

and write

ψ =

Ψ (8.201)

in plumes and basal boundary layer. Suppose that



∞

φdΨ =C

, −C

, (8.202)

give the values of plume buoyancy in the left and right plumes (C

> 0).

Note that if we choose T

such that θ ≡ 0 in the core, then the value of φ at the

base is undetermined, let us say φ = φ

at z = 0. (This just reﬂects the fact that

the internal core temperature T

is not exactly known.) Thus we have to determine

and C

The core ﬂow problem for

ψ =δ

ψ is the following:

∇

ψ =0, (8.203)

with

ψ =0,

on x =0,

ψ =0,

=−C

on x =a,

ψ =0,

=0onz =0,

ψ =0,

=0onz =1.

(8.204)

Given the solution to this, then the anti-clockwise boundary velocity on the edges

OABC (see Fig. 8.4)isu =∂

ψ/∂n, where n is the inward normal, and if we deﬁne

504 8 Mantle Convection

τ (no longer the stress!) by

τ =



uds, (8.205)

where s measures distance anti-clockwise along OABC from O, then φ satisﬁes

∂φ

∂τ

∂

∂Ψ

in 0 <Ψ <∞, 0 <τ<τ

, (8.206)

where Ψ = δ

ψ is the boundary layer coordinate, and if τ = 0,τ

,τ

at the

points O,A,B,C, then suitable boundary conditions are

φ →0asΨ →∞,

=0onΨ =0, 0 <τ<τ

and τ

<τ <τ

φ =φ

on Ψ =0,τ

<τ <τ

(8.207)

Outer Thermal Layer At this point we must enquire what happens to the up-

welling plume after it impinges the lid at x = 0. We have already described the

delamination layer, in which ψ ∼ 1/δ, θ ∼ 1 and 1 − νs − z ∼ δ. In the plume,

however, ψ ∼1/δ

and θ ∼δ

/δ. The plume turns the corner and is carried across

in an outer thermal layer below the delamination layer. In this layer, we write (from

the delamination layer variable)

ζ =





1/2

η, (8.208)

which is appropriate if S

∞

=0in(8.193).

We then ﬁnd that at leading order,

−s



∂P

∂η

=−

∂T

∂η

−

∂P

∂η

=φ,

=−N

∂

∂η

=−2Ns



∂

∂η

∗

=|T

N =

∗n−1

−Ψ

=0,

(8.209)

provided ε (δ

/δ)

1/2

, which we assume.

8.5 Variable Viscosity 505

It ﬁrstly follows that φ ≈ φ(Ψ), and the plume buoyancy structure is advected

across the top surface unchanged. Secondly, S =−T

and Ψ satisfy the equations

∂

∂η

=|S|

n−1

∂S

∂η



φ(Ψ),

(8.210)

with matching conditions

Ψ =Ψ

=0,S=S

∞

on η =0,

S →0asη →∞.

(8.211)

The extra condition determines the value of S

∞

, which is then used to solve the

delamination layer equations (8.191).

Plume Circulation Returning now to the model for φ in the plumes, (8.206) and

(8.207), we see that if the upwelling plume value at x = 0, τ =τ

is φ

(Ψ ), then



∞

(Ψ ) dΨ =−C

. (8.212)

The change in the buoyancy



∞

φdΨ across the top surface is (C

), but this

is manifested in the delamination layer (where Ψ ∼δ

/δ 1). Therefore the initial

condition for the plume at x =a, τ =0, may be approximately represented as

φ|

τ =0

=φ

(Ψ ) −(C

)δ

(Ψ ), (8.213)

where the half-range delta function δ

(Ψ ) is zero for Ψ>0, and



∞

(Ψ ) dΨ =1.

In addition, we have the heat ﬂux conditions

=−



∂φ

∂Ψ



Ψ =0

dτ (8.214)

at the base, and similarly at the top



Γdx=(1 −T

)



. (8.215)

We can simplify the model as follows. Deﬁne

γ =

, 2C =C

. (8.216)

Given s,(8.215) deﬁnes C. The core stream function can be written as

ψ =C

1/2

, (8.217)

506 8 Mantle Convection

where ψ

(γ is merely a labelling superscript, not an exponent) depends only on γ ,

the relevant plume conditions being

=2γ on x =0,

=2(1 −γ) on x =a.

(8.218)

Next we deﬁne σ by

φ =2Cσ −

√

πτ

exp



−

4τ



, (8.219)

so that σ satisﬁes

=σ

ΨΨ

in 0 <Ψ <∞, 0 <τ<τ

, (8.220)

with boundary conditions

σ →0asΨ →∞,

=0onΨ =0, 0 <τ<τ

σ =α +

√

πτ

on Ψ =0,τ

<τ <τ

=0onΨ =0,τ

<τ <τ

(8.221)

where

=2αC, (8.222)

and we require, from (8.212)–(8.214),

σ |

τ =0

=σ |

τ =τ

−

√

πτ

exp



−

4τ



, (8.223)

−



∂σ

∂Ψ



Ψ =0

dτ =1, (8.224)



∞

σ |

τ =0

dΨ =γ. (8.225)

Now since

ψ =C

1/2

, it follows that

τ =C

1/2



∂ψ

∂n

ds, (8.226)

so we deﬁne

τ =C

1/2

ˆτ, Ψ =C

1/4

χ, σ =Σ(χ, ˆτ)/C

1/4

, (8.227)

8.5 Variable Viscosity 507

and τ

1/2

ˆτ

,etc.Σ satisﬁes

ˆτ

=Σ

χχ

, (8.228)

with

Σ →0asχ →∞,

=0onχ =0, 0 < ˆτ<ˆτ

Σ =ˆα +

√

π ˆτ

on χ =0, ˆτ

< ˆτ<ˆτ

=0onχ =0, ˆτ

< ˆτ<ˆτ

(8.229)

where ˆα =αC

1/4

, and the constraints

Σ|

ˆτ =ˆτ

=Σ|

ˆτ =0



π ˆτ

exp



−

4 ˆτ



, (8.230)

−



ˆτ

∂Σ

∂χ



χ=0

d ˆτ =1, (8.231)



∞

Σ|

τ =0

dχ =γ. (8.232)

The value of ˆα is chosen so that (8.231) is satisﬁed. Given an initial function

(χ) =Σ|

ˆτ =0

,wesolve(8.229) with (8.231) till ˆτ =ˆτ

, which determines Σ|

ˆτ =ˆτ

as a linear functional of Σ

.(8.230) is then a linear inhomogeneous integral equa-

tion, with (so we suppose) a unique solution, which depends on γ , since ˆτ

, ˆτ

and

ˆτ

depend on γ .(8.232) then provides an equation for γ . Note that the plume head

value of φ

(Ψ ) used in (8.210) is given by

(Ψ ) =2C

3/4



Ψ/C

1/4



. (8.233)

It remains to determine s and thus C,via(8.215):

2C =(1 −T

)



. (8.234)

This appears to require that s can be fully determined through the solution of

(8.191). We now consider whether this can be true.

Delamination Layer: Similarity Solution The delamination equations (8.191)

form a ﬁfth order set of differential equations. One can argue that the boundary

conditions in (8.192) and (8.193) actually constitute six conditions, since the implied

conditions Ψ →0asζ →−∞and θ

→0asζ →∞are understood. If that is the

case, then we can expect that s



in (8.191) is a nonlinear eigenvalue for the model,

and thus (since also Γ depends on s), that the solution will implicitly determine a

differential equation of the form s



=Υ(s)for s.

508 8 Mantle Convection

If that is the case, then the question arises what is the appropriate boundary con-

dition for s? One possibility is that s(0) =0; another would be that s



(0) =0, which

would determine s(0) as the root (if it exists) of Υ = 0. Without detailed matching

of the delamination layer to the corner ﬂow, it is difﬁcult to be more precise.

Some insight into this question can be gained if we suppose that S

∞

= 0in

(8.193). This would follow if the function φ

(Ψ ) used in the solution of (8.210)

is identically zero. In this case, Eqs. (8.191) have a similarity solution, and this is

given by

ξ =Γ(x)ζ, θ =g(ξ), T

−s



h(ξ), Ψ =









n−1

f(ξ).

(8.235)

The functions f, g, h then satisfy the equations



=g,



=|h|

n−1





=0,

(8.236)

with

g(∞) =h(∞) =0,f(−∞) =0,g



(−∞) =1, (8.237)

providing s satisﬁes the equation











n−1





Γ. (8.238)

Given Γ , this equation determines the location of the (unknown) lithosphere base.

To solve (8.236), we use a shooting method for the equation with B

=1:



=G,



=|H |

n−1



+FG



=0,

(8.239)

with

F =F



=0,G



=c, G =−cM, H =H

at ξ =−M,

(8.240)

and M is chosen to be large (M = 20 is adequate). The values of H

and c are

adjusted via Newton iteration at ξ = M until G = H =0 there. Once a solution is

found, then f, g, h are determined by

g(ξ) =G(ξ/c), h(ξ ) =cH (ξ /c), f (ξ) =c

n+2

F(ξ/c), (8.241)

8.5 Variable Viscosity 509

Table 8.1 Va lu es o f B

1 0.087

24.74 ×10

−3

31.39 ×10

−4

3.5 1.98 ×10

−5

42.53 ×10

−6

53.15 ×10

−8

and

n+3

. (8.242)

AsshowninTable8.1, B

is small and decreases rapidly as n increases.

Lithosphere Base We substitute (8.199)into(8.238), so that



n+2

n





(1 −T

)

n+3

(8.243)

(since Γ>0, and thus we can take s



> 0). In order to complete the solution for

s(x), we require two boundary conditions. It is not entirely obvious where these are

to come from. One reasonable choice would appear to be s



(0) =0, which ensures

that Ψ → 0asx → 0 (as would appear necessary, not only because of the vertical

boundary at x =0, but also because Ψ ∼δ

/δ 1 in the plume). If s =s

at x =0,

then the solution of (8.243)as



n+2



n+2

−s

n+2



n+1

=bx, (8.244)

where

b =



n +1

n +2



(1 −T

)

n+3



n+1

. (8.245)

More likely, however, is that we should choose s(0) = 0, on the basis that the

similarity variable is ξ =

(1−T

)ζ

, and identiﬁcation of x →0 with ζ →∞, as usual

in similarity solutions, requires this condition. In that case, the solution for s is

s =



(2n +3)bx

n +1



n+1

2n+3

, (8.246)

and the determination of the solution is complete.