We present a conceptual and numerical approach to model processes in the Earth's interior that involve multiple phases that simultaneously interact thermally, mechanically and chemically. The approach is truly multiphase in the sense that each dynamic phase is explicitly modelled with an individual set of mass, momentum, energy and chemical mass balance equations coupled via interfacial interaction terms. It is also truly multi-component in the sense that the compositions of the system and its constituent thermodynamic phases are expressed by a full set of fundamental chemical components (e.g. SiO$_2$, Al$_2$O$_3$, MgO, etc) rather than proxies. In contrast to previous approaches these chemical components evolve, react with, and partition into, different phases with different physical properties according to an internally-consistent thermodynamic model. This enables a thermodynamically-consistent coupling of the governing set of balance equations. Interfacial processes such as surface tensions and/or surface energy contributions to the dynamics and energetics of the system are also taken into account. The model presented here describes the evolution of systems governed by Multi-Phase Multi-Component Reactive Transport (MPMCRT) based on Ensemble Averaging and Classical Irreversible Thermodynamics principles. This novel approach provides a flexible platform to study the dynamics and non-linear feedbacks occurring within various natural systems at different scales. This notably includes major-and trace-element transport, diffusion-controlled trace-element re-equilibration or rheological changes associated with melt generation and migration in the Earth's mantle.
This is a pre-copyedited, author-produced PDF of an article accepted for publication in
Geophysical journal international following peer review. The version of record Oliveira, B.,
Afonso, J., Zlotnik, S., Diez, P. Numerical modelling of multi-phase multi-component reactive transport in the Earth's interior. "Geophysical journal international", 1 Gener 2018, vol. 212, núm. 1, p. 345-388 is available online at: https://doi.org/10.1093/gji/ggx399.
Practical applications of controlled-source electromagnetic (EM) modelling require solutions for multiple sources at several frequencies, thus leading to a dramatic increase of the computational cost. In this paper, we present an approach using block Krylov subspace solvers that are iterative methods especially designed for problems with multiple right-hand sides (RHS). Their main advantage is the shared subspace for approximate solutions, hence, these methods are expected to converge in less iterations than the corresponding standard solver applied to each linear system. Block solvers also share the same preconditioner, which is constructed only once. Simultaneously computed block operations have better utilization of cache due to the less frequent access to the system matrix. In this paper, we implement two different block solvers for sparse matrices resulting from the finite-difference and the finite-element discretizations, discuss the computational cost of the algorithms and study their dependence on the number of RHS given at once. The effectiveness of the proposed methods is demonstrated on two EM survey scenarios, including a large marine model. As the results of the simulations show, when a powerful preconditioning is employed, block methods are faster than standard iterative techniques in terms of both iterations and time.
One thousand seven hundred and eighty-six vertical-component, short-period observations of microearthquake codas from regional earthquakes recorded by 17 stations belonging to the National Seismological Network of Colombia were used to estimate seismic wave attenuation in Colombia. Local magnitudes range from 2.9 to 6.0 and only events occurring at hypocentral distances up to 255 km were considered for the analysis. The frequencies of interest lay between 1 and 19 Hz and the analysis was performed for each seismic station separately. Coda-wave attenuation (Q-1c) was estimated by means of a single-scattering method whereas the separation of intrinsic absorption (Q-1i) and scattering attenuation (Q-1s) from total attenuation (Q-1t) was performed using a multiple lapse time-window analysis based on the hypothesis of multiple isotropic scattering and uniform distribution of scatterers. A regionalization of the estimated Q0 (Qc at 1 Hz) values was performed and a contour map of seismic coda attenuation in Colombia is presented, where four zones with significant variations of attenuation related to different geological and tectonic characteristics can be observed. The highest attenuation is linked to the central and western regions (Q0 around 50 and 56) whereas a lower attenuation (Q0 around 69 and 67) is assigned to the northern and eastern regions. Results show that the Q-1 values are frequency dependent in the considered frequency range, and are approximated by a least-square fit to the power law Q-1(f) =Q-10(f/f0)-?. The exponents of the frequency dependence law ranged from ?= 0.65 to 1.01 for Q-1c, ?= 0.62 to 1.78 for Q-1i, ¿= 0.28 to 1.49 for Q-1s, and ¿= 0.53 to 1.67 for Q-1t. On the other hand, intrinsic absorption is found to dominate over scattering in the attenuation process for most of the stations and frequency bands analysed. Some discrepancies have been observed between the theoretical model and the observations for some frequency bands which indicate that it would be necessary to consider models for depth-dependent velocity structure and/or non-isotropic scattering patterns.
Caselles, J.; Canas, J.; Pujades, L.G.; Lana, F.J.; Badal, J.; Corchete, V.; Navarro Bernal, Manuel Geophysical journal international Vol. 126, num. 3, p. 269-281 DOI: 10.1111/j.1365-246X.1997.tb04505.x Data de publicació: 1997-06 Article en revista
Attenuation coefficients of Rayleigh waves propagating across the Iberian Peninsula, and quality factor models were obtained for different two-station great-circle paths. The data used were long-period wave trains contained in the seismograms provided by the ILIHA array installed in the peninsula in the 1980s. More than 300 seismograms were analysed to obtain 14 sets of similar paths, consisting of either paths totally included in one of the four pure tectonic regions or mixed paths crossing several regions. We used the two-station method to determine attenuation parameters (¿R) and then an inversion method to derive Qß-1 models. The attenuation coefficients of seismic energy in the 10-120 s range are in the range 2.8 × 10-3-0.1 × 10-3 km-1, and correlate well with both previous results obtained from analogue WWSSN and most seismogeotectonic characteristics of Iberia. The Qß models obtained also agree with a few models for Iberia determined earlier, and they show significant variations from one region to another. As a general feature, the ¿-values are higher in the active seismotectonic regions (typical maximum of 2.6 × 10-3 km-1) than in the more stable ones (typical maximum of 1.4 × 10-3). Moreover, the maximum of Qß-1 in the active regions (145 × 10-3 -120 × 10-3) is higher than in the stable ones (40 × 10-3 -100 × 10-3). In general, both the regionalized attenuation coefficients of Rayleigh waves and the quality factor models exhibit a good agreement with previously published Lg wave attenuation and coda Q distribution patterns in Iberia.
Pujades, L.G.; Ugalde, A.; Canas, J.; Navarro Bernal, Manuel; Badal, J.; Corchete, V. Geophysical journal international Vol. 129, num. 1, p. 281-291 DOI: 10.1111/j.1365-246X.1997.tb01581.x Data de publicació: 1997-01 Article en revista
The anelastic attenuation in the Almeria Basin (southeastern Iberian Peninsula) is investigated by using seismic data collected during the summer of 1991. A multiple-lapse time-window analysis is applied to high-frequency seismograms corresponding to 20 shallow seismic events with low magnitudes (m= 2.5) and distances less than 71 km, recorded at six short-period seismographic stations. We have constructed corrected geometrical spreading and normalized energy-distance curves for the region over the frequency bands 1-2, 2-4, 4-8, 8-14 and 14-20 Hz. A theoretical model for body-wave energy propagation in a randomly heterogeneous medium has been employed to interpret the observations. Two parameters describe the medium in this model: the scattering attenuation coefficient ¿S= kQ-1S and the intrinsic attenuation coefficient ¿1= kQ-11, where k is the wavenumber and Q-11 and Q–1S are the intrinsic and scattering attenuation respectively. This model assumes that scattering is isotropic, including all orders of multiple scattering, and predicting the spatial and temporal energy distribution of seismic energy. A least-squares fitting procedure has been used to find the best estimates of the model parameters. The analysis of the spectral amplitude decay of coda waves has provided coda Q–1C values at the same frequency bands. The results obtained show that Q–11, Q–1S and Q–1C decrease with increasing frequency; for frequencies lower than 3 Hz scattering attenuation is stronger than intrinsic absorption and coda Q–1C takes values between intrinsic and total attenuation, being very close to Q–1T. Q–1S is more frequency-dependent than Q–1I; for frequencies greater than 3 Hz intrinsic absorption is the dominant attenuation effect and Q–11 and Q–1S have significant frequency dependence. In order to correlate the results obtained with the major geological and tectonic features of the region, a geotectonic framework for the area is provided and the predominant frequency decay in coda waves is analysed in order to obtain the coda Q frequency dependence following a power law QC= Q0(f/f0)¿, where f0 is a reference frequency. In this way we have obtained regionalized values of coda Q at 1 Hz (Q0). Finally, a first-order approach has allowed us to obtain intrinsic and scattering quality factors from the obtained Q0 and ¿ values, leading us to obtain tentative distributions of QI. QS and Q0 at 1 Hz for the area. The derived intrinsic and scattering quality-factor distributions are in good agreement with the tectonic history and the main geological features of the region. Large scattering and intrinsic attenuation (QS~ 80, QI~ 100) are found in the sedimentary Neogene and Quaternary basin, while scattering is the dominant effect in the old Palaeozoic rocks of the mountains (QS~ 200, QI~ 1000). Intrinsic Q shows a higher sensitivity to the geological characteristics than scattering Q.
Badal, J.; Corchete, V.; Payo, G.; Pujades, L.G.; Canas, J. Geophysical journal international Vol. 124, num. 1, p. 591-611 DOI: 10.1111/j.1365-246X.1996.tb07039.x Data de publicació: 1996-02 Article en revista
A detailed dispersion analysis of Rayleigh waves propagating across the Iberian Peninsula is carried out. The starting data are high-quality long-period data recorded at the broad-band NARS stations installed in (he Iberian Peninsula during the ILIHA project. We apply methods to obtain a correct selection of data and subsequent two-station surface-wave velocity measurements. A total of 64 teleseismic events recorded by the NARS array and 143 seismic paths have been studied. Several techniques which provide a significant improvement in the signal-to-noise ratio are employed to remove higher-mode interference efficiently and improve the isolation of the fundamental-mode Rayleigh wave from the seismograms. Thereafter, the interstation Rayleigh wave phase and group velocities are determined. We perform simultaneous inversion of phase-and group-velocity dispersion data by means of the stochastic inverse operator, and lest the reliability of the results by computing resolving kernels and also by forward modelling. A regionalization procedure based on the Backus-Gilbert approach for linear inversion of traveltime data is applied.
Both the inversion results and the contoured shear-wave velocity panoramas display the main features of the deep structure of Iberia. We find a subcrustal low-velocity channel which extends over practically the whole peninsular area and spans a depth interval of approximately 40–50 km; it exhibits velocities of between 4.30 and 4.50 km s-1. At depths of 66–81 km, we find the highest velocities in the lithosphere, which reach values of 4.85 km s-1 in many cases. The low-velocity channel of the asthenosphere spans a large depth interval of approximately 80-180 km; it shows the lowest velocity values computed by us. We find velocities decreasing with depth, which are of the order of 4.25–4.36 km s-1 for the first 40 km and of the order of 4.00–4.25 km s-1 for the rest. The upper mantle under the asthenosphere exhibits high velocities, which range between 4.62 and 4.82 km s-1 in most cases.
The shear-wave velocity structure of the Iberian subcrustal lithosphere and asthenosphere is mapped at 11 depth intervals from 24 to 201 km. At the top of the mantle, relatively low velocities span the Ebro Valley and also the southern third of the peninsula. Low velocities appear in the south-southwest quadrant, and high velocities occur over the Hercynian basement. At greater lithospheric depths, very low velocities extending over the whole peninsula suggest a low-velocity channel of non-uniform lateral structure, where a reduced zone to the west of the Iberian plateau shows relatively high velocities. At the greatest lithospheric depths, the whole Iberian block is fairly homogeneous laterally. The asthenosphere shows a notable lateral heterogeneity as well. We distinguish two parts: the upper asthenosphere, a 40 km thick layer with predominant velocity values of 4.25 km s or more; and the lower asthenosphere, a 60 km thick layer with velocity values generally below 4.25 km s-1. The upper asthenosphere seems to be less laterally heterogeneous than the lower asthenosphere. The lower asthenosphere exhibits a more pronounced negative velocity gradient than the upper asthenosphere.
Badal, J.; Corchete, V.; Payo, G.; Serón, F.; Canas, J.; Pujades, L.G. Geophysical journal international Vol. 108, num. 1, p. 71-88 DOI: 10.1111/j.1365-246X.1992.tb00839.x Data de publicació: 1992-01 Article en revista
A rigorous study of velocity dispersion of surface waves generated by teleseismic events propagating across the Iberian Peninsula and traversing main geological units, has been carried out from a set of selected analogue data, as digital records have only become available recently. Dispersed seismic signals have been obtained over a period of 16 years, between 1967 and 1982, at the five Iberian stations having long-period instruments. In our study, we have considered many earthquakes thus obtaining a fairly good path coverage of most of the peninsula for two-station Rayleigh wave velocity measurements. In all cases, the approach azimuths of the wavefronts were carefully checked. Several digital filtering techniques have been employed to remove the effects of multipathing and modal contamination, and to isolate the fundamental mode from Rayleigh wavetrains. Thus, we have obtained good estimates for both phase and group velocities. A time-variable filter has reduced the influence of noise and removed higher mode interference. Multiple filtering is then used to compute group velocity. Frequency-domain Wiener deconvolution is used to compute the interstation phase velocity. The determined average Rayleigh wave velocities reveal differences in the propagation conditions of the seismic energy across the peninsula. A mapping of velocities for various periods of reference, together with a mapping of errors in velocity, are the basis for obtaining the Rayleigh wave velocity distribution in the peninsula. Theoretical 2-D layered earth models are obtained by joint inversion of phase and group velocity dispersion curves using the stochastic inverse operator. In our inversion scheme, we use velocities corrected for anelastic effects. Finally, a 3-D mapping of S velocity is performed. This study shows important regional features of the deep structure of Iberia; we see small lateral inhomogeneities and also two low-velocity layers: one with shear velocities usually ranging from 4.23 to 4.31 km s-1 directly under the Moho, and another, the asthenosphere, with a negative velocity gradient for depths between 81 and 181 km, terminated at the bottom by a sharp discontinuity.
A method based on the coda attenuation law: Q=Q0(f/f0)v leads to the determination of the lateral variation of coda-Q in the southern part of the Iberian Peninsula using seismograms belonging to the seismological network of the Cartuja Observatory, located in Granada. The lateral variation of Q0 (Q value corresponding to a reference frequency f0 of 1 Hz) and its frequency dependence for the 1 to 5 Hz frequency range are, in general, in agreement with coda-Q values for frequencies less than about 1 Hz, previously determined in the region under study.
To determine the coda-Q values analytical functions have been used to fit the magnification curves of the vertical component short-period seismographs belonging to the Cartuja network. The problem is solved by using least-squares techniques and non-linear inversion. The determined coda-Q0 values and its frequency dependence correlate well with several known geophysical parameters in the southern part of the Iberian Peninsula.
Badal, J.; Corchete, V.; Pujades, L.G.; Payo, G.; Canas, J.; Seron, F. Geophysical journal international Vol. 100, num. 2, p. 193-202 DOI: 10.1111/j.1365-246X.1990.tb02479.x Data de publicació: 1990-01 Article en revista
Several filtering techniques have been used to remove the effects of multipathing and modal contamination, and to isolate the fundamental mode from Rayleigh wavetrains. Group velocity data are obtained by means of the multiple-filter technique. A time-variable filter has allowed the influence of noise as well as the interference produced by higher modes to be removed. Multiple filtering is then used again to compute group velocities at each station. The interstation group velocity for the fundamental mode Rayleigh wave is estimated according to the velocities at two stations. Frequency-domain Wiener deconvolution is used to compute the phase velocity between two stations. The well-known three-station method is applied to correct the distances travelled by the waves across the array and therefore to determine interstation phase and group velocities in a more accurate manner. On the other hand, lateral refraction at the Atlantic continental edge of the Peninsula is also studied. Phase velocities are corrected for the anelastic effect. Inversion of the interstation Rayleigh wave phase velocities is then made in accordance with generalized inversion theory to obtain theoretical 2-D layered earth models. In this paper, these methods are applied to Rayleigh waves generated by teleseismic events propagating across the Iberian Peninsula and recorded at WWSSN stations. As a consequence, new and principal features for the Iberian lithosphere-asthenosphere system are obtained. A very interesting feature of the the Iberian lithosphere was found—a low-velocity layer directly under the Moho, between 39 and 64 km depth, with shear velocities ranging from 4.12 to 4.37 km s-1. The Iberian asthenosphere, which lies between 100 and 180km depth, is not an homogeneous layer and shows a negative velocity gradient from top to bottom together with a sudden increase in shear velocity beneath the low-velocity zone.
Pujades, L.G.; Canas, J.; Egozcue, J. J.; Puigvi, M.; Gallart, J.; Lana, X.; Pous, J.; Casas, A. Geophysical journal international Vol. 100, num. 2, p. 285-302 DOI: 10.1111/j.1365-246X.1990.tb02486.x Data de publicació: 1990-01 Article en revista
Several attenuation studies have established a frequency dependence law of the anelastic attenuation factor Q in the form Q = Q0(f/f0)v for the approximate 1–10 Hz frequency range. We propose a method that leads to the determination of Q0, which is a function of the reference frequency f0, and the real exponent v with a single station. To carry out the problem we determine a set of master curves as a function of v. We discuss the method, and the different features of the master curves, when it is applied to the complicated regions of the Iberian Peninsula and to several instruments with different responses. Using this new method and the seismographic stations available in the Iberian Peninsula we have mapped iso-Q0 lines, at a reference frequency of 1 Hz, applying inversion methods. The Q0 values determined for Iberia vary between about 100 and about 600. Values close to 100 correspond to the southern part of Iberia. In general, Q0 values increase from south to north with values about 600 near the NW part of Iberia. The Pyrenees Mountains and adjacent areas present Q0 values between about 200 and about 350. These results suggest a strong Q0 lateral variation in Iberia. A considerable frequency dependence of coda-Q has also been determined. The v values vary between 0.3 and 0.8. The Q0 values obtained in the Iberian Peninsula show very good agreement with several Q0 values obtained in other regions of the world. Comparison between the iso-Q0 lines and other geophysical parameters, like regional variations of Pn velocities, heat flow, isoseismal intensity distribution and crustal thickness, indicates that lower Q0 values are associated with higher isoseismal intensity attenuation, higher heat flow, lower Pn velocities and thinner crust.
Badal, J.; Canas, J.; Corchete, V.; Pujades, L.G.; Serón, F. Geophysical journal international Vol. 103, num. 1, p. 135-146 DOI: 10.1111/j.1365-246X.1990.tb01758.x Data de publicació: 1990-01 Article en revista
Coda wave analysis is used to obtain frequency-dependent coda-Q values for different seismic zones of the Iberian area. Seventeen source regions around the Geophysical Observatory of Toledo and some four seismic events per region have been considered in this study. We have used an iterative Fourier analysis technique to see the variation of the frequency along the coda, also taking into account the instrument response. We have applied a suitable criterion to select the predominant frequency every 5 s along the coda. The variation of the frequency with time for each region is averaged with a second-degree polynomial, which is compared to master curves obtained directly from the response of the seismograph system, in order to determine the elastic quality factor Q. It has been observed that the frequency-time curves thus obtained are better explained if Q is considered as an exponential function of the peak frequency. The main result of this work is a set of 1 Hz Q values with a clear indication that frequency dependence of Q exists, although the bandwidth from which our conclusions are reached is only 0.5–1 Hz. The coda-Q values obtained for the tectonically most stable areas (north Spain) appear somewhat higher, Q > 300, than those corresponding to the seismic active zones (south Spain), Q < 250. Thus, a clear relationship is established between Q values and the two major tectonic provinces in Iberia. These results may be helpful for seismic risk and earthquake engineering purposes.