- Email: [email protected]
A geochemical model of mid-ocean ridge magma chambers
Earth and Planetary Scrence Letters, 60 (1982) 93- 104 Elsevier Scientific Publishing Company, Amsterdam - Printed
A geochemical
93 in The Netherlands
model of mid-ocean
ridge magma chambers
Dave Robson and J.R. Cann Department of Geolog)?, The University, Newcastle upon Tyne NE1 7R U (England)
Received November 26. 198 I Revised manuscript received March 29. 1982
A computer model of mid-ocean ridge basalt generation using trace element geochemistry has been developed. The model simulates a periodically replenished, continually cooled and fractionated magma chamber. with periodic lava extrusion. Primitive basalts from the ocean floor are used to generate likely evolution paths for the magma chamber. The steady state variant of this model has led to the isolation of several variables which crittcally affect the basalt composition. Although the fraction of cumulates is an important parameter, other variables such as the volume of incoming magma batches, their frequency, and the volume of the mixing cell. play a critical part especially on slow-spreading ridges. The growing magma chamber model uses random number generators to simulate the initiation and growth of a chamber. This modei predicts a rapid increase in incompatible element concentrations, immediately after chamber initiation on a fast-spreading ridge. This would occur in situations such as propagating rifts and may help in the understanding of ferrobasalt generation.
1. Introduction
In this paper a computer model is presented which simulates the geochemical initiation and evolution of long-lived magma chambers at midocean ridge crests. Residence of basaltic magma within such chambers during its ascent to the ocean floor has recently become an important element in models of mid-ocean ridge basalt generation. The model explores the consequences of such residence in relation to basalt geochemistry. Early workers on ocean floor basalts [l-5] concluded that they could represent primary partial melts from the mantle extracted and erupted essentially unmodified. Geophysical results showing a lack of deep seismicity and presence of high heat flow in this region of continuous upwelling and spreading of mantle, suggest that the asthenosphere comes within a few kilometres of the ocean floor at mid-ocean ridge crests. Therefore eruption of unmodified mantle melts might be expected. O’Hara [6], however, pointed out petrological reasons why oceanic basalts are unlikely to be primary 0012-821X/82/0000-0000/%02.75
0 1982 Elsevier Scientific
melts, and other authors have confirmed and extended his conclusions (e.g. [7-91). The existence of a magma reservoir beneath mid-ocean ridge crests was suggested on theoretical grounds by Cann [ IO,1 I] and from examination of the Troodos complex in Cyprus by Greenbaum [ 121. Geophysical work on the East Pacific rise [ 13- 151 has revealed the existence of a low-velocity zone beneath the ridge axis. This could represent a magma reservoir although any similar structure beneath the slower-spreading FAMOUS area of the Mid-Atlantic Ridge [ 16,171 has escaped detection. The petrography of mid-ocean ridge basalts also points to the operation of magma chamber processes. In many mid-ocean ridge basalts two generations of phenocrysts exist. These are earlier irregular megucrysrs [3,4] of basic plagioclase and forsteritic olivine (Foss-Fo,,) (An ,*- An,,) with later euhedral less basic phenocrysts in equilibrium with the groundmass ([18-201, etc.). It seems likely that the resultant basalt is a mixture of primitive and more evolved melt, the mixing Publishing
Company
94
taking place in a magma chamber undergoing active crystal fractionation. Work on ophiolite complexes, which are now generally regarded as fragments of ancient ocean crust indicate the existence of large magma reservoirs and the operation of magma chamber processes, e.g. Troodos, Cyprus [21], Oman [22]. Such a magma chamber can either exist as a closed system or as a periodically replenished open system reservoir [23], continually accreting ocean crust. Bearing in mind the situation existing at an ocean ridge with continually uprising asthenosphere providing a reasonably plentiful supply of melt, the latter situation appears to be the most likely. The spreading rate will have a fundamental effect on the magma chamber. Thermal modelling [24-263 indicates the existence of a magma chamber beneath ridge crests spreading faster than 0.45-0.9 cm/yr half rate. As spreading rate increases this chamber expands until at a spreading rate of 6.0 cm/yr it will underlie 10 km of ocean floor either side of the spreading ridge [26]. A decrease in the thickness of the sheeted dyke unit as spreading rate increases is also indicated. Below the critical spreading rate (0.45-0.9 cm/yr) no chamber can exist on thermal grounds and since the heat flux is also low the supply of melt will be reduced. In this situation the infinite leak model of Nisbet and Fowler [27] would seem likely, where basaltic liquid fills cracks propagated from the base of the crust. The model presented in this paper is based on the infinite onion model of Cann [ 10,111. This model implies the existence of a high-level magma chamber, from which ocean crust is constantly peeled off, beneath all constructive margins spreading faster than the critical limits.
2. The model: structure The following sections of this paper refer to the structure and results of the model. The assumptions made are those used in the model, based on evidence from ocean ridges. The model is designed to simulate the geochemical evolution of mid-ocean ridge basalts. Evolution of the magma chamber liquids is monitored
by changes in trace element concentrations, rather than using major elements. This allows the use of bulk solid-liquid distribution coefficients (D. Appendix 1, equation ( 1)) to provide a good approximation to the complex crystal-liquid processes. The use of major elements would require a much more cumbersome approach involving phase equilibria. Among the trace elements, incompatible elements prove particularly useful. These are elements for which D is much less than 1, being strongly rejected by crystalline phases present. In ocean floor basalts some incompatible elements such as Zr, Nb, Y are relatively immobile during secondary processes (e.g. [28,29]) and are particularly useful in determining basalt petrogenesis. Complementary compatible elements such as Cr and Ni are also of great value. In the majority of the modelling presented here incompatible elements will be represented by an element with a bulk solid-liquid distribution coefficient of 0.1 (labelled as Zr for convenience) and compatible elements by one with a distribution coefficient of 4 or 5 (labelled as Ni). The structure of ocean crust assumed by the model is that shown in Fig. 1. This is based on geophysical studies on ocean crust (e.g. [30,31]) and on ophiolite complexes (e.g. [32,33]). From the point of view of this geochemical model, the crust can be divided into three units; a lava unit composed of material erupted periodically from the magma chamber, a dyke unit formed from liquid samples derived continuously from the chamber during spreading (sheeted dykes and isotropic gaband a cumulate unit continuously bros), fractionated from the magma chamber (the cumulate gabbros and ultramafics). As mentioned previously the model assumes the existence of a permanent to semi-permanent open system magma chamber beneath all active ocean ridges except those spreading extremely slowly (Fig. 1). The chamber is periodically replenished from below, whilst continually cooling and fractionating to produce oceanic crust. The dykes and cumulates together are assumed to be a constant thickness (usually 5 km) but the lava thickness is allowed to vary for reasons discussed later. The chamber size is directly proportional to the spreading rate [24-261, thus spreading ridges have
95
Magma chamber deflated volume V B Fig. 1. Diagrammatic
representation
Blob of primitive volume VC of a fast-spreading
Lava thickness
TL
Dyke thickness
TD
Cumulate thickness
TC
melt
ocean ridge showing
large continuous chambers while slow-spreading ridges have small discontinuous semi-permanent The material supplying the magma reservoirs. chamber consists of blobs composed of primitive basalt, generated by partial melting of the mantle at depth. The degree of partial melting will have a marked effect on the chemistry of the resultant liquid, but the assumption is made, for the purposes of this modelling, that each blob is produced by a similar degree of partial melting of a homogeneous mantle source region. Although melt production may be continuous, blobs will only rise up buoyantly through the asthenosphere on reaching
_
some of the parameters
referred
to in the text.
a certain minimum size, and therefore blob emplacement into the high level magma chamber is essentially periodic. This is supported by evidence from ophiolite complexes, particularly the cumulates in the Samail ophiolite, Oman [22]. The rate of supply must be proportional to spreading rate to maintain the approximately constant thickness of oceanic crust produced at spreading centres (see Appendix 1, equation (2)). In periods between the addition of blobs the chamber undergoes continual open system surface equilibrium (Rayleigh) fractionation (see Appendix 1, equation (3)) gradually producing liquids
96
more enriched in incompatible elements (i.e. more evolved). As the blob enters, it can either mix with the chamber contents or lie on the chamber floor forming a temporary pool [34]. Since the Rayleigh number is thought to be approximately 10” [26] mixing is assumed to occur rapidly along the length of mixing cell (J in Fig. 1) as active convection will be taking place. If pooling did occur, producing for example olivine cumulate, mixing would still occur after a relatively short time [35]. The effects on the model would be slight, resulting in a small reduction in the effective blob volume to be mixed and a similar reduction in the concentration within the magma chamber liquid of elements compatible with the cumulate phase. The model can correct for such a situation. The blob will mix with the chamber contents according to a simple mass continuity equation (see Appendix 1, equation (4)). The chamber liquids will be evolved (i.e. incompatible element enriched) and mixing of a primitive (i.e. incompatible element depleted) blob with this liquid will result in its composition being reset back to an intermediate composition more depleted in incompatible elements. These chemical changes will be reflected in the dyke unit which is formed continuously with time. In effect the dyke unit will act as a tape recorder for the geochemical evolution of the chamber liquids. The dyke intruded immediately after a blob has entered and mixed will have the most primitive composition for that cycle, and is termed the first dyke composition. After fractionation and immediately before the next blob enters the chamber and mixes, the dyke intruded will have the most evolved composition for that cycle and is termed the lust dyke composition. The magnitude of this sawtooth between first and last dyke compositions will increase as the chamber volume/blob volume ratio decreases, that is as the spreading rate decreases (assuming chamber size is directly proportional to spreading rate and blob volume is constant). Thus fast-spreading ridges should produce dyke units in general with less variation than slow-spreading ridges (Fig. 2). The magnitude of this sawtooth provides dimensional information on some of the parameters involved, and detailed work on ophiolite complex dyke units may extract this information, and pro-
Slow
spreading
ridge
First
dyke
Last
dyke
Time -
Fast
I
spreading
ridge
First dyke
Last
dyke
Time
-
Fig. 2. Variation in dyke composition with time under ideal conditions. Equal-sized blobs are injected at regular intervals into a fully developed magma chamber. The larger sawtooth at slow-spreading rates is related to the small size of the magma chamber, and the lower overall level of Zr to the smaller ratio between cumulates and dykes.
vide a useful indicator of paleospreading rates. The magma chamber is assumed to exist in isostatic equilibrium with the crust above. As the blob enters the chamber, the magma pressure will increase [36] and expand the chamber upwards and outwards upsetting the equilibrium. If the pressure is great enough exceeding the lithostatic pressure and retention characteristics of the chamber roof, then material will be expelled by decompression. This material, in the form of lava, will be extruded until the pressure is such that it can be successfully contained within the chamber [36]. As a result of this on any given ridge a small chamber will erupt less lava than a larger chamber (assuming entry of uniform sized blobs), since the retention pressure will be exceeded by a smaller amount. A chamber in the process of growing will therefore
97
expel less of its contents as lava than would a fully grown chamber. This excess material which is retained will cause further chamber growth. Eventually, the chamber will reach a quasi-steady size, which depends on the rate of cooling and the spreading rate, at which inputs and outputs are, on balance, equal over the long term. In this situation the lava thickness will also reach a steady state at the same time. The extrusion of lava is therefore mainly periodic, as a result of blob input. This is in contrast to the continuous abstraction of dykes and cumulates. As mixing occurs rapidly, then the majority of lava will be first dyke in character. The average lava composition is consequently more primitive than that of the dyke unit. Some material, which was at the very apex of the chamber when the blob entered may be pushed out and extruded without mixing. This will result in the extrusion of some evolved last dyke lavas and give a more primitive aspect to the dyke unit sawtooth. Thus the lava unit will not show all of the compositions developed in the magma chamber, but may show the extremes. The dykes and cumulates are the only units to record the complete geochemical evolution of the chamber.
3. The model magma chamber: typical history The typical history of a chamber from initiation to freeze up could be expected to proceed as follows: The first blob will rise up and occupy the empty area beneath the ridge crest (assuming the previous chamber has dried up completely). The pressure increase will be minimal and only a small amount of lava will be extruded. This lava will be primitive, virtually unmodified mantle melt, and this will be one of the few situations in which such melt can be erupted. The blob, or small magma chamber, as it is now will begin to freeze and crystallize, and the chamber liquids will become’ gradually richer in incompatible elements. Since the chamber will be small, the liquids will become extremely evolved after only a short time. This crystallization will continue either until the chamber dries up completely producing trondhjemite, or until another blob enters. This initial period will be the most critical for the magma chamber. If a
blob enters in time it will mix with the small volume of highly evolved material in the chamber, resetting the composition of the liquid to something between that of the chamber and of the blob. The extent of this reset will be large since the chamber volume/blob volume ratio will be small. The pressure increase will force out lava in greater quantity than previously. Fractionation and crystallization will continue until the next blob enters. This cycle of events will continue indefinitely, unless the supply of blobs to the chamber falters, causing the chamber to shrink and even dry up. This will occur more often on a slow spreading ridge where the chamber will be very small, and the blobs may be infrequent. During freezing up of the chamber, with few blobs entering, and the chamber volume decreasing, progressively less lava will be extruded. Thus although trondhjemitic compositions may exist at depth, very little will be extruded on the surface. If freeze up does not occur the chamber will grow, gradually producing more lava on entry of an equal volume blob. After a finite time the chamber will reach its steady state size. At this point the mean thickness of lava will be extruded, inputs will be balanced by outputs and no further chamber growth will occur under normal circumstances. This can be termed the steady state, and at this time the chamber liquids will have evolved to the point where first dyke and last dyke compositions will be constant for every cycle (assuming constant blob supply). These can be termed the steadv state composition, and it is useful to look at the effect of critical variables upon these.
4. The steady state model: results The results presented here are based on the magma chamber at steady state. The steady state compositions can be calculated directly (see Appendix 1, equations (5) and (6)) without the need to run the model until steady state is achieved. These results would apply to a well established magma chamber with reliable blob supply. As mentioned previously Zr is used on the label for an incompatible element (D = 0.1) and Ni for a compatible element (D = 4-5).
98
duce the difference between last dyke and first dyke Zr concentrations. In Fig. 3 it can also be seen how increasing the dyke unit thickness causes a decrease in the difference between last dyke and first dyke compositions. If as suggested by thermal modelling [24-261 the dyke unit thickness increases as spreading rate decreases. then this effect would counteract that of decreasing the chamber volume/blob volume ratio (i.e. decreasing the spreading rate). An increase in lava thickness has a similar effect but to a different degree, due to the non-continuous production of this unit. Fig.4 shows the variation of Zr with Ni. This incompatible/compatible element plot shows how increasing chamber volume (at constant blob volume of 0.4 km’) causes an increase in the Zr content of the first dykes and a complementary decrease in the Ni content. The last dykes show the opposite, that is a decrease in Zr and an increase in Ni. Increasing the Ni distribution coefficient causes the fields to tend towards the horizontal. If pooling and crystallization of, for example, olivine cumulate, takes place prior to mixing [35] the effects will be as follows. Since Ni is
Fig. 3 shows the effect of increasing the chamber volume/blob volume ratio on the difference between the Zr content of last dykes and first dykes. As the chamber volume/blob volume ratio increases (as spreading rate increases, assuming constant blob volume) the difference between the last dyke and first dyke compositions decreases. After the volume of the chamber exceeds approximately 50 times that of the blob no further observable decrease occurs and the first dykes are equal in composition to the last dykes. These results are similar to those produced by O’Hara and Mathews [37]. Increasing the chamber volume/blob volume ratio leads to a corresponding Zr increase in the mean dyke composition. Thus fast-spreading ridges will have more Zr-rich dyke units on average, but with less variation than would slow-spreading ridges. This may be useful as a palaeospreading rate indicator for use on ophiolite complexes, and detailed studies of the sheeted dyke unit in these areas should reveal information about parameters such as this ratio. The effect of trickling very small blobs into a large chamber almost continuously will also be to re-
CB-
Deflated
chamber
vol./Blob
vol.
Fig. 3. Difference between the concentration of Zr (D =0.1) in the last dykes (C,) and the first dykes (CA) related to changing ratio between deflated chamber volume (V,) and blob volume (V,) at steady state. The crustal thickness used is 6 km. The dotted lines indicate a lava thickness of 0 km while the full lines indicate a lava thickness of 1 km. The numbers refer to the dyke thickness (in km).
99
.^
.__ I””
I”
18
Ni (ppm)
Fig. 4. Relation between Zr (D =0.1) and Ni (D =4.0) at steady state for differing values of deflated chamber volume (I/a) and dyke thickness (To). The blob volume used is 0.4 km3 and the lava thickness 1 km. Dyke thicknesses of 0, 1, 2, 3 and 4 km are shown by the dotted lines. Ve values of 0.1, 1.O and 10 km3 are shown with the full lines. A parental composition of 40 ppm Zr. 250 ppm Ni is assumed.
strongly partitioned into olivine, Ni will be removed from the liquid into the cumulate, and the concentration of Ni in the remaining liquid reduced. The liquid will now precipitate less olivine on crystallization, and therefore the bulk Ni distribution coefficient in the magma chamber will be reduced. The effect of this will be to incline the fields toward the vertical. From the steady state modelling it can be seen how the chamber volume/blob volume ratio, dyke thickness, lava thickness and distribution coefficient are the critical variables in determining the steady state composition.
5. The growing magma chamber model: results In this version of the model a more natural approach is taken. The chamber is allowed to grow from zero to steady state size. Random number generators produce blobs of differing size and at varying intervals around some defined means, and
in this way a more realistic set of results can be produced. More variables need to be introduced including the degree of randomness, length of mixing cell (i.e. the length of ridge crest under which mixing of the magma chamber with the blob takes place and products from this chamber extend), and functions relating lava production to chamber size. The assumption is made that steady state chamber size increases in proportion to spreading rate. Fig. 5 shows the geochemical evolution of the first dykes on a slow-spreading (1 .O cm/yr) ocean ridge, evolving over 1.0 Ma. The random number generators produce an intermittent rate of blob supply, causing freezing up of the chamber to occur. The elemental concentrations fluctuate extensively as the chamber falters in growth and steady state conditions are not necessarily permanent. This would be the situation beneath, for example, the Mid-Atlantic Ridge.
d.
0.2
0.4
TIME (MILLION
0.6
0.8
1.0
YRS)
Fig. 5. Geochemical evolution of a slow-spreading ocean ridge with time from chamber initiation showing the variation in relative enrichment for various elements. The spreading rate used as 1.O cm/yr, and the lines represent elements with distribution coefficients of 0.01, 0.1, 0.5, 1.0, 3.0, 6.0 (from top to bottom). The respective steady state compositions are shown by the arrows.
100
In Fig. 6 the spreading rate has been increased to 3 cm/yr, equivalent to a medium-velocity ridge such as the Galapagos Spreading Centre. A peak in the incompatible element concentrations is developed, immediately following the primitive, incompatible-element-depleted values of chamber initiation. The incompatible element concentrations then fall slowly and eventually reach steady state a significant time after the birth of the chamber. This peak is termed an overshoot peak. A much smaller undershoot peak occurs in the compatible elements. As spreading rate increases the magnitude of the overshoot peak increases. Fig. 7 shows the evolution of a magma chamber in terms of first dyke compositions (i.e. lavas) on a 6 cm/yr ridge (fast-spreading part of the East Pacific Rise). The overshoot peak is now quite large, and steady state valves are not achieved for some 5 Ma. A maximum enrichment over steady state (for an element with D = 0.01) of 40% occurs some
J
I
'O-A.0
0.4
0.2 TIME
-
*/
‘Oo.0 -11
C
I
0.2
0.4 TIME
(MILLION
0.6
0.6
1.0
YRS)
Fig. 6. Geochemical evolution of a medium-spreading ocean ridge with time from chamber initiation, showing the variation in relative enrichment for various elements. The spreading rate used is 3.0 cm/yr and the lines represent elements with distribution coefficients of 0.01,0.1,0.5, 1.0, 3.0 and 6.0 (from top to bottom). The respective steady state compositions are shown by the arrows. Notice the small overshoot peak.
0.6 (MILLION
0.6
1.0
YRS)
Fig. 7. Geochemical evolution of a fast-spreading ocean ridge with time from chamber initiation showing the variation in relative enrichment for various elements. The spreading rate used is 6.0 cm/yr and the lines represent elements with distribution coefficients of 0.01. 0.1. 0.5. 1.0. 3.0 and 6.0 (from top to bottom). The respective steady state compositions are shown by the arrows. Notice the large overshoot peak.
0.013 Ma after chamber initiation. It is therefore possible to produce high-evolved, incompatibleelement-enriched basalts (ferrobasalts) on fastspreading ridges, soon after a new chamber has formed. The steady state basalts on these ridges may be much more depleted in incompatible elements. The effect of increasing the dyke thickness is generally to lower the profile in terms of incompatible element concentrations. Increasing the lava thickness causes a similar lowering but also results in a more pronounced overshoot peak. The overshoot peak is a result of allowing material to remain in the chamber during its early stages of evolution, when the chamber volume is small, and the amount of liquid erupted low. A fraction of liquid is removed as lava, this fraction being proportional to the ratio between the chamber volume and the steady state chamber volume. Most of the liquid remaining from the lava alloca-
101
tion is kept within the chamber to allow the chamber to grow. The magnitude of the overshoot peak depends on the amount of magma retained to build the magma chamber volume. In the model, lava eruption takes place even when the magma chamber is very small. If magma is retained in larger amounts in the chamber to allow it to grow faster, by reducing the output from it or by increasing the input in relation to eruption, then the overshoot peak is reduced in size. and under some conditions may vanish. However the conditions leading to this state of affairs are considered geologically unrealistic.
6. Conclusion Since many of the variables applied in the computer model are taken from estimates of similar variables in the real world, conclusions can be drawn as to the effects of these variables on midocean ridge basalt geochemistry. The steady state modelling shows that a greater range of steady state compositions can exist on slow-spreading ridges than can on fast-spreading ridges. The fast-spreading ridges will have on average more evolved basalts, but with less extremes of composition. The more highly evolved character of fast-spreading ridge basalts has been suggested by many authors (e.g. [38]). Slow-spreading ridges will probably be underlain by small, often semipermanent magma chambers, whereas fast-spreading ridges will have a large permanent magma chamber. Steady state compositions are reached after few input cycles on slow-spreading ridges if blob supply is adequate. More commonly, especially at very slow spreading rates, intermittent blob supply will cause the chamber to falter and even dry up. Spreading must then be taken up by other means such as faulting, until the next chamber establishes itself. The model indicates that ophiolites of slowspreading ridge origin will have a larger proportion of extreme differentiates in their stratigraphy than would fast-spreading ridges, and the dyke unit will be on average more evolved than the lavas. The lavas may also show more highly evolved material at the base of large eruptive units (i.e. last
dyke lava pushed out before mixing). Intense study of these areas may indicate whether mixing is rapid or not. Fast-spreading ridges, allowed to develop fully will eventually achieve steady state compositions. A considerable time may elapse between chamber birth and steady state. The growing magma chamber model shows how the overshoot peak is developed, and how this peak becomes more pronounced as spreading rate increases. Thus on fastspreading ridges many of the basalts sampled may not be steady state at all, but may be overshoot peak basalts, produced by newly formed chambers, associated with ridge crest jumps or propagating rift tips. Work in the Pacific [39,40] shows that similar profiles exist around propagating rift zones with primitive basalts at the very tip followed by ferrobasalts and finally more magnesian basalts (steady state?) behind the active rift. The model presented here would be consistent with a growing magma chamber associated with the propagation of such a rift. There may be no necessity to invoke differing source areas or even varying dyke unit thickness to explain anomalous areas of fast-spreading ocean crust. For example, the Costa Rica Rift (spreading rate 3 cm/yr) is erupting magnesian basalts [41], not dissimilar from parts of the Mid-Atlantic Ridge. Nearby the Galapagos rift is erupting ferrobasalts [41-431, yet the spreading rate is almost identical. The main difference between these two areas are the presence of structures indicating the movement of spreading centres (propagating rifts) in the Galapagos area, whereas the Costa Rica Rift area shows no evidence of such structures. The two areas could be tentatively linked, the Costa Rica Rift basalts representing the steady state compositions, and the ferrobasalts of the Galapagos area, the overshoot peak compositions. Further refinement of the model, especially with relation to specific areas, should prove useful and help unravel the problems of oceanic basalt petrogenesis.
Acknowledgements We would like to thank Andy Adamson and Ros Strens for discussion and critically reading the
102
manuscript. Dave Robson would also like to thank his wife for support and help in drafting of the diagrams. This study was supported by a N.E.R.C. grant.
that the Rayleigh fractionation equation itself does not change, but the distribution coefficient must be expressed differently [44]. In this situation: f=
Ma M,+M,-M,
However,
Appendix 1 -Fundamental
equations
a new distribution
coefficient
must be defined:
D.M,+M,
D’=
Mc+M, Defnitions V J H PC Ro MO MB MA ML MD MC
Then: C,=C-f’D’_”
half spreading rate (cm/yr) mixing cell length (km) thickness of ocean crust (km) density of ocean crust (kg/m3) rate of blob supply mass of blob (kg) mass of chamber before entry of blob (kg) mass of chamber after entry of blob (kg) mass of lava (kg) mass of dykes (kg) mass of cumulates (kg)
and for any trace element
Mixing
CL CS CO CB
G
Steady state equations At steady state the first and last dyke compositions can be calculated as follows. combining equations (3) and (4):
e
bulk solid-liquid distribution coefficient for element e open system Rayleigh distribution coefficient for element e concentration of element e in liquid concentration of element e in solid concentration of element e in blob (primitive) concentration of element e in the magma chamber before mixing (last dyke) concentration of element e in the magma chamber after mixing (first dyke)
Bulk solid-liquid distribution
coefficient
= Co.Mo+(C,.f’D’-‘)).MB
(1)
Considering a single cycle of injection of a blob, extraction lava, fractionation of cumulates and continuous abstraction dykes then:
of of
Blob supply
(2) (Rayleigh fractionation)
A M,+"fo
therefore: Co.Mo
c, =
MB + MO - M,( fcD’-
in an open system
The conventional Rayleigh fractionation equation does not apply in the case of open systems, where liquid fractions are continuously removed (in the form of dykes) as well as cumulate. The removal of the liquid fraction is similar to the abstraction of interstitial liquid in cumulates. It can be shown
“)
(5)
and: C, = CAf(D’-‘)
Appendix 2-The
D = C,/C,
Surfore equilibrium
equation
(4)
c D D’
(3)
(6)
program
The program was written in the FORTRAN IV language and run on the Northumbrian Universities IBM 370/168 mainframe computer under the University of Michigan Terminal operating system. The only special package required is a random number generator. In this particular program the National Algorithms group (NAG) subroutines were used. Input to the program consists of primitive basalt trace element compositions, spreading rate, crustal thickness, lava and dyke thickness, mean blob volume, degree of randomness, length of mixing cell and the elemental distribution coefficients. The random number generators produce log normally distributed numbers, these are restrained by inputting the arithmetic mean and standard deviation of the distribution. Steady state compositions and dimensional parameters are calculated initially. The program calculates the first and last dyke compositions, chamber volume and lava thickness in each cycle. Fig. A-l shows a simplified flow diagram of the ocean crust simulation program.
103
READ
INPUT
PARAMETERS
f INITIALIZE VARIABLES
f -
GENERATE
BLOB
1 CALCULATE TO NEXT
TIME BLOB
I
9
IO 11 12
13
REMOVE
14
LAVA
15
1 FRACTIONATE
CHAMBER
REMOVING AND
DYKES
CUMULATES
16
1 STOP
Fig. A-l. program.
Simplified
flow diagram
17
of ocean
crust
simulation
References 1 A.E.J. Engel and C.G. Engel, Composition of basalts from the Mid-Atlantic Ridge, Science 144 (1964) 1330-1333. 2 A.E.J. Engel, C.G. Engel and R.G. Havens, Chemical characteristics of oceanic basalts and the upper mantle, Geol. Sot. Am. Bull. 76 (1965) 719-734.
18
19
20
I.D. Muir, C.E. Tilley and J.H. Scoon, Basalts from the northern part of the rift zone of the Mid-Atlantic Ridge, J. Petrol. 5 (1964) 409-434. I.D. Muir, C.E. Tilley and J.H. Scoon, Basalts from the northern part of the Mid-Atlantic Ridge, II. The Atlantis collection from near 30°N, J. Petrol. 7 (1966) 193-201. G.D. Nicholls, A.J. Nalwalk and E.E. Hays, The nature and composition of rock samples dredged from the Mid-Atlantic Ridge between 22ON and 52ON, Mar. Geol. 1 (1964) 333343. M.J. O’Hara, Are ocean floor basalts primary magma?, Nature 220 (1968) 683-686. R.W. Kay, N.J. Hubbard and P.W. Gast, Chemical characteristics and origin of ocean ridge volcanic rocks, J. Geophys. Res. 75 (1970) 1585-1613. C.H. Donaldson and R.W. Brown, Refractory megacrysts and magnesium-rich melt inclusions within spine1 in oceanic tholeiites: indicators of magma mixing and parental magma composition, Earth Planet. Sci. Lett. 37 (1977) 81-89. H. Sato. Nickel content of basaltic magmas: identification of primary magmas and a measure of the degree of olivine fractionation, Lithos IO (1977) 113-120. J.R. Cann, New model for the structure of the ocean crust, Nature 226 (1970) 928-930. J.R. Cann, A model for oceanic crustal structure developed, Geophys. J.R. Astron. Sot. 39 (1974) 169-187. D. Greenbaum, Magmatic processes at ocean ridges and evidence from the Troodos Massif, Cyprus, Nature 238 (1972) 18-21. J. Orcutt, B. Kennett, D. LeRoy and W. Prothero, A low velocity zone underlying a fast spreading rise crest, Nature 256 (1975) 475-476. B.R. Rosendahl, Evolution of oceanic crust, 2. Constraints, implications and inferences, J. Geophys. Res. 81 (1976) 5305-5314. B.R. Rosendahl, R.W. Raitt, L.M. Dorman, L.D. Bibee, D.M. Hussong and G.H. Sutton, Evolution of oceanic crust, I. A physical model of the East Pacific Rise crest derived from seismic refraction data, J. Geophys. Res. 81 (1976) 5294-530s. C.M.R. Fowler, The crustal structure of the Mid-Atlantic Ridge crest at 37”N, Geophys. J.R. Astron. Sot. 47 (1976) 459-49 1. C.M.R. Fowler, The Mid-Atlantic Ridge: structure at 45’N, Geophys. J.R. Astron. Sot. 54 (1978) 167-183. D.P. Blanchard, J.M. Rhodes, M.A. Dungan, K.V. Rodgers, C.H. Donaldson, J.C. Brannon, J.W. Jacobs and E.K. Gibson, The chemistry and petrology of basalts from Leg 37 of the Deep Sea Drilling Project, J. Geophys. Res. 81 (1976) 423 I-4246. M.A. Dungan and J.M. Rhodes, Residual glasses and melt inclusions in basalts from D.S.D.P. Legs 45 and 46: evidence for magma mixing, Contrib. Mineral. Petrol. 67 (1978) 417-431. M.F.J. Flower, W. Ohnmacht, H.-U. Schmincke, I.L. Gibson, P.T. Robinson and R. Parker, Petrology and geochemistry of basalts from Hole 396B, Leg 46, in L. Dmitriev, J.
104
21
22
23
24 25
26
27
28 29
30
31
32 33
34
Heirtzler et al., Initial Reports of the Deep Sea Drilling Project, 46 (U.S. Government Printing Office, Washington, D.C., 1978) 179-213. CR. Allen, The petrology of a portion of the Troodos plutonic complex, Cyprus Ph.D. Thesis, University of Cambridge (1975) unpublished. J.D. Smewing, Mixing characteristics and compositional differences in mantle-derived melts beneath spreading axes: evidence from cyclically layered rocks in the ophiolite of North Oman, J. Geophys. Res. 86 (1981) 2645-2659. M.J. O’Hara, Geochemical evolution during fractional crystallization of a periodically refilled magma chamber, Nature 266 (1977) 503-507. N.H. Sleep, Formation of oceanic crust, some thermal constraints, J. Geophys. Res. 80 (1975) 4037-4042. N.J. Kusznir and M.H.P. Bott, A thermal study of the formation of oceanic crust, Geophys. J.R. Astron. Sot. 47 (1976) 83-95. N.J. Kusznir, Thermal evolution of the oceanic crust; its dependence on spreading rate and effect on crustal structure, Geophys. J.R. Astron. Sot. 61 (1980) 167-181. E.G. Nisbet and C.M.R. Fowler, The Mid-Atlantic Ridge at 37 and 45’N: some geophysical and petrological constraints, Geophys. J.R. Astron. Sot. 54 (1978) 631-660. J.R. Cann, Rb, Sr, Y, Zr and Nb in some ocean floor basaltic rocks, Earth Planet. Sci. Lett. 10 (1970) 7-l 1. D.A. Wood, I.L. Gibson and R.N. Thompson, Elemental mobility during zeohte facies metamorphism of the Tertiary basalts of eastern Iceland, Contrib. Mineral. Petrol. 55 (1976) 241-254. W.J. Ludwig, J.E. Nafe and CL. Drake, Seismic refraction, The Sea, Vol. 4, A.E. Maxwell, ed. (Wiley-Interscience, New York, N.Y., 1970) 53-84. G.L. Maynard, Crustal layer of seismic velocity 6.9-7.6 kilometres per second under the deep oceans, Science 168 (1970) 120-121. E.H. Moores and E.D. Jackson, Ophiolites and oceanic crust, Nature 250 (1974) 136-139. M.H. Salisbury and N.I. Christensen, The seismic velocity structure of a traverse through the Bay of Islands ophiolite complex, Newfoundland, an exposure of oceanic crust and upper mantle, J. Geophys. Res. 83 (1978) 805-817. R.S.J. Sparks, P. Meyer and H. Sigurdsson, Density varia-
35
36 37
38 39
40
41
42
43
44
tion amongst mid-ocean ridge basalts: implications for magma mixing and the scarcity of primitive lavas, Earth Planet. Sci. Lett. 46 (1980) 419-430. H.E. Huppert and R.S.J. Sparks. Restrictions on the composition of mid-ocean ridge basalts: a fluid dynamical investigation. Nature 286 (1980) 46-48. S. Blake, Volcanism and the dynamics of open magma chambers, Nature 289 (1981) 783-785. M.J. G’Hara and R.E. Mathews, Geochemical evolution in an advancing, periodically replenished, periodically tapped, continuously fractionated magma chamber, J. Geol. Sot. London 138 (1981) 237-277. E.G. Nisbet and J.A. Pearce, TiO, and a possible guide to past oceanic spreading rates, Nature 246 (1973) 468-469. J.M. Sinton, D.M. Christie, D.S. Wilson. R. Hey, Petrological consequences of rift propagation on oceanic spreading ridges, Abstr., Am. Geophys. Union Chapman Conf. on The Generation of the Oceanic Lithosphere, Virginia, U.S.A.. April (1981). J.R. Delaney, H.P. Johnson, J.L. Karsten, The Juan de Fuca ridge hot spot-propagating rift system: new tectonic, geochemical and magnetic data, J. Geophys. Res. 86 (198 1) 747-750. R.N. Anderson, D.A. Clague, K.D. Klitgord, M. Marshall and R.K. Nishimori, Magnetic and petrological variations along the Galapagos spreading centre and their relation to the Galapagos melting anomaly, Geol. Sot. Am. Bull. 86 (1975) 683-694. D.P. Mattey and I.D. Muir, Geochemistry and mineralogy of basalts from the Galapagos spreading centre, Deep Sea Drilling Project, Leg 54, in: B.R. Rosendahl, R. Hekinian et al.. Initial Reports of the Deep Sea Drilling Project, 54 (U.S. Government Printing Office, Washington, D.C.. 1980) 755-770. J.H. Natland and W.G. Melson, Composition of basaltic glasses from near the East Pacific Rise and Siqueros fracture zone, near 9”N, in: B.R. Rosendahl, R. Hekinian et al., Initial Reports of the Deep Sea Drilling Project, 54 (U.S. Government Printing Office, Washington, D.C., 1980) 705723. J.R. Cann, Rayleigh fractionation with continuous removal of liquid, Earth Planet. Sci. Lett. 60 (1982) 114-I 16 (this issue).
A geochemical
93 in The Netherlands
model of mid-ocean
ridge magma chambers
Dave Robson and J.R. Cann Department of Geolog)?, The University, Newcastle upon Tyne NE1 7R U (England)
Received November 26. 198 I Revised manuscript received March 29. 1982
A computer model of mid-ocean ridge basalt generation using trace element geochemistry has been developed. The model simulates a periodically replenished, continually cooled and fractionated magma chamber. with periodic lava extrusion. Primitive basalts from the ocean floor are used to generate likely evolution paths for the magma chamber. The steady state variant of this model has led to the isolation of several variables which crittcally affect the basalt composition. Although the fraction of cumulates is an important parameter, other variables such as the volume of incoming magma batches, their frequency, and the volume of the mixing cell. play a critical part especially on slow-spreading ridges. The growing magma chamber model uses random number generators to simulate the initiation and growth of a chamber. This modei predicts a rapid increase in incompatible element concentrations, immediately after chamber initiation on a fast-spreading ridge. This would occur in situations such as propagating rifts and may help in the understanding of ferrobasalt generation.
1. Introduction
In this paper a computer model is presented which simulates the geochemical initiation and evolution of long-lived magma chambers at midocean ridge crests. Residence of basaltic magma within such chambers during its ascent to the ocean floor has recently become an important element in models of mid-ocean ridge basalt generation. The model explores the consequences of such residence in relation to basalt geochemistry. Early workers on ocean floor basalts [l-5] concluded that they could represent primary partial melts from the mantle extracted and erupted essentially unmodified. Geophysical results showing a lack of deep seismicity and presence of high heat flow in this region of continuous upwelling and spreading of mantle, suggest that the asthenosphere comes within a few kilometres of the ocean floor at mid-ocean ridge crests. Therefore eruption of unmodified mantle melts might be expected. O’Hara [6], however, pointed out petrological reasons why oceanic basalts are unlikely to be primary 0012-821X/82/0000-0000/%02.75
0 1982 Elsevier Scientific
melts, and other authors have confirmed and extended his conclusions (e.g. [7-91). The existence of a magma reservoir beneath mid-ocean ridge crests was suggested on theoretical grounds by Cann [ IO,1 I] and from examination of the Troodos complex in Cyprus by Greenbaum [ 121. Geophysical work on the East Pacific rise [ 13- 151 has revealed the existence of a low-velocity zone beneath the ridge axis. This could represent a magma reservoir although any similar structure beneath the slower-spreading FAMOUS area of the Mid-Atlantic Ridge [ 16,171 has escaped detection. The petrography of mid-ocean ridge basalts also points to the operation of magma chamber processes. In many mid-ocean ridge basalts two generations of phenocrysts exist. These are earlier irregular megucrysrs [3,4] of basic plagioclase and forsteritic olivine (Foss-Fo,,) (An ,*- An,,) with later euhedral less basic phenocrysts in equilibrium with the groundmass ([18-201, etc.). It seems likely that the resultant basalt is a mixture of primitive and more evolved melt, the mixing Publishing
Company
94
taking place in a magma chamber undergoing active crystal fractionation. Work on ophiolite complexes, which are now generally regarded as fragments of ancient ocean crust indicate the existence of large magma reservoirs and the operation of magma chamber processes, e.g. Troodos, Cyprus [21], Oman [22]. Such a magma chamber can either exist as a closed system or as a periodically replenished open system reservoir [23], continually accreting ocean crust. Bearing in mind the situation existing at an ocean ridge with continually uprising asthenosphere providing a reasonably plentiful supply of melt, the latter situation appears to be the most likely. The spreading rate will have a fundamental effect on the magma chamber. Thermal modelling [24-263 indicates the existence of a magma chamber beneath ridge crests spreading faster than 0.45-0.9 cm/yr half rate. As spreading rate increases this chamber expands until at a spreading rate of 6.0 cm/yr it will underlie 10 km of ocean floor either side of the spreading ridge [26]. A decrease in the thickness of the sheeted dyke unit as spreading rate increases is also indicated. Below the critical spreading rate (0.45-0.9 cm/yr) no chamber can exist on thermal grounds and since the heat flux is also low the supply of melt will be reduced. In this situation the infinite leak model of Nisbet and Fowler [27] would seem likely, where basaltic liquid fills cracks propagated from the base of the crust. The model presented in this paper is based on the infinite onion model of Cann [ 10,111. This model implies the existence of a high-level magma chamber, from which ocean crust is constantly peeled off, beneath all constructive margins spreading faster than the critical limits.
2. The model: structure The following sections of this paper refer to the structure and results of the model. The assumptions made are those used in the model, based on evidence from ocean ridges. The model is designed to simulate the geochemical evolution of mid-ocean ridge basalts. Evolution of the magma chamber liquids is monitored
by changes in trace element concentrations, rather than using major elements. This allows the use of bulk solid-liquid distribution coefficients (D. Appendix 1, equation ( 1)) to provide a good approximation to the complex crystal-liquid processes. The use of major elements would require a much more cumbersome approach involving phase equilibria. Among the trace elements, incompatible elements prove particularly useful. These are elements for which D is much less than 1, being strongly rejected by crystalline phases present. In ocean floor basalts some incompatible elements such as Zr, Nb, Y are relatively immobile during secondary processes (e.g. [28,29]) and are particularly useful in determining basalt petrogenesis. Complementary compatible elements such as Cr and Ni are also of great value. In the majority of the modelling presented here incompatible elements will be represented by an element with a bulk solid-liquid distribution coefficient of 0.1 (labelled as Zr for convenience) and compatible elements by one with a distribution coefficient of 4 or 5 (labelled as Ni). The structure of ocean crust assumed by the model is that shown in Fig. 1. This is based on geophysical studies on ocean crust (e.g. [30,31]) and on ophiolite complexes (e.g. [32,33]). From the point of view of this geochemical model, the crust can be divided into three units; a lava unit composed of material erupted periodically from the magma chamber, a dyke unit formed from liquid samples derived continuously from the chamber during spreading (sheeted dykes and isotropic gaband a cumulate unit continuously bros), fractionated from the magma chamber (the cumulate gabbros and ultramafics). As mentioned previously the model assumes the existence of a permanent to semi-permanent open system magma chamber beneath all active ocean ridges except those spreading extremely slowly (Fig. 1). The chamber is periodically replenished from below, whilst continually cooling and fractionating to produce oceanic crust. The dykes and cumulates together are assumed to be a constant thickness (usually 5 km) but the lava thickness is allowed to vary for reasons discussed later. The chamber size is directly proportional to the spreading rate [24-261, thus spreading ridges have
95
Magma chamber deflated volume V B Fig. 1. Diagrammatic
representation
Blob of primitive volume VC of a fast-spreading
Lava thickness
TL
Dyke thickness
TD
Cumulate thickness
TC
melt
ocean ridge showing
large continuous chambers while slow-spreading ridges have small discontinuous semi-permanent The material supplying the magma reservoirs. chamber consists of blobs composed of primitive basalt, generated by partial melting of the mantle at depth. The degree of partial melting will have a marked effect on the chemistry of the resultant liquid, but the assumption is made, for the purposes of this modelling, that each blob is produced by a similar degree of partial melting of a homogeneous mantle source region. Although melt production may be continuous, blobs will only rise up buoyantly through the asthenosphere on reaching
_
some of the parameters
referred
to in the text.
a certain minimum size, and therefore blob emplacement into the high level magma chamber is essentially periodic. This is supported by evidence from ophiolite complexes, particularly the cumulates in the Samail ophiolite, Oman [22]. The rate of supply must be proportional to spreading rate to maintain the approximately constant thickness of oceanic crust produced at spreading centres (see Appendix 1, equation (2)). In periods between the addition of blobs the chamber undergoes continual open system surface equilibrium (Rayleigh) fractionation (see Appendix 1, equation (3)) gradually producing liquids
96
more enriched in incompatible elements (i.e. more evolved). As the blob enters, it can either mix with the chamber contents or lie on the chamber floor forming a temporary pool [34]. Since the Rayleigh number is thought to be approximately 10” [26] mixing is assumed to occur rapidly along the length of mixing cell (J in Fig. 1) as active convection will be taking place. If pooling did occur, producing for example olivine cumulate, mixing would still occur after a relatively short time [35]. The effects on the model would be slight, resulting in a small reduction in the effective blob volume to be mixed and a similar reduction in the concentration within the magma chamber liquid of elements compatible with the cumulate phase. The model can correct for such a situation. The blob will mix with the chamber contents according to a simple mass continuity equation (see Appendix 1, equation (4)). The chamber liquids will be evolved (i.e. incompatible element enriched) and mixing of a primitive (i.e. incompatible element depleted) blob with this liquid will result in its composition being reset back to an intermediate composition more depleted in incompatible elements. These chemical changes will be reflected in the dyke unit which is formed continuously with time. In effect the dyke unit will act as a tape recorder for the geochemical evolution of the chamber liquids. The dyke intruded immediately after a blob has entered and mixed will have the most primitive composition for that cycle, and is termed the first dyke composition. After fractionation and immediately before the next blob enters the chamber and mixes, the dyke intruded will have the most evolved composition for that cycle and is termed the lust dyke composition. The magnitude of this sawtooth between first and last dyke compositions will increase as the chamber volume/blob volume ratio decreases, that is as the spreading rate decreases (assuming chamber size is directly proportional to spreading rate and blob volume is constant). Thus fast-spreading ridges should produce dyke units in general with less variation than slow-spreading ridges (Fig. 2). The magnitude of this sawtooth provides dimensional information on some of the parameters involved, and detailed work on ophiolite complex dyke units may extract this information, and pro-
Slow
spreading
ridge
First
dyke
Last
dyke
Time -
Fast
I
spreading
ridge
First dyke
Last
dyke
Time
-
Fig. 2. Variation in dyke composition with time under ideal conditions. Equal-sized blobs are injected at regular intervals into a fully developed magma chamber. The larger sawtooth at slow-spreading rates is related to the small size of the magma chamber, and the lower overall level of Zr to the smaller ratio between cumulates and dykes.
vide a useful indicator of paleospreading rates. The magma chamber is assumed to exist in isostatic equilibrium with the crust above. As the blob enters the chamber, the magma pressure will increase [36] and expand the chamber upwards and outwards upsetting the equilibrium. If the pressure is great enough exceeding the lithostatic pressure and retention characteristics of the chamber roof, then material will be expelled by decompression. This material, in the form of lava, will be extruded until the pressure is such that it can be successfully contained within the chamber [36]. As a result of this on any given ridge a small chamber will erupt less lava than a larger chamber (assuming entry of uniform sized blobs), since the retention pressure will be exceeded by a smaller amount. A chamber in the process of growing will therefore
97
expel less of its contents as lava than would a fully grown chamber. This excess material which is retained will cause further chamber growth. Eventually, the chamber will reach a quasi-steady size, which depends on the rate of cooling and the spreading rate, at which inputs and outputs are, on balance, equal over the long term. In this situation the lava thickness will also reach a steady state at the same time. The extrusion of lava is therefore mainly periodic, as a result of blob input. This is in contrast to the continuous abstraction of dykes and cumulates. As mixing occurs rapidly, then the majority of lava will be first dyke in character. The average lava composition is consequently more primitive than that of the dyke unit. Some material, which was at the very apex of the chamber when the blob entered may be pushed out and extruded without mixing. This will result in the extrusion of some evolved last dyke lavas and give a more primitive aspect to the dyke unit sawtooth. Thus the lava unit will not show all of the compositions developed in the magma chamber, but may show the extremes. The dykes and cumulates are the only units to record the complete geochemical evolution of the chamber.
3. The model magma chamber: typical history The typical history of a chamber from initiation to freeze up could be expected to proceed as follows: The first blob will rise up and occupy the empty area beneath the ridge crest (assuming the previous chamber has dried up completely). The pressure increase will be minimal and only a small amount of lava will be extruded. This lava will be primitive, virtually unmodified mantle melt, and this will be one of the few situations in which such melt can be erupted. The blob, or small magma chamber, as it is now will begin to freeze and crystallize, and the chamber liquids will become’ gradually richer in incompatible elements. Since the chamber will be small, the liquids will become extremely evolved after only a short time. This crystallization will continue either until the chamber dries up completely producing trondhjemite, or until another blob enters. This initial period will be the most critical for the magma chamber. If a
blob enters in time it will mix with the small volume of highly evolved material in the chamber, resetting the composition of the liquid to something between that of the chamber and of the blob. The extent of this reset will be large since the chamber volume/blob volume ratio will be small. The pressure increase will force out lava in greater quantity than previously. Fractionation and crystallization will continue until the next blob enters. This cycle of events will continue indefinitely, unless the supply of blobs to the chamber falters, causing the chamber to shrink and even dry up. This will occur more often on a slow spreading ridge where the chamber will be very small, and the blobs may be infrequent. During freezing up of the chamber, with few blobs entering, and the chamber volume decreasing, progressively less lava will be extruded. Thus although trondhjemitic compositions may exist at depth, very little will be extruded on the surface. If freeze up does not occur the chamber will grow, gradually producing more lava on entry of an equal volume blob. After a finite time the chamber will reach its steady state size. At this point the mean thickness of lava will be extruded, inputs will be balanced by outputs and no further chamber growth will occur under normal circumstances. This can be termed the steady state, and at this time the chamber liquids will have evolved to the point where first dyke and last dyke compositions will be constant for every cycle (assuming constant blob supply). These can be termed the steadv state composition, and it is useful to look at the effect of critical variables upon these.
4. The steady state model: results The results presented here are based on the magma chamber at steady state. The steady state compositions can be calculated directly (see Appendix 1, equations (5) and (6)) without the need to run the model until steady state is achieved. These results would apply to a well established magma chamber with reliable blob supply. As mentioned previously Zr is used on the label for an incompatible element (D = 0.1) and Ni for a compatible element (D = 4-5).
98
duce the difference between last dyke and first dyke Zr concentrations. In Fig. 3 it can also be seen how increasing the dyke unit thickness causes a decrease in the difference between last dyke and first dyke compositions. If as suggested by thermal modelling [24-261 the dyke unit thickness increases as spreading rate decreases. then this effect would counteract that of decreasing the chamber volume/blob volume ratio (i.e. decreasing the spreading rate). An increase in lava thickness has a similar effect but to a different degree, due to the non-continuous production of this unit. Fig.4 shows the variation of Zr with Ni. This incompatible/compatible element plot shows how increasing chamber volume (at constant blob volume of 0.4 km’) causes an increase in the Zr content of the first dykes and a complementary decrease in the Ni content. The last dykes show the opposite, that is a decrease in Zr and an increase in Ni. Increasing the Ni distribution coefficient causes the fields to tend towards the horizontal. If pooling and crystallization of, for example, olivine cumulate, takes place prior to mixing [35] the effects will be as follows. Since Ni is
Fig. 3 shows the effect of increasing the chamber volume/blob volume ratio on the difference between the Zr content of last dykes and first dykes. As the chamber volume/blob volume ratio increases (as spreading rate increases, assuming constant blob volume) the difference between the last dyke and first dyke compositions decreases. After the volume of the chamber exceeds approximately 50 times that of the blob no further observable decrease occurs and the first dykes are equal in composition to the last dykes. These results are similar to those produced by O’Hara and Mathews [37]. Increasing the chamber volume/blob volume ratio leads to a corresponding Zr increase in the mean dyke composition. Thus fast-spreading ridges will have more Zr-rich dyke units on average, but with less variation than would slow-spreading ridges. This may be useful as a palaeospreading rate indicator for use on ophiolite complexes, and detailed studies of the sheeted dyke unit in these areas should reveal information about parameters such as this ratio. The effect of trickling very small blobs into a large chamber almost continuously will also be to re-
CB-
Deflated
chamber
vol./Blob
vol.
Fig. 3. Difference between the concentration of Zr (D =0.1) in the last dykes (C,) and the first dykes (CA) related to changing ratio between deflated chamber volume (V,) and blob volume (V,) at steady state. The crustal thickness used is 6 km. The dotted lines indicate a lava thickness of 0 km while the full lines indicate a lava thickness of 1 km. The numbers refer to the dyke thickness (in km).
99
.^
.__ I””
I”
18
Ni (ppm)
Fig. 4. Relation between Zr (D =0.1) and Ni (D =4.0) at steady state for differing values of deflated chamber volume (I/a) and dyke thickness (To). The blob volume used is 0.4 km3 and the lava thickness 1 km. Dyke thicknesses of 0, 1, 2, 3 and 4 km are shown by the dotted lines. Ve values of 0.1, 1.O and 10 km3 are shown with the full lines. A parental composition of 40 ppm Zr. 250 ppm Ni is assumed.
strongly partitioned into olivine, Ni will be removed from the liquid into the cumulate, and the concentration of Ni in the remaining liquid reduced. The liquid will now precipitate less olivine on crystallization, and therefore the bulk Ni distribution coefficient in the magma chamber will be reduced. The effect of this will be to incline the fields toward the vertical. From the steady state modelling it can be seen how the chamber volume/blob volume ratio, dyke thickness, lava thickness and distribution coefficient are the critical variables in determining the steady state composition.
5. The growing magma chamber model: results In this version of the model a more natural approach is taken. The chamber is allowed to grow from zero to steady state size. Random number generators produce blobs of differing size and at varying intervals around some defined means, and
in this way a more realistic set of results can be produced. More variables need to be introduced including the degree of randomness, length of mixing cell (i.e. the length of ridge crest under which mixing of the magma chamber with the blob takes place and products from this chamber extend), and functions relating lava production to chamber size. The assumption is made that steady state chamber size increases in proportion to spreading rate. Fig. 5 shows the geochemical evolution of the first dykes on a slow-spreading (1 .O cm/yr) ocean ridge, evolving over 1.0 Ma. The random number generators produce an intermittent rate of blob supply, causing freezing up of the chamber to occur. The elemental concentrations fluctuate extensively as the chamber falters in growth and steady state conditions are not necessarily permanent. This would be the situation beneath, for example, the Mid-Atlantic Ridge.
d.
0.2
0.4
TIME (MILLION
0.6
0.8
1.0
YRS)
Fig. 5. Geochemical evolution of a slow-spreading ocean ridge with time from chamber initiation showing the variation in relative enrichment for various elements. The spreading rate used as 1.O cm/yr, and the lines represent elements with distribution coefficients of 0.01, 0.1, 0.5, 1.0, 3.0, 6.0 (from top to bottom). The respective steady state compositions are shown by the arrows.
100
In Fig. 6 the spreading rate has been increased to 3 cm/yr, equivalent to a medium-velocity ridge such as the Galapagos Spreading Centre. A peak in the incompatible element concentrations is developed, immediately following the primitive, incompatible-element-depleted values of chamber initiation. The incompatible element concentrations then fall slowly and eventually reach steady state a significant time after the birth of the chamber. This peak is termed an overshoot peak. A much smaller undershoot peak occurs in the compatible elements. As spreading rate increases the magnitude of the overshoot peak increases. Fig. 7 shows the evolution of a magma chamber in terms of first dyke compositions (i.e. lavas) on a 6 cm/yr ridge (fast-spreading part of the East Pacific Rise). The overshoot peak is now quite large, and steady state valves are not achieved for some 5 Ma. A maximum enrichment over steady state (for an element with D = 0.01) of 40% occurs some
J
I
'O-A.0
0.4
0.2 TIME
-
*/
‘Oo.0 -11
C
I
0.2
0.4 TIME
(MILLION
0.6
0.6
1.0
YRS)
Fig. 6. Geochemical evolution of a medium-spreading ocean ridge with time from chamber initiation, showing the variation in relative enrichment for various elements. The spreading rate used is 3.0 cm/yr and the lines represent elements with distribution coefficients of 0.01,0.1,0.5, 1.0, 3.0 and 6.0 (from top to bottom). The respective steady state compositions are shown by the arrows. Notice the small overshoot peak.
0.6 (MILLION
0.6
1.0
YRS)
Fig. 7. Geochemical evolution of a fast-spreading ocean ridge with time from chamber initiation showing the variation in relative enrichment for various elements. The spreading rate used is 6.0 cm/yr and the lines represent elements with distribution coefficients of 0.01. 0.1. 0.5. 1.0. 3.0 and 6.0 (from top to bottom). The respective steady state compositions are shown by the arrows. Notice the large overshoot peak.
0.013 Ma after chamber initiation. It is therefore possible to produce high-evolved, incompatibleelement-enriched basalts (ferrobasalts) on fastspreading ridges, soon after a new chamber has formed. The steady state basalts on these ridges may be much more depleted in incompatible elements. The effect of increasing the dyke thickness is generally to lower the profile in terms of incompatible element concentrations. Increasing the lava thickness causes a similar lowering but also results in a more pronounced overshoot peak. The overshoot peak is a result of allowing material to remain in the chamber during its early stages of evolution, when the chamber volume is small, and the amount of liquid erupted low. A fraction of liquid is removed as lava, this fraction being proportional to the ratio between the chamber volume and the steady state chamber volume. Most of the liquid remaining from the lava alloca-
101
tion is kept within the chamber to allow the chamber to grow. The magnitude of the overshoot peak depends on the amount of magma retained to build the magma chamber volume. In the model, lava eruption takes place even when the magma chamber is very small. If magma is retained in larger amounts in the chamber to allow it to grow faster, by reducing the output from it or by increasing the input in relation to eruption, then the overshoot peak is reduced in size. and under some conditions may vanish. However the conditions leading to this state of affairs are considered geologically unrealistic.
6. Conclusion Since many of the variables applied in the computer model are taken from estimates of similar variables in the real world, conclusions can be drawn as to the effects of these variables on midocean ridge basalt geochemistry. The steady state modelling shows that a greater range of steady state compositions can exist on slow-spreading ridges than can on fast-spreading ridges. The fast-spreading ridges will have on average more evolved basalts, but with less extremes of composition. The more highly evolved character of fast-spreading ridge basalts has been suggested by many authors (e.g. [38]). Slow-spreading ridges will probably be underlain by small, often semipermanent magma chambers, whereas fast-spreading ridges will have a large permanent magma chamber. Steady state compositions are reached after few input cycles on slow-spreading ridges if blob supply is adequate. More commonly, especially at very slow spreading rates, intermittent blob supply will cause the chamber to falter and even dry up. Spreading must then be taken up by other means such as faulting, until the next chamber establishes itself. The model indicates that ophiolites of slowspreading ridge origin will have a larger proportion of extreme differentiates in their stratigraphy than would fast-spreading ridges, and the dyke unit will be on average more evolved than the lavas. The lavas may also show more highly evolved material at the base of large eruptive units (i.e. last
dyke lava pushed out before mixing). Intense study of these areas may indicate whether mixing is rapid or not. Fast-spreading ridges, allowed to develop fully will eventually achieve steady state compositions. A considerable time may elapse between chamber birth and steady state. The growing magma chamber model shows how the overshoot peak is developed, and how this peak becomes more pronounced as spreading rate increases. Thus on fastspreading ridges many of the basalts sampled may not be steady state at all, but may be overshoot peak basalts, produced by newly formed chambers, associated with ridge crest jumps or propagating rift tips. Work in the Pacific [39,40] shows that similar profiles exist around propagating rift zones with primitive basalts at the very tip followed by ferrobasalts and finally more magnesian basalts (steady state?) behind the active rift. The model presented here would be consistent with a growing magma chamber associated with the propagation of such a rift. There may be no necessity to invoke differing source areas or even varying dyke unit thickness to explain anomalous areas of fast-spreading ocean crust. For example, the Costa Rica Rift (spreading rate 3 cm/yr) is erupting magnesian basalts [41], not dissimilar from parts of the Mid-Atlantic Ridge. Nearby the Galapagos rift is erupting ferrobasalts [41-431, yet the spreading rate is almost identical. The main difference between these two areas are the presence of structures indicating the movement of spreading centres (propagating rifts) in the Galapagos area, whereas the Costa Rica Rift area shows no evidence of such structures. The two areas could be tentatively linked, the Costa Rica Rift basalts representing the steady state compositions, and the ferrobasalts of the Galapagos area, the overshoot peak compositions. Further refinement of the model, especially with relation to specific areas, should prove useful and help unravel the problems of oceanic basalt petrogenesis.
Acknowledgements We would like to thank Andy Adamson and Ros Strens for discussion and critically reading the
102
manuscript. Dave Robson would also like to thank his wife for support and help in drafting of the diagrams. This study was supported by a N.E.R.C. grant.
that the Rayleigh fractionation equation itself does not change, but the distribution coefficient must be expressed differently [44]. In this situation: f=
Ma M,+M,-M,
However,
Appendix 1 -Fundamental
equations
a new distribution
coefficient
must be defined:
D.M,+M,
D’=
Mc+M, Defnitions V J H PC Ro MO MB MA ML MD MC
Then: C,=C-f’D’_”
half spreading rate (cm/yr) mixing cell length (km) thickness of ocean crust (km) density of ocean crust (kg/m3) rate of blob supply mass of blob (kg) mass of chamber before entry of blob (kg) mass of chamber after entry of blob (kg) mass of lava (kg) mass of dykes (kg) mass of cumulates (kg)
and for any trace element
Mixing
CL CS CO CB
G
Steady state equations At steady state the first and last dyke compositions can be calculated as follows. combining equations (3) and (4):
e
bulk solid-liquid distribution coefficient for element e open system Rayleigh distribution coefficient for element e concentration of element e in liquid concentration of element e in solid concentration of element e in blob (primitive) concentration of element e in the magma chamber before mixing (last dyke) concentration of element e in the magma chamber after mixing (first dyke)
Bulk solid-liquid distribution
coefficient
= Co.Mo+(C,.f’D’-‘)).MB
(1)
Considering a single cycle of injection of a blob, extraction lava, fractionation of cumulates and continuous abstraction dykes then:
of of
Blob supply
(2) (Rayleigh fractionation)
A M,+"fo
therefore: Co.Mo
c, =
MB + MO - M,( fcD’-
in an open system
The conventional Rayleigh fractionation equation does not apply in the case of open systems, where liquid fractions are continuously removed (in the form of dykes) as well as cumulate. The removal of the liquid fraction is similar to the abstraction of interstitial liquid in cumulates. It can be shown
“)
(5)
and: C, = CAf(D’-‘)
Appendix 2-The
D = C,/C,
Surfore equilibrium
equation
(4)
c D D’
(3)
(6)
program
The program was written in the FORTRAN IV language and run on the Northumbrian Universities IBM 370/168 mainframe computer under the University of Michigan Terminal operating system. The only special package required is a random number generator. In this particular program the National Algorithms group (NAG) subroutines were used. Input to the program consists of primitive basalt trace element compositions, spreading rate, crustal thickness, lava and dyke thickness, mean blob volume, degree of randomness, length of mixing cell and the elemental distribution coefficients. The random number generators produce log normally distributed numbers, these are restrained by inputting the arithmetic mean and standard deviation of the distribution. Steady state compositions and dimensional parameters are calculated initially. The program calculates the first and last dyke compositions, chamber volume and lava thickness in each cycle. Fig. A-l shows a simplified flow diagram of the ocean crust simulation program.
103
READ
INPUT
PARAMETERS
f INITIALIZE VARIABLES
f -
GENERATE
BLOB
1 CALCULATE TO NEXT
TIME BLOB
I
9
IO 11 12
13
REMOVE
14
LAVA
15
1 FRACTIONATE
CHAMBER
REMOVING AND
DYKES
CUMULATES
16
1 STOP
Fig. A-l. program.
Simplified
flow diagram
17
of ocean
crust
simulation
References 1 A.E.J. Engel and C.G. Engel, Composition of basalts from the Mid-Atlantic Ridge, Science 144 (1964) 1330-1333. 2 A.E.J. Engel, C.G. Engel and R.G. Havens, Chemical characteristics of oceanic basalts and the upper mantle, Geol. Sot. Am. Bull. 76 (1965) 719-734.
18
19
20
I.D. Muir, C.E. Tilley and J.H. Scoon, Basalts from the northern part of the rift zone of the Mid-Atlantic Ridge, J. Petrol. 5 (1964) 409-434. I.D. Muir, C.E. Tilley and J.H. Scoon, Basalts from the northern part of the Mid-Atlantic Ridge, II. The Atlantis collection from near 30°N, J. Petrol. 7 (1966) 193-201. G.D. Nicholls, A.J. Nalwalk and E.E. Hays, The nature and composition of rock samples dredged from the Mid-Atlantic Ridge between 22ON and 52ON, Mar. Geol. 1 (1964) 333343. M.J. O’Hara, Are ocean floor basalts primary magma?, Nature 220 (1968) 683-686. R.W. Kay, N.J. Hubbard and P.W. Gast, Chemical characteristics and origin of ocean ridge volcanic rocks, J. Geophys. Res. 75 (1970) 1585-1613. C.H. Donaldson and R.W. Brown, Refractory megacrysts and magnesium-rich melt inclusions within spine1 in oceanic tholeiites: indicators of magma mixing and parental magma composition, Earth Planet. Sci. Lett. 37 (1977) 81-89. H. Sato. Nickel content of basaltic magmas: identification of primary magmas and a measure of the degree of olivine fractionation, Lithos IO (1977) 113-120. J.R. Cann, New model for the structure of the ocean crust, Nature 226 (1970) 928-930. J.R. Cann, A model for oceanic crustal structure developed, Geophys. J.R. Astron. Sot. 39 (1974) 169-187. D. Greenbaum, Magmatic processes at ocean ridges and evidence from the Troodos Massif, Cyprus, Nature 238 (1972) 18-21. J. Orcutt, B. Kennett, D. LeRoy and W. Prothero, A low velocity zone underlying a fast spreading rise crest, Nature 256 (1975) 475-476. B.R. Rosendahl, Evolution of oceanic crust, 2. Constraints, implications and inferences, J. Geophys. Res. 81 (1976) 5305-5314. B.R. Rosendahl, R.W. Raitt, L.M. Dorman, L.D. Bibee, D.M. Hussong and G.H. Sutton, Evolution of oceanic crust, I. A physical model of the East Pacific Rise crest derived from seismic refraction data, J. Geophys. Res. 81 (1976) 5294-530s. C.M.R. Fowler, The crustal structure of the Mid-Atlantic Ridge crest at 37”N, Geophys. J.R. Astron. Sot. 47 (1976) 459-49 1. C.M.R. Fowler, The Mid-Atlantic Ridge: structure at 45’N, Geophys. J.R. Astron. Sot. 54 (1978) 167-183. D.P. Blanchard, J.M. Rhodes, M.A. Dungan, K.V. Rodgers, C.H. Donaldson, J.C. Brannon, J.W. Jacobs and E.K. Gibson, The chemistry and petrology of basalts from Leg 37 of the Deep Sea Drilling Project, J. Geophys. Res. 81 (1976) 423 I-4246. M.A. Dungan and J.M. Rhodes, Residual glasses and melt inclusions in basalts from D.S.D.P. Legs 45 and 46: evidence for magma mixing, Contrib. Mineral. Petrol. 67 (1978) 417-431. M.F.J. Flower, W. Ohnmacht, H.-U. Schmincke, I.L. Gibson, P.T. Robinson and R. Parker, Petrology and geochemistry of basalts from Hole 396B, Leg 46, in L. Dmitriev, J.
104
21
22
23
24 25
26
27
28 29
30
31
32 33
34
Heirtzler et al., Initial Reports of the Deep Sea Drilling Project, 46 (U.S. Government Printing Office, Washington, D.C., 1978) 179-213. CR. Allen, The petrology of a portion of the Troodos plutonic complex, Cyprus Ph.D. Thesis, University of Cambridge (1975) unpublished. J.D. Smewing, Mixing characteristics and compositional differences in mantle-derived melts beneath spreading axes: evidence from cyclically layered rocks in the ophiolite of North Oman, J. Geophys. Res. 86 (1981) 2645-2659. M.J. O’Hara, Geochemical evolution during fractional crystallization of a periodically refilled magma chamber, Nature 266 (1977) 503-507. N.H. Sleep, Formation of oceanic crust, some thermal constraints, J. Geophys. Res. 80 (1975) 4037-4042. N.J. Kusznir and M.H.P. Bott, A thermal study of the formation of oceanic crust, Geophys. J.R. Astron. Sot. 47 (1976) 83-95. N.J. Kusznir, Thermal evolution of the oceanic crust; its dependence on spreading rate and effect on crustal structure, Geophys. J.R. Astron. Sot. 61 (1980) 167-181. E.G. Nisbet and C.M.R. Fowler, The Mid-Atlantic Ridge at 37 and 45’N: some geophysical and petrological constraints, Geophys. J.R. Astron. Sot. 54 (1978) 631-660. J.R. Cann, Rb, Sr, Y, Zr and Nb in some ocean floor basaltic rocks, Earth Planet. Sci. Lett. 10 (1970) 7-l 1. D.A. Wood, I.L. Gibson and R.N. Thompson, Elemental mobility during zeohte facies metamorphism of the Tertiary basalts of eastern Iceland, Contrib. Mineral. Petrol. 55 (1976) 241-254. W.J. Ludwig, J.E. Nafe and CL. Drake, Seismic refraction, The Sea, Vol. 4, A.E. Maxwell, ed. (Wiley-Interscience, New York, N.Y., 1970) 53-84. G.L. Maynard, Crustal layer of seismic velocity 6.9-7.6 kilometres per second under the deep oceans, Science 168 (1970) 120-121. E.H. Moores and E.D. Jackson, Ophiolites and oceanic crust, Nature 250 (1974) 136-139. M.H. Salisbury and N.I. Christensen, The seismic velocity structure of a traverse through the Bay of Islands ophiolite complex, Newfoundland, an exposure of oceanic crust and upper mantle, J. Geophys. Res. 83 (1978) 805-817. R.S.J. Sparks, P. Meyer and H. Sigurdsson, Density varia-
35
36 37
38 39
40
41
42
43
44
tion amongst mid-ocean ridge basalts: implications for magma mixing and the scarcity of primitive lavas, Earth Planet. Sci. Lett. 46 (1980) 419-430. H.E. Huppert and R.S.J. Sparks. Restrictions on the composition of mid-ocean ridge basalts: a fluid dynamical investigation. Nature 286 (1980) 46-48. S. Blake, Volcanism and the dynamics of open magma chambers, Nature 289 (1981) 783-785. M.J. G’Hara and R.E. Mathews, Geochemical evolution in an advancing, periodically replenished, periodically tapped, continuously fractionated magma chamber, J. Geol. Sot. London 138 (1981) 237-277. E.G. Nisbet and J.A. Pearce, TiO, and a possible guide to past oceanic spreading rates, Nature 246 (1973) 468-469. J.M. Sinton, D.M. Christie, D.S. Wilson. R. Hey, Petrological consequences of rift propagation on oceanic spreading ridges, Abstr., Am. Geophys. Union Chapman Conf. on The Generation of the Oceanic Lithosphere, Virginia, U.S.A.. April (1981). J.R. Delaney, H.P. Johnson, J.L. Karsten, The Juan de Fuca ridge hot spot-propagating rift system: new tectonic, geochemical and magnetic data, J. Geophys. Res. 86 (198 1) 747-750. R.N. Anderson, D.A. Clague, K.D. Klitgord, M. Marshall and R.K. Nishimori, Magnetic and petrological variations along the Galapagos spreading centre and their relation to the Galapagos melting anomaly, Geol. Sot. Am. Bull. 86 (1975) 683-694. D.P. Mattey and I.D. Muir, Geochemistry and mineralogy of basalts from the Galapagos spreading centre, Deep Sea Drilling Project, Leg 54, in: B.R. Rosendahl, R. Hekinian et al.. Initial Reports of the Deep Sea Drilling Project, 54 (U.S. Government Printing Office, Washington, D.C.. 1980) 755-770. J.H. Natland and W.G. Melson, Composition of basaltic glasses from near the East Pacific Rise and Siqueros fracture zone, near 9”N, in: B.R. Rosendahl, R. Hekinian et al., Initial Reports of the Deep Sea Drilling Project, 54 (U.S. Government Printing Office, Washington, D.C., 1980) 705723. J.R. Cann, Rayleigh fractionation with continuous removal of liquid, Earth Planet. Sci. Lett. 60 (1982) 114-I 16 (this issue).