Modeling of the geothermal history of the Okhotsk Sea lithosphere

Modeling of the geothermal history of the Okhotsk Sea lithosphere

Tectonophysics, 164 (1989) 189-201 Elsevier Science Publishers 189 B.V., Amsterdam - Printed in The Netherlands Modeling of the geothermal hist...

985KB Sizes 3 Downloads 135 Views

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