# Yamauchi & Takei anelastic shear wave velocity-temperature conversion benchmark#

*This section was contributed by Fred Richards.*

This benchmark tests the implementation of the anelastic shear wave
velocity-to-temperature conversion derived by Yamauchi and Takei [2016]
based on forced-oscillation experiments conducted at seismic
frequencies on polycrystalline borneol. This anelasticity parameterization has
been calibrated against a range of observational constraints on upper mantle
temperature, attenuation and viscosity structure, using the surface wave
tomography model of Priestley *et al.* [2012] to constrain
shear wave velocity (\(V_S\)) and the plate model of McKenzie *et al.* [2005]
to estimate lithospheric temperature structure. The resulting
\(V_S\)-to-temperature conversion accurately accounts for the strongly
non-linear temperature dependence of \(V_S\) at near-solidus conditions and is
therefore especially useful for initializing models with accurate temperature
structure in the upper \(\sim\) 400 km of the mantle. This benchmark is
located in the folder benchmarks/yamauchi_takei_2016_anelasticity.

The parameterization of Yamauchi and Takei [2016] defines \(V_S\) as

where \(\rho\) is the density and \(J_1\) and \(J_2\) represent real and imaginary components of the complex compliance, \(J^*\), which is a quantity describing the sinusoidal strain resulting from the application of a unit sinusoidal stress. \(J_1\) represents the strain amplitude in phase with the driving stress, whilst the \(J_2\) component is \(\frac{\pi}{2}\) out of phase, resulting in dissipation. Density is calculated using the expression

where \(\rho_{0} = 3291~\text{ kg m}^{-3}\) and \(\alpha = 3.59 \times 10^{-5}~\text{ K}^{-1}\) are the density and thermal expansivity corresponding to \(T_{0} = 873~\text{ K}\), \(P\) = pressure and \(K = 115.2~\text{ GPa}\) is the bulk modulus. \(J_1\) and \(J_2\) are expressed as

where \(A_B = 0.664\) and \(\alpha_B = 0.38\) represent the amplitude and slope of background stress relaxation and \(J_U\) is the unrelaxed compliance. Parameters \(A_P\) and \(\sigma_P\) represent the amplitude and width of a high frequency relaxation peak superimposed on this background trend such that

and

where \(T^{\prime}\) is the homologous temperature (\(\frac{T}{T_{s}}\)) with \(T\) the temperature and \(T_s\) the solidus temperature, both in Kelvin. \(\phi_m\) is the melt fraction and \(\beta(\phi_m)\) describes the direct poroelastic effect of melt (assumed to be negligible under upper mantle conditions). For this case, \(J_U\) is the inverse of the unrelaxed shear modulus, \(\mu_{U}(P,T)\), such that

where \(\mu_U^0\) is the unrelaxed shear modulus at surface pressure-temperature conditions, the differential terms are assumed to be constant and the pressure, \(P\), in GPa is linearly related to the depth, \(z\), in km by \(\frac{z}{30}\). The normalised shear wave period, \(\tau_S^{\prime}\), in (51) is equal to \(\frac{\tau_S}{2\pi\tau_M}\), where \(\tau_S = 100~\text{ s}\) is the shear wave period and \(\tau_M = \frac{\eta}{\mu_U}\) is the normalised Maxwell relaxation timescale. \(\tau_P^{\prime}\) represents the normalised shear-wave period associated with the center of the high frequency relaxation peak, assumed to be \(6 \times 10^{-5}\). The shear viscosity, \(\eta\), is

where \(d\) is the grain size, \(m\) the grain size exponent (assumed to be 3), \(R\) the gas constant, \(E_a\) the activation energy and \(V_a\) the activation volume. Subscripts \([X]_r\) refer to reference values, assumed to be \(d_r = d = 1~\text{ mm}\), \(P_r = 1.5~\text{ GPa}\) and \(T_r = 1473~\text{ K}\) for the upper mantle. \(A_{\eta}\) represents the extra reduction of viscosity due to an increase in \(E_a\) near the solidus, expressed as

where \(T^{\prime}_\eta\) is the homologous temperature above which activation energy becomes \(E_a + \Delta E_a\), and \(\gamma = 5\) is the factor of additional reduction. \(\lambda\phi\) describes the direct effect of melt on viscosity, assumed to be negligible here. The solidus temperature, \(T_s\), is fixed to a value of 1599 K at 50 km equivalent to a dry peridotite solidus [Hirschmann, 2000] and linearly increases below this depth according to

where \(\frac {\partial T_s}{\partial z}\) is the solidus gradient.

In this benchmark, a 2D input ASCII file containing \(V_S\) specified at 50 and 75 km depth is read in and converted to temperature using an initial temperature model that implements the \(V_{S}(P,T)\) formulation detailed above. The default parameters governing the relationship between \(V_S\) and temperature are set to the values calibrated by Yamauchi and Takei [2016], where \(\mu_U^0 = 72.45~\text{ GPa}\), \(\frac{\partial{\mu_U}}{\partial{T}} = -0.01094~\text{ GPa K}^{-1}\), \(\frac{\partial{\mu_U}}{\partial{P}} = 1.987\), \(\eta_r = 6.22 \times 10^{21}~\text{ Pa s}\), \(E_a = 452.5~\text{ kJ mol}^{-1}\), \(V_a = 7.913 \times 10^{-6}~\text{ m}^{3}~\text{ mol}^{-1}\) and \(\frac{\partial T_s}{\partial z} = 1.018~\text{ K km}^{-1}\). As \(V_S\) is a complex function of temperature, a Brent minimization algorithm is used to find optimal values. Fig. 224 shows that the implementation of this parameterization can accurately recreate the results shown by Yamauchi and Takei [2016] in their Fig. 20.

The parameter file and initial temperature model for this benchmark can be found at benchmarks/yamauchi_takei_2016_anelasticity/yamauchi_takei_2016_anelasticity.prm and benchmarks/yamauchi_takei_2016_anelasticity/anelasticity_temperature.cc. Code to recreate Fig. 224 is provided in benchmarks/yamauchi_takei_2016_anelasticity/plot_output.