16% in the P –T region between the experimentally determined garnet-in
and plagioclase-out phase transitions in basalt^56. The physical param-
eters for the numerical experiments are presented in Extended Data
Table 2^45 ,^46 ,^48 ,^57 –^59.
Hydration model. During the evolution of the model, water is liber-
ated from the subducted oceanic crust as a consequence of dehydration
reactions and compaction (see Sizova et al.^23 ). We assume incomplete
hydration of the mantle wedge as a consequence of the channelization
of slab-derived fluids. To account for this behaviour, we arbitrarily
assign 2 wt.% H 2 O as an upper limit for mantle wedge hydration. The
hydrated mantle is subdivided into two parts: (i) an upper, colder ser-
pentinized zone, and (ii) a lower, warmer, hydrated, but not serpenti-
nized zone. The stable mineralogical composition and water content
were obtained by free energy minimization^60 as a function of pressure
and temperature from thermodynamic data^61. We would like to note
that owing to existing uncertainties of the high-pressure mineralogy
of subducted crust and hydrated mantle in the mantle transition zone
and of the amount of water in the nominally anhydrous minerals, we
cannot very accurately predict the amount of water that is released
from subducting slabs in the mantle transition zone. This amount is
however high enough to trigger the hydrous diapirism processes in
the mantle transition zone^26 ,^27 that are documented in our models.
Melting model. Melting of the mantle and the crust, as well as melt
extraction and percolation across the crust–mantle interface and to
the surface is implemented in a simplified manner^62. According to our
model, mafic magma added to the crust is balanced by melt production
and extraction in the source mantle/crustal region. Melt percolation
is not modelled directly and is considered to be nearly instantane-
ous^63 ,^64 compared to the long computational time step. Lagrangian
markers track the amount of melt extracted during the evolution of
each experiment. The total amount of melt, M, for every marker takes
into account the amount of previously extracted melt and is calculated
as in Vogt et al.^62
MM=−0e∑mMxt
where ∑mMext is the total melt fraction extracted during the previous
m extraction episodes. At the beginning all mantle rocks are assumed
to be melt-depleted (∑=mMMext0) in accordance with the assumed
adiabatic mantle potential temperature profile under a mid-ocean
ridge. The rock is considered to be non-molten (refractory) when the
extracted melt fraction is larger than the standard one (that is, when
∑>mMMext0). If M > 0 for a given marker, the melt fraction Mext = M is
extracted and ∑mMext is updated. In our experiments, all melt is
extracted from a partially molten rock marker after the volume of
melt reaches a threshold of 2 vol.%. As a result, the volume of melt in
the partially molten mantle is always <2 vol.%. Melt extracted from
the mantle rises and adds to the crust as hot melt intrusions (plutons,
80% of all melts) and to the surface as volcanic rocks (20% of all
melts^62 ,^65 ).
One key component of our numerical model is that it takes into
account decreases in the mantle density related to melt extraction.
This density reduction is a well known phenomenon that has been
studied both on the basis of natural data^6 ,^66 and experimentally (for
example, see Schutt and Lesher^67 ). In our models, standard density
of the melt-depleted mantle is therefore corrected for the degree of
depletion as Gerya et al.^38 :
ρρ0(depl)0=(1−0. (^04) ∑mMext)
where ρ0(depl) is the standard density of depleted solid mantle and
∑mMext is the degree of melt extraction that changes with time.
The volumetric degree of melting M 0 in partially molten rocks is com-
puted as in Vogt et al.^62. For the mantle, we use the P–T–H 2 O-dependent
melting model for peridotite as in Katz et al.^47. For crustal rocks,
we assume that the degree of both hydrous and dry melting is a
piecewise-linear function of P–T (for example, ref.^40 ):
M
TT
TT
TTTTT
TT
0<
−
− <<
1>
0 ,
solidus
solidus
liquidus solidus solidusliquidus
liquidus
where Tsolidus and Tliquidus are, respectively, the solidus temperature (wet
and dry solidi are used for the hydrated and dry mantle, respectively)
and dry liquidus temperature at a given pressure and rock composition
(Extended Data Table 2).
The effect of latent heat due to equilibrium melting/crystallization
is included implicitly by increasing the effective heat capacity CP,eff and
thermal expansion coefficient αeff of the partially molten/solidified
material (see ref.^68 ):
CCH
M
T
=+
∂
P,effPL∂ Pconstant
and
ααρ
H
T
M
P
=+
∂
eff ∂ T
L
constant
where CP is the heat capacity and HL is the latent heat of melting.
Visco-plastic rheological model
We employ a visco-plastic rheology with experimentally calibrated
flow laws. The effective creep viscosities of rocks are represented as
a function of temperature and stress by experimentally determined
flow laws (Extended Data Table 2). The viscosity for dislocation creep
depends on strain rate, pressure and temperature and is defined in
terms of deformation invariants^45 as follows:
ηε̇ A EVP
nRT
creep=(II)e(1−)nn/ D−1/n xp + ,
where εIİ̇=1/2εεijij̇ is the square root of second invariant of the strain
rate tensor ε̇ij and AD, E, V and n are the experimentally determined flow
law parameters, the material constant, the activation energy, the acti-
vation volume and the stress exponent, respectively (Extended Data
Table 2). We use the same ductile rheology for both solid and partially
molten rocks, since the maximum melt fraction in our experiments is
limited to 2 vol.% by melt extraction. Melt extraction and related water
removal from the mantle results in the increase of the mantle viscos-
ity^69. Therefore, in our models, we use strong dry olivine rheology for
both lithospheric and asthenospheric mantle (except for mantle rocks
hydrated by subduction, Extended Data Table 2).
The ductile rheology of rocks is combined with a brittle/plastic rhe-
ology to yield an effective visco-plastic rheology. For this purpose the
Drucker–Prager yield criterion^45 is implemented by limiting creep
viscosity, ηcreep, as follows:
η
cPγ
ε
≤
- 2
creep ,
Iİ
where c is the compressive strength at P = 0, γ is the effective internal
friction coefficient, which includes effects of fluid and melt weakening,
γ = γdryλfluid and γ = γdryλmelt, and γdry is the internal friction coefficient
of dry rocks:
λ
P
fluid=1−P
fluid
solid