Tectonophysics,
164 (1989) 189-201
Elsevier Science Publishers
189
B.V., Amsterdam
- Printed
in The Netherlands
Modeling of the geothermal history of the Okhotsk Sea lithosphere I.K. TUEZOV ’ and V.D. EPANESHNIKOV
*
’ Institute of Tectonics and Geophysics, U.S..% R. Academy of Sciences, Khabarousk (U.S.S.R.) 2 Khabarovsk Polyiechnic Institute, Ministry of Education, Khabarovsk (U.S.S.R.) (Received
November
2, 1987; revised version accepted
July 18, 1988)
Abstract Tuezov, I.K. and Epaneshnikov, term&k,
L. Rybach
and
V.D., 1989. Modeling E.R.
Decker
of the geothermal
(Editors),
Heat
Flow
and
history
of the Okhotsk
Lithosphere
Sea lithosphere.
In: V.
Tectonophysics, 164:
Structure.
189-201. A method
of determination
flow data on the Earth’s Sea lithosphere
surface
is carried
south of the Okhotsk
of evolution
Modeling
out. It is determined
Sea and 30-35 melting
melting
that the partial
is mapped.
zone within
of the evolution
m.y. ago in the north.
zone at the initial stage of its development layer above the partial
of the partial
is presented.
melting
The location
The decrease
the lithosphere
of the partial zone began
melting developing
with time in the thickness
kali
basaltic
volcanism
stage. At this initial rifts and zones of epiplatform
oro-
tion on the basis of geological and geophysical data and the other involves mathematical modeling of the evolution of the lithosphere in the active regions. Analysis of the temporal evolution of the composition of the mantle and crustal xenoliths in rocks has been carried out by Gliko and (1987). The following stages of magmatic
of the granitic
melting (dioritic)
0 1989 Elsevier Science
Publishers
is characteristic
stage,
thinning
of this
of the litho-
modeling of the thermophysical processes occurring in active regions has been carried out. The first mathematical approach involves deep temperature estimations carried out on the basis of the solution of a one-dimensional heat conduction equation Cermak,
(e.g., Smirnov and 1982). In this approach,
rence of the upper boundary is determined in the active
evolution are distinguished. The generation of magma occurs at depths of about 100 km, the temperatures being 1250 + 100” C at initial stage of activation. Melts strongly unsaturated with respect to silica predominate. The process then penetrates to depths of 50-60 km and occurs at temperatures of 1150 “-1200 o C. Widespread al0040-1951/89/$03.50
m.y. ago in the
of the partial
sphere by 50-60 km has already occurred. The duration of this activation process which results in lithospheric thinning and rifting may cover a time period of 30 m.y. Works involving mathematical
genesis are the evidence of activation of stable regions of the Earth. Two approaches to investigating this activation may be made. One projects the temporal succession of the lithospheric activa-
eruptive Grachev
40-45
of the upper boundary
heat
the Okhotsk
zone is recorded.
Introduction Continental
using the observed zone within
B.V.
Sugrobov, a shallow
1980, occur-
of the asthenosphere regions; this was the
basis of further investigations. The second mathematical approach includes construction of two-dimensional steady-state models of the crust and mantle inhomogeneities which cause the heat flow anomalies seen at the Earths surface (Horai, 1976; Lubimova et al., 1983; Tuezov and Epaneshnikov, 1987). These investigations revealed
190
the relationship
between
the shape and the depth
of occurrence
of the inhomogeneities
of the
flow
heat
anomalies.
and the shape
A third
approach
includes modeling of the thermophysical processes based on the solution of a two-dimensional non-steady-state method problem and
conduction
of lithospheric
in which mantle
Models of the partial
thinning
melting
below, and the thermal
(1979)
was per-
of the lithosphere heat flow from
models of mechanical
of the lithosphere
sub-
by the asthenosphere
have been examined by Zorin and Lepina (1985). In order to construct the partial melting model it was necessary to increase the heat flow values to 2 300 mW/m2 at the base of the lithosphere. These values should lead to intense melting of the Earths crust, such an extensive process being proved by the data on volcanism. The authors came to the conclusion that the model of mechanical substitution is preferable. The data obtained on the velocity
history
of substitu-
tion (- 3 km/m.y.) do not contradict the data on the rheology of the lithosphere and asthenosphere (Artyushkov, 1979; Bott, 1981) and allow the consistency between the available geological and geophysical data and the heat flow data to be maintained. The shape and the time of existence of the intrusive body which provides the positive heat flow anomaly on the surface were determined by Kutas and Tsvyaschchenko (1984) and Golubev and Osokina (1980). Modeling of the area beneath the BaikaI rift zone has been carried out. It is shown that to obtain the heat flow anomaly of the observed shape, it is necessary to assume the presence of an intrusive body 6-8 km wide over a period of - 2.5-2.8 m.y. below the Baikal rift. The dependence of the shape of the positive heat flow anomaly on the shape, the depth and the time of existence of the intrusive body located below the anomaly is revealed in this paper. Intrusive bodies with increased temperatures, within the limits of which partial melting zones may exist, play a key role in the formation of positive heat flow anomalies. This paper presents
development
of the method
determination
taking
of thermal
into account
the ob-
served heat flow data. The method of calculation
The
two-phase
heat flow was increased.
due to the effect of the increased stitution
equation.
of the Stephan-type
was worked by Gliko and Efimov,
modeling
formed
heat
of solution
the further
As has already ing zones positive
play
been pointed a key role
heat flow anomalies.
termination be framed evolution
of the thermal
The external
boundary
of
of de-
can therefore
of determination
melting
of the
zones with time.
(Fig. 1) consists
B, = Bi + Bf + B: + Bf.
parts
melt-
The problem history
in the problem of the partial
out, partial
in the formation
of four
Bi
Here,
is the
upper part of the boundary (the land surface and the sea bottom), Bf and B: are the lateral sides of the boundary and B: is the lower part of the boundary at a depth of 400 km. B, is the internal boundary, containing the closed zone inside B,. Boundary B, contours the partial melting zone, and its location may vary with time. The melting set at B, boundary The two-dimensional
temperature is automatically at a corresponding depth.
non-steady-state heat conduction equation with variable coefficients is solved inside the boundaries B, and
B,, which contour
aax 4x3 z> +& [
aT(x
-
[
K(X, z)
z,
ix
Here, tivity,
1
aT(x;Zzy “I +H(x,
=p(x, z)c(x, z) referring
t)
the area Qt:
aqx,z,
for (x,
at
z)
to Q ; t varies 0 to tH .
x, z are coordinates, H is heat generation,
The boundary
and
conditions
z, t) = 0°C
for (x,
(1)
K is thermal conducp is density, C is heat
capacity, T is temperature which eq. (1) is solved.
T(x,
t)
z)
t is the time for
at B, are as follows: z) referring
to B:;
t variesOtot,; aT(x,
z, t) 3X
(Bf, T(x,
B:);
=
0°C
for (x,
z) referring
to
t varies 0 to t,;
I, t) = 1630 “C for (x, z) referring
t varies 0 to t n
.
to Bf; (2)
191
I
0
0: T-0 "C
I
I
I
I
x , WY
I
b
5-l
R
R 400
Fig. 1. Test modeling
1
\
of the evolution
a. Heat flow. b. The partial
of the partial
I
’ B’<
melting
T=l63O'C zone. I-present
day;
2-15
m.y. ago; 3-30
melting zone. Bi, Bf, B: and 8: are the parts of the external 51 is the area of the solution
In these conditions at the upper part of the boundary the temperature was set at O”C, while at the lower part of the boundary a temperature of the olivine-spine1 phase transition equal to
of the heat conduction
boundary,
m.y. ago; 4-45
B, is the internal
m.y. ago.
boundary
and
equation.
1630” C was set (Zharkov, 1983). The heat flow through the lateral sides of the boundary was assumed to be zero. The conditions at boundary B, were constant
192
t~ou~out
the time
corresponded point
O-t,.
interval
to the present
distant
time, while zero is a
boundary
equation.
This
observed
heat flow is a result
initial
3, we obtain:
approximation
steadiness
l”(x, z,
t) = G(x,
referring
z, t) for (x,
.z)
z, t) is the unknown
The initial
conditions
(3) function.
are as follows:
zone and
obtained
(4) where
T,(x,
z) is the solution
by
upper
boundary
depth
of 8-10
internal
boundary melting
a depth of 100 We denote for any G(x, observed heat
corresponding
to the
km. the solution of the problem (l-4) z, t) by T(x, z, t; G), while the flow through B: is denoted by
Q(& z). The problem the partial
conditions
zone in the shape of a thin plate at
the
steady-state
melting
calculations,
of the uplift
of the partial
melting
of the
zone to a
km.
numerically
the problem
by the method
differences
form grid (40 X 40 points) it was solved numerically
(l-5)
was
of successive
Arsenin, by the
1974). locally
scheme on a non-uni(Samarsky,
1977), and
for the boundary
condi-
tions (2 and 3) and the initial conditions (4) using the modified elimination method for a tridiagonal matrix (Samarsky, 1987). The finite”differences analogue of the mized function (5) was assumed to be:
mini-
(6) of determining
melting
the evolution
of
zone is then to find the func-
tion G(x, z, t) for (x, z) referring to B, and t is then varying 0 to tH; the following function minimized :
..f(c;) = iB,)/ K(X,
regime. As an of a weak non-
approximations (Tikhonov and Equation (1) was approximated one-dimensional
of the two-dimen-
sional steady-state heat conduction equation with the variable coefficients for the external (2) and partial
of the stower ther-
to be solved is provided
also by the restriction
solved
to D
that the
of the shape of the partial
In such an arrangement,
T(x, z, 0) = Tc(x, Z) for (x, 2) referring
us to suggest
the setting
of the problem
by the setting
to B,; t varies 0 to t,
where G(x,
fact allows
mal effects in the quasi-steady-state
to t,.
in the past relative
At the internal
t,
Moment
z)
aT(x* “,‘Zt”; G,
1
Here, Qc, Q:’ and AQ,?’ are the calculated heat flow, the observed heat flow and the error in the observed
heat flow measurement respectively. The was determined on the basis of the value Qf solution of eq. (1) and appeared to be the function of the boundary conditions G( x, z, t) at the internal boundary 3,. Summing up was carried out on all the grid
Here, J(G) is the curvilinear integral along the boundary Bi, the integer is a square of the calculated and the observed heat flow differences at the of the moment t,, while dS is the differential integration curve Bi. In such an arrangement, the problem is referred to a number of ill-determined problems, which thus permit only appro~mated solutions; many solutions usually exist. The wide range of the evolutionary heat conduction problems has been investigated by Latters and Lions (1970), and their solvability is proved. The observed heat flow is rather well described by the solution of the steady-state heat conduction
boundary
points
along
the upper
part
of the
B,.
Thus, the problem was to search for the function G(x, z, t), using the method of successive approximations, which describes evolution of the shape of the partial melting zone in the time interval O-t, as well as in the interval duration of O-t, in which the function (6) was minimized. The initial appro~mation of the function G(x, z, tf was taken from the solution of the steady-state heat conduction equation. Function (6) was minimized by the multiple numerical solution of eq. (1) for the time interval O-t,. The shape of the partial melting zone was systematically
varied
in the course
of solution
in
193
order to provide calculated
the best consistency
and
the
observed
flows
the
at the
The problem intervals between
The time interval
possible
to reach
the calculated
was considered development
during
the best
which it
consistency
and the observed
of the partial
to be
the
melting
The arbitrary
time of
zone.
melting
solution
which
Evolu-
zone was was
being
was verified evolution
by test calcula-
of the partial
melt-
ing zone from the thin plate at a depth of 100 km
A map compiled for the calculation Sea lithosphere. mW/m2)
The
flow was given
in
heat flow
t, (the present day). steady-state heat conduction
prob-
lem was first solved in accordance with the method chosen. As an initial approximation, the obtained shape of the partial melting non-steady-state calculations the evolution of the zone.
zone was set by the while determining
Figure 2, shows the dependence of the minimized function (6) on the duration of the time interval O-t, for which the inverse problem was solved. It is seen that the minimum of function (6) corresponds
to the time interval
of 45 m.y., i.e.
one can clearly distinguish the time of evolution of the partial melting zone as approximately coinciding with the weak non-steadiness. Figure 3 shows the determined
evolution
of the
partial melting zone. Comparison with Fig. 1 demonstrates that it correlates well with the given evolution of the zone.
Fig. 2. Test modeling. Dependence of the function minimized on the time of evolution of the partial melting zone.
flow
minima
Kuri-Kamchatka
(30-50
Okhotsk
Sea
area, and along the
arc (Fig. 4). Heat flow maxima,
with values of loo-125 mW/m’, are observed east of Sakhalin, in the Kuril Basin, and in the central
map was compiled
heat
as the observed
heat
the Sikhote-Alin
(Fig.
test calculations
(1980) was used
occur in the northwestern
part of the Okhotsk
1). The calculated
by Smirnov
of the heat field of the Okhotsk
was given in the time interval O-t, (45 m.y. ago-present) and the direct problem was solved
at the moment The inverse
zone in the
Results
and adjacent
searched for. The above method
further
melting
Sea lithosphere.
heat flows
to be the most probable
tion of the shape of the partial assumed
was used for determin-
of the partial
was then solved for the other time
O-t,.
appeared
Thus the above method ing the evolution Okhotsk
t,.
moment
tions.
between
heat
Sea.
The heat flow data
on the basis
were recorded
of which this
in the computer
memory and plotted on the calculated profiles by linear interpolation. Figure 4 shows the calculated profiles
along
simplifies
with
the heat
the comparison
flow contours.
This
with the results
of the
calculations made with the heat flow. The sedimentary, dioritic and basaltic
layers
were plotted on the profiles and the mantle was divided into two layers down to 70 km and below 70 km. The thermophysical properties of the rocks (the thermal conductivity (W/mK), taking into account the mean temperatures, and the heat generation (pW/m3)) state calculations 1987).
were the same as in the steady(Tuezov and Epaneshnikov,
The values
for thermal
conductivity
and
heat generation,
respectively,
are as follows:
ments-1.6
1.2; dioritic
layer-l.6
and
and
basaltic
layer-l.4
granitic
and layer-l.5
0.9;
sedi1.0;
and 0.4; mantle down to 70 km-3.0 and 0.05; and mantle below 70 km-5.0 and 0.005. Within the mantle layers, the initial values of heat generation are indicated. These appeared to be the significant parameters in the steady-state calculations and were specified for every profile. The melting versus depth curve was taken from Yoder (1965) the values being as follows: 0 km-llOO”C, 50 km-1200°C and 100 km1360° C. The calculations were carried out to a depth of 400 km for all the profiles (the olivine-spine1 phase transition depth) (Zharkov, 1983).
194
60-
80 -
too -
Fig. 3. Test modeling. The determined evolution of the partial melting zone. For explanation see Fig. 1.
Figure 5 shows the steady-state temperature field on profile X-X. The partial melting zone was determined while solving the inverse steady-state heat conduction problem. It is seen that the shape of the zone does not correlate well with the loca-
tion of the Moho. It should be emphasized that the granitic layer is absent on profile X-X while the thickness of the dioritic layer is less than 5 km above the partial melting zone. In general, the thicknesses of the dioritic and granitic layers above
195
T
1.x
156”
48 O_
the partial melting zone increase to the north of profile X-X and, to a lesser degree, to the south of it. The problem (l-5) was solved for each profile, i.e. the time interval O-T, was determined for
each profile where the best consistency between the calculated and the observed heat flows was achieved at the moment t, (the present day). Figures 6-8 show the evolution of the upper boundary of the partial melting zone . Figure 9
TABLE 1 Velocities
at which the upper boundary
Time interval
Duration (m.y.)
(m.y. ago)
of the partial melting zone ascended Displacement
Displacement
Displacement
interval (km)
amplitude (km)
velocity (cm/yr)
loo-50
50
1.00
30
0.60
10
0.20
10
0.04
40-35
5
35-30
5
50-20
5
20-10
30-25 25-
0
25
lo-
0
196
100
1-JJ .T*--*\.
3
.A’
.
c
.
A.
.
’
‘\.
*--.-._.
_ km
\
\
300
--
_---
1600
----_
---
--
--
--
1
km
Fig. 5. The steady-state 3 = water;
4 = sedimentary
temperature layer;
field for profile
X-X.
(Fig. 4). I = the observed
heat
5 = dioritic layer; 6 = basaltic layer; 7 = Moho; 8 = additional of 70 km; 9 = partial melting zone; 10 = geotherm.
flow;
2 = the calculated
boundary
in the mantle
heat
flow;
at a depth
197
Fig. 6. Depth
contours
to the upper boundary (thicker
contours)
of the partial
melting
of the upper boundary
shows the location of the upper boundary of the zone determined by the steady-state calculations. Figure boundary
6 shows the location of the upper of the partial melting zone 35 m.y. ago.
Movement of this upper boundary, beginning 40-45 m.y. ago, resulted in the formation of two
zone 35 m.y. ago. 1 = the calculated
of the partial
melting
zone; 3 = isobaths
profiles;
2 = isolines
(km)
(m).
Sea, the upper boundary of the partial melting zone ascended to depths of 50-70 km at this time. Figure
8 shows
the
location
boundary Expansion
of the partial melting of the already-formed
the zone,
formation
of
the
upper
zone 25 m.y. ago. uplifted parts of
of the new areas within
the
major uplifts of the boundary, and one small uplift between the island of Sakhalin and the Kuriles. The upper boundary of the partial melt-
central part of the Okhotsk Sea and in the south of Kamchatka and a general rising up of the upper boundary to depths of 50 km were characteristic
ing zone ascended from a depth of 100 km to 40 km at the central part of the uplift of the boundary during a time span of 5-10 m.y. Figure 7 shows the location of the upper boundary of the partial melting zone 30 m.y. ago. Further upward movement of the upper boundary
of this period. Figure 9 shows the location of the upper boundary of the partial melting zone determined by the steady-state calculations. This stage may be assumed to be the present state of the upper boundary of the zone.
to depths of 20 km occurred up to this period. In the remainder (the major part) of the Okhotsk
Table 1 shows the velocities at which the upper boundary of the partial melting zone ascended. It
198
Fig. 7. Depth contours
to the upper boundary
of the partial
may be seen that the migration velocity of the upper boundary within the depth range of 100-10 km changes by a factor of 25 times in the time span of 40 m.y. The greatest values, equal to 1.0 cm/y, correspond to the initial stage of partial melting zone development, while the lowest values, 0.04 cm/y, correspond to the final stage. The difference in the velocity of the upward movement of the upper boundary of the partial melting zone is probably caused by the dominant role of convective heat transfer through fractures in the initial stage of development, while in the final stage heat transfer was mainly by conductive processes. Comparison of the location of the upper boundary of the partial melting zone 25 m.y. ago with the location determined by the steady-state
melting zone 30 m.y. ago. For explanation,
see Fig. 6.
calculations shows that the zone had nearly completed formation 25 m-y. ago. Further, it corresponds to the character of the tectonic setting during the Late Cenozoic and to the recent stages of development of the Okhotsk Sea. The tectonic movements have stabilized during this period of time. However, the geothermal and tectonic processes have not yet completely ceased -the formation of the Kuril Basin is the evidence for this. In contrast to the Okhotsk Sea, the Ku~l-K~chatka island arc remained essentially active, in terms of geo~erm~sm, in the Late Cenozoic and Quaternary. The active volcanic processes testify to this activity. The results of modeling on profile X-X show that the partial melting zone has existed for 45 m.y., while on profiles 8-8 and 5-5 the zone
199
Fig. 8. Depth contours to the upper boundary of the partial melting zone 25 m.y. ago. For explanation, see Fig. 6.
demonstrates a 35 m.y.-old history. The thickness of the dioritic layer above the zone is -C5 km in profile X-X, while the total thickness of the granitic and dioritic layers on profiles 8-8 and 5-5 reaches lo-15 km above the partial melting zone. We may conclude that there is a decrease in thickness of the granitic (dioritic) layer over the partial melting zone with time. It should be emphasized that the present calculations are not yet concerned with the problems of the energy which produces the partial melting zone, and the evolution of the zone is examined without any regard to the causes of it. Additionally, the evolution of the zone may cause changes in the location of seismic boundaries and in the Earth’s surface. These problems are also not examined in this paper.
Conclusions
(1) On the basis of presently available heat flow data and data on the structure of the crust and upper mantle, a method of the determination of the evolution of the partial melting zone in the Okhotsk Sea is worked, approximating to a weak non-steadiness. (2) The vast partial melting zone is located by numerical modeling. It began developing 40-45 m.y. ago. The spatial location of the stages of development of the upper boundary of the zone is mapped. (3) A correlation is recorded between the time of existence of the partial melting zone and reduction in thickness of the granitic (dioritic) layers above the zone.
200
Fig. 9. Depth contours to the present upper boundary of the partial melting zone determined by the steady-state calculations.
For
explanation, see Fig. 6.
Bott,
Acknowledgements
M.H.P.,
1981. Crustal
doming and the mechanism
continental rifting. Tectonophysics,
The authors are delighted to thank Dr. V. cerm&k for the efforts made in critically reviewing the manuscript and Dr. T. Nagao for valuable comments. Thanks are also due to N.N. Kovriga, M.G. Shiryaeva, I.E. Avdeeva, A.N. Zagorovich., I.G. Goncharova, S.G. Povyakalo, L.G. Telegina, N.A. Kozhevnikova, Latushkina for
N.P. Lisogurskaya and translation and technical
L.A. assis-
term&k, V., 1982. Geotermischeskaya moshchnosti
litosfery na tenitotii
of
73: l-8. model litosfery i karta SSSR.
Izv. Akad. Nauk
SSSR, Ser. Fiz. Zemli, 1: 25-38. Gliko,
A.O. and Efimov,
A.B.,
1979. 0
razdela faz pri vozmushchenii
dvizhenii granitzy
granichnykh
usloviy. Izv.
Akad. Nauk SSSR, Ser. Fiz. Zemli, 3: 3-14. Gliko, A.O. and Grachev, A.F.,
1987. 0 prirode glubinnykh
processov, obuslavlivayushchikh latformennogo
razvitiye oblastey vnutrip-
magmatizma i kontinentalnykh
riftov. Dokl.
Akad. Nauk SSSR, 295, 1: 64-67. Golubev,
V.A.
and
Osokina,
A.V.,
1980.
Raspredeleniye
teplovogo potoka i priroda ego lokalnykh anomaliy v ra-
tance.
ione ozera Baikal. Izv. Akad. Nauk SSSR, Ser. Fiz. Zemli, 4: 63-75.
References Artyushkov, PP.
Horai, E.V., 1979. Geodinamika.
Nauka, Moscow, 327
K.,
1976.
Heat
flow anomaly
associated
with dike
intrusion. J. Geophys. Res., 81 (5): 894-898. Kutas, K.I. and Tsvyashchenko,
V.A., 1984. Ob ispolzovanii
201
kharakternykh anomaliy. Latters,
tochek
Geofiz.
R. and Lions, Zh.L.,
i ego Prilozheniya. Chislennye Moscow,
Modeli
Moscow, Samarsky,
Poley
O.I., 1983. Nauka,
Raznostnykh
Skhem.
Nauka,
v Chislennye
Metody.
Nauka,
Kartogr.,
Raionov. and Sugrobov,
Vulkanoi.
glubinnykh Seismol.,
Glavnoye
Territorii
Upravleniye
SSSR Geod.
150 pp.
v Kurilo-Kamchatskoi
Otzenka
Okhotskogo
morya.
V.M.,
i Aleutskoi
temperatur 2: 3-18.
1980. Zemnoi
teplovoi
provintziyakh.
i moshcbnost
Moscow,
Resheniya
223 pp.
V.D., 1987. Chislennoye teplovogo
potoka
modlitosfery
Izv. Alcad. Nauk SSSR, Ser. Fiz. Zemh.
G.C. and Tilli, K.A.. 1965. ~ois~o~deniye
3.
litosfery.
V.N.,
Nauka, Zorin,
Moscow,
Ya.B.
statsionarnogo
Zharkov,
283 pp.
i Sopredelnykh
I.K. and Epaneshnikov,
vykh magm.
656 pp.
Potoka
V.Ya., 1974. Metody Nauka,
7: 94-100.
Zemli.
Teplovogo
Zadach.
elirovaniye
Yoder,
1977. Teoriya
Ya.B., 1980. Karta
potok
Tuezov,
Kvaziobrashcheniya
V.M. and Parfenyuk,
A.N. and Arsenin,
Nekorrektnykh
SSR, 6 (4): 72-73.
336 pp.
Teplovykh
A.A., 1987. Vvedeniye
Moscow,
Smimov,
1970. Metod
Tikhonov,
geotermicheskikh
Ukr.
125 pp. A.A.,
Samarsky,
Nauk
Mir, Moscow,
E.A., Luboshitz,
Lubimova,
Smimov,
dlya interpretatzii
J. Akad.
Nauka,
Moscow,
1983. Vnutrenneye
Moscow,
Zemli
i Planet.
S.V., 1985. Geotermischeskye
astenosfemykh
riftovyrni
zonami.
(Editors),
Geotermicheskiye
Kazakhstane.
Stroeniye
415 pp.
Yu. A. and Lepina,
eli razvitiya
bazalto-
247 pp.
vystupov
In: A.V. Shcherbakov
Nauka,
and V.I. Dvorov
Issledovaniya
Moscow,
mod-
pod kontinentalnymi
pp. 187-199.
v Srednei
Azii i