<< . .

. 23
( : 51)



. . >>

(second term on left side), and the applied load 2
are solving for only half of the pro¬le (i.e. the
q(x) (right side).
other half of the pro¬le supports the other half
of the load). Substituting Eq. (5.11) into Eq. (5.10)
5.2 Methods for 1D problems gives
V ± 3 ’x/± x x
w(x) = cos + sin
e (5.12)
± ±
5.2.1 Displacement under line loading 8D
for x ≥ 0 (Figure 5.2b).
In some cases the topographic load is nar-
row enough compared to the ¬‚exural wave- In the case of de¬‚ection of a diving board,
length of the lithosphere that the load can be we found that the solution is characterized by
112 FLEXURAL ISOSTASY


thickness of 40 km, for example, the ¬‚exural pa-
200
rameter under crustal loading is ≈ 130 km, while
for ice loading it is approximately 90 km.
crust
150
(km) 5.2.2 Series solutions and Fourier ¬ltering
The ¬‚exure equation is a linear equation. As such,
100
ice we can use Fourier series and superposition to de-
velop general solutions for the de¬‚ection of the
50 lithosphere under any topographic load, given
the solution for de¬‚ection under a periodic load.
Let™s assume a periodic topography with ampli-
0
60
20 40 80 tude h 0 and wavelength »:
0
T e (km)
2π x
h(x) = h 0 sin (5.16)
Fig 5.3 Plot of ¬‚exural parameter ± as a function of elastic »
thickness Te for crustal and ice loads, assuming E = 70 GPa
The vertical load corresponding to this topogra-
and ν = 0.25.
phy is given by
a cubic polynomial with no characteristic length 2π x
qa (x) = ρc gh 0 sin (5.17)
»
scale. In the case of ¬‚exural-isostasy, however, the
de¬‚ection is characterized by the length scale ±.
and hence the ¬‚exural-isostatic equation we must
Also, the zone of depression is accompanied by
solve is
a smaller zone of uplift ¬‚anking the margin of
d4 w
the load. This uplift zone is the forebulge. The + (ρm ’ ρc )gw = qa (x)
D (5.18)
dx 4
location of the forebulge is determined by dif-
ferentiating w with respect to x and setting the The displacement of the lithosphere will follow
result equal to zero: the same periodic form as the loading -- the mag-
nitude of the de¬‚ection will depend sensitively on
dw 2w 0 ’x/± x
=’ sin = 0
e (5.13) the wavelength of the load and ¬‚exural param-
± ±
dx
eter of the lithosphere, but there will be no off-
Therefore,
set or shift between the load and the response.
xb = ± sin’1 (0) = π± (5.14) Therefore, we look for solutions of the form
The amplitude of the forebulge can be found by 2π x
w = w 0 sin (5.19)
»
substituting xb into Eq. (5.12) to obtain
w b = ’w 0 e’π = ’0.0432w 0 (5.15) where w 0 is yet to be determined. Substituting
Eq. (5.19) into Eq. (5.18) gives
Therefore, the height of the forebulge is equal
h0
to about 4% of the maximum de¬‚ection of the
w0 = (5.20)
ρm 2π 4
’1+
lithosphere beneath the load. This calculation as- D
ρc ρc g »
sumes that the lithosphere is entirely rigid. In
For small wavelengths, i.e. where
some cases, the loading causes the lithosphere
to break under the stresses of the load. In these 1
D 4
»
cases, the sin term in Eq. (5.12) becomes zero and 2π (5.21)
ρc g
the width of the zone of subsidence becomes nar-
rower (Turcotte and Schubert, 2002). the denominator in Eq. (5.20) is very large
The relationship between ± and T e is plotted and hence w 0 h 0 . This means that small-
in Figure 5.3 for loading by crustal rocks (i.e. wavelength topographic variations cause very lit-
ρ ≈ 600 kg/m3 ) and ice (i.e. ρ ≈ 2300 kg/m3 ), as- tle displacement of the crust and we say that the
suming E = 70 GPa and ν = 0.25. For an elastic crust is rigid for loads of this scale. At large scales,
5.3 METHODS FOR 2D PROBLEMS 113


where (so called because it ¬lters out low wavenumbers)
before the inverse Fourier transform is applied to
1
D 4
bring the function back to real space.
» 2π (5.22)
ρc g Appendix 4 provides code that implements
the Fourier ¬ltering technique for any input to-
Eq. (5.20) gives
pographic load. Figure 5.4 illustrates example
ρc
w0 = h0 output of that program for an east--west pro¬le
(5.23)
ρm ’ ρ c
through the central Andes. The east--west topo-
graphic pro¬le of the central Andes is character-
which is equivalent to the case of pure isostatic
ized by a high-elevation plateau (the Altiplano-
balance given by Eq. (5.2). In the large-scale limit,
Puna basin), high peaks on the margins of the
therefore, ¬‚exure plays no role and the topo-
Altiplano-Puna, and an active fold thrust belt to
graphic load is said to be fully compensated. At
the east of the Altiplano-Puna. In this example,
smaller scales, where the de¬‚ection is smaller
the topographic load of the Andes is considered
than the Airy isostatic limit, the topography is
to be all topography above 1 km in elevation. Fig-
said to be only partially compensated.
ure 5.4b plots the modeled de¬‚ection of the crust
The solution for periodic topography is use-
for ρc = 2.7 g/cm3 , ρm = 3.3 g/cm3 , and three hy-
ful for determining the degree of compensation
pothetical values of ±. The depth of total de¬‚ec-
of a topographic load. The real power of the pe-
tion beneath the Andes (≈ 15 km) depends lin-
riodic solution, however, is that it can be used
early on the relative density contrast between the
in conjunction with the Fourier series approach
crust and mantle, but only weakly on the ¬‚ex-
to determine the response of the lithosphere to
ural parameter for wide topographic loads. Low
any load using superposition. If the topographic
values of the ¬‚exural parameter predict narrow
load within a domain of length L is expressed as
foreland sedimentary basins, while higher val-
a Fourier series:
ues predict broad basins (Figure 5.4b--5.4c). The

πx
h(x) = height of the ¬‚exural forebulge (Figure 5.4c) also
an sin (5.24)
L increases with higher values of the ¬‚exural pa-
n=1

rameter. Gravity modeling has been used to con-
with coef¬cients
strain the width of the Andean foreland basin by
L
2 nπ x taking advantage of the density contrast between
an = h(x ) sin dx (5.25)
L L the basin in¬ll and the crust below. Those con-
0

straints, coupled with forward modeling of ¬‚ex-
then, according to Eq. (5.20), the de¬‚ection
ural pro¬les, suggest that ± = 100 km is a good
caused by the load is given by
estimate for the ¬‚exural parameter of the central
∞ nπ x
Andes.
sin
w(x) = L
an (5.26)
ρm nπ 4
’1+ D
n=1 ρc ρc g L


Equation (5.26) can be used to obtain analytic ex-
5.3 Methods for 2D problems
pressions for the ¬‚exural-isostatic displacement
underneath loads with simple (rectangular, tri-
5.3.1 Integral method
angular, etc.) loads.
Most real topography, however, does not lend Thus far we have considered solutions in 1D only.
itself to simple analytic forms. Equation (5.26) These solutions are very useful for long moun-
can be implemented numerically for any topo- tain chains or other landforms where the loading
graphic load, however, using the Fourier ¬lter- is relatively constant in one direction. In many
ing technique. In this technique, the Fast Fourier cases, however, (e.g. an isolated seamount or curv-
Transform (FFT) is used to quickly calculate the ing mountain belt) a full 2D solution is needed.
Fourier coef¬cients an for an input load function. In such cases, 2D methods analogous to those we
Then, each coef¬cient is multiplied by the ¬lter developed for 1D are available.
114 FLEXURAL ISOSTASY



(a) 5

4

h3
(km)
(c) 0.4
2
w
0.2
1 (km)
0
0.0
(b) 0
± = 50 km
75 km
0.2
100 km
w5
0.4
(km) ± = 50 km
1200
800 1000 1400
75 km
100 km x (km)
10


15
600 800 1000 1200
400
200
0
x (km)
Solutions corresponding to a spatially dis-
Fig 5.4 Flexural model of the response of the lithosphere
to the distributed load of the central Andes. (a) Topographic tributed load can be obtained by representing
pro¬le of the central Andes. (b) De¬‚ection of the lithosphere the load as a collection of many point sources
beneath the central Andes for three hypothetical values of
and integrating the response to each point. For
± = 50, 75, 100 km. (c) Close-up plot of the ¬‚exural pro¬le
a spatially distributed topographic load given by
of the foreland basin and forebulge, illustrating the
h(x, y), the de¬‚ection is given by:
dependence of basin width and forebulge amplitude on ±.

∞ ∞
ρc
First we consider the solution to the 2D w(r ) =
ρm ’ ρc
¬‚exural-isostatic equation: 0 0

2
— kei (x ’ x )2 + (y ’ y )2
D ∇ 4 w + ρgw = q ±
(5.27)
— h(x , y )dx dy (5.29)
corresponding to a point load V at r = 0. The
solution is the Kelvin function kei:

V ±2 This approach is broadly analogous to the series
2r
w(r ) = kei (5.28)
πD ± solution method of Eq. (5.26).
Few analytic solutions can be obtained using
The shape of the kei function is broadly similar Eq. (5.29), both because topographic loads rarely
to Eq. (5.12), i.e. it has a deep zone of subsidence conform to analytic expressions and because the
beneath the load and a smaller forebulge ¬‚ank- kei function is dif¬cult to integrate. Appendix 4
ing the load. provides code for series approximations to the kei
5.3 METHODS FOR 2D PROBLEMS 115


ral parameter, the de¬‚ection follows the disk-like
shape of the load except at the edges where ¬‚ex-
0 0.2
ural bending stresses limit the curvature of the
w/h de¬‚ection pro¬le.
1 R/a =
1.0
4.0
5.3.2 Fourier ¬ltering
2 In the previous section we showed how the de¬‚ec-
tion beneath spatially distributed loads could be
calculated by dividing the load into many small
3
point sources and summing the contributions of
each point. This is a powerful approach, but it
4
can be computationally challenging when calcu-
’1 0 1
lating the response at high spatial resolution over
x/a
a large region. Such calculations can require eval-
uating the kei function 104 --106 times for each of
Fig 5.5 De¬‚ection of disk-shaped loads for various values
104 --106 points. The Fourier-¬ltering method pro-
of R /±.
vides a faster way by taking advantage of the
computational speed of the Fast Fourier Trans-
function that can be used to integrate Eq. (5.29)
form algorithm.
numerically for any topographic load.
In order to implement the Fourier ¬ltering al-
One analytic expression is available, how-
gorithm in 2D, the FFT of the 2D load h(x, y) over
ever, that nicely illustrates the response of the
a domain of width L x and length L y is computed
lithosphere to loads of different size. The re-
to obtain a set of Fourier coef¬cients an,m . Each
sponse of the lithosphere to a disk-shaped load
Fourier coef¬cient is then multiplied by the fac-
depends on the ratio of the disk radius R to the
tor
¬‚exural parameter ±. The solution is given by
⎛ ⎛ ⎞4 ⎞’1
(Hetenyi, 1979): 2 2
⎝ ρm ’ 1 + D ⎝π n m ⎠⎠
+
hρc R R r
ρc ρc g
w(r ) = Lx Ly
ker ber
ρm ’ ρ c ± ± ±
(5.32)
R R r
’ kei +1
bei (5.30)
± ± ± Then, the inverse Fourier transform is applied to
the ¬ltered data set.
for r < R and
Appendix 4 provides a code for implementing
hρc R R r
w(r ) = the Fourier ¬ltering method in 2D. Figure 5.6 il-
ber ker
ρm ’ ρ c ± ± ±
lustrates example output from this program for
R R r
the central Andes for ± = 50, 75, and 100 km. Fig-
’ bei kei (5.31)
± ± ±
ures 5.6b--5.6d are grayscale plots of the absolute
for r > R . value of the de¬‚ection on a logarithmic scale. As
the value of ± increases, the pattern of de¬‚ection
Figure 5.5 plots Eq. (5.31) for loads of differ-
ent radii, expressed as a ratio of the ¬‚exural pa- becomes smoother. In each grayscale image the
rameter ±. For cases in which the load radius forebulge is represented by the ¬rst ring around
is much less than the ¬‚exural parameter, the the central zone of subsidence. Subsequent rings
load is only partially compensated. For cases in are the backbulge and other smaller zones of al-
which the radius is equal to the ¬‚exural parame- ternating uplift and subsidence with increasing
ter, the maximum de¬‚ection is equal to the Airy distance from the load. Figure 5.6e illustrates the
isostatic limit (i.e. Eq. (5.2)), but the maximum de- variations in forebulge amplitude in the vicinity
¬‚ection is attained only beneath the center of the of the Bolivian Orocline using a contour map. The
load. For load radii much larger than the ¬‚exu- forebulge amplitude varies from 160 m to 400 m
116 FLEXURAL ISOSTASY




basin environments are also attractive environ-
Fig 5.6 2D ¬‚exural-isostatic response of the lithosphere to
ments for studying the linkages of climate, tec-
the topographic load of the central Andes (a) for ± = 50 (b),
tonics, and erosion, because the tectonic param-
75 (c), and 100 km (d). In (b)“(d), the grayscale images
eters controlling foreland basin geometry can be
represent the absolute value of the de¬‚ection on a
unusually well constrained. Speci¬cally, thrust-
logarithmic scale. This color scale clearly shows that the
width of the forebulge increases with increasing ±. (e) A belt migration rates can be determined using
contour map of the forebulge in the vicinity of the Bolivian a mass-balance framework that relates crustal
Orocline. This example shows that the height of the shortening to thrust-belt migration (DeCelles and
forebulge changes dramatically (from a minimum of 160 m to
DeCelles, 2001) and basin subsidence is often
a maximum of 400 m) in the vicinity of a curved load.
well characterized using ¬‚exural modeling, grav-
ity measurements, and seismic stratigraphy (e.g.
within a distance of ≈ 100 km, illustrating the
Watts et al., 1995). Sediment transport in ¬‚uvial
importance of 2D effects in the vicinity of a bend
channels has a linear relationship with channel
in the topographic load.
slope, resulting in a diffusion model for the basin
topographic pro¬le (Begin et al., 1981; Paola et al.,
1992).
5.4 Modeling of foreland basin
Foreland basin evolution can be modeled as
geometry a combination of thrust-belt migration, ¬‚exural
subsidence, and diffusive ¬‚uvial sediment trans-
port (Flemings and Jordan, 1989). Modeling this
Flexure plays an important role in large-scale
system requires solving the diffusion equation
landform evolution because, as erosion and de-
with a moving boundary (the thrust front, where
position change the topographic load, ¬‚exural
sediment enters the basin) and moving sink
isostasy causes uplift and subsidence responses
term (the ¬‚exural depression). Typically, mov-
that vary as a function of spatial scale and litho-
ing boundaries and sinks in diffusion problems
spheric rigidity. In this section we explore the
require numerical approaches. However, in this
interaction between ¬‚exural isostasy and sedi-
case the problem can be solved analytically by
mentation in a foreland basin adjacent to a mi-
working in the frame of reference of the migrat-
grating thrust belt. Foreland basins provide a po-
ing thrust front and assuming steady-state con-
tentially useful record of tectonics and erosion
ditions. Conceptually, basins are in steady state
of the adjacent range, but our ability to inter-
when sediment progradation occurs at the same
pret that record uniquely is limited. Foreland

<< . .

. 23
( : 51)



. . >>