УДК 550.37
LARGE-SCALE THREE-DIMENSIONAL INVERSION OF EARTHSCOPE MT DATA USING THE INTEGRAL EQUATION METHOD © 2010 г. Michael S. Zhdanov, Alisa Green, Alexander Gribenko, Martin Cuma
Department of Geology and Geophysics, University of Utah Received January 25, 2010
In this paper we apply 3D inversion to MT data collected in the Northwestern United States as a part of the EarthScope project. By the end of 2009 MT data had been collected from 262 stations located throughout Oregon, Washington, Idaho, and most of Montana and Wyoming. We used data from 139 MT stations in this analysis. We developed fully parallelized rigorous 3D MT inversion software based on the integral equation method with variable background conductivity. We also implemented a receiver footprint approach which considerably reduced the computational resources needed to invert the large volumes of data covering vast areas. The data set used in the inversion was obtained through the Incorporated Research Institutions for Seismology (IRIS). The inversion domain was divided into 2.7 M cells. The inverted electrical conductivity distribution agrees reasonably well with geological features of the region.
INTRODUCTION
EarthScope is the United States Earth Science Program to explore the structure and evolution of the North American continent and understand processes controlling earthquakes and volcanoes. A major part of the EarthScope project is the Transportable Array which will be deployed over the next decade over the entire continental United States. This Transportable Array will provide an unparalleled means to study the geology
ofthe United States. through seismology and other geophysical data. EMScope is the magnetotelluric component of the National Science Foundation's EarthScope USArray program. The MT Transportable Array comprises shorter-period investigations at hundreds of sites in the continental United States. By the end of 2009 MT data had been collected from 262 stations located throughout Oregon, Washington, Idaho, and most of Montana and Wyoming, 139 of which were used in this analysis (Fig. 1).
WAA04 j» f
ЭРI
» WAAÏÛ f щ I DA11
й
200 km
Fig. 1. The map of the Northwestern United States with the locations of the EarthScope MT stations shown by yellow pins.
The unique geological structure of the Northwestern United States is very important both for the study of its geodynamical history and for understanding the physical processes controlling earthquakes and volcanic eruptions [Bishop, 2003]. For such a complex region, definitive structural interpretations based purely on seismological observations may not be sufficient for reliable study ofthe deep earth interior. Conductivity in the subsurface plays a significant role in determining subsurface tectonic activities because ofits sensitivity to temperature, the presence of interstitial fluids, melts, volatiles, and bulk composition. For example, many geological provinces merge in Oregon from the Columbia River basalts in the north to the Basin and Range in the south [Baldwin, 1981]. It is a tectonically active region with the subducting Juan de Fuca plate and volcanically interesting from the effects of the North American Plate traveling over the potential plume currently located beneath Yellowstone.
During recent years, magnetotelluric (MT) methods have experienced a rapid development. Significant improvements have been made both in MT data acquisition systems and in the quality of processing and analyzing this data. Modern MT surveys provide high quality data from an array of densely distributed MT stations, which contain unique information about the geoelectrical structure of geological formations. However, development of a truly three-dimensional (3D) inversion method still represents a very challenging numerical and practical problem. The reasons are twofold. First, 3D forward modeling is a highly complicated and time-consuming mathematical problem itself, especially for large-scale geoelectrical models. Second, the inversion of MT data is an unstable and nonunique problem. One should use regularization methods and physical constraints to obtain a stable and geologically meaningful solution of the inverse problem. There are several algorithms available now for 3D MT inversion. Some of these algorithms are based on rigorous forward modeling [Newman and Alumbaugh, 1997; Sasaki, 2004; Mackie and Watts, 2004; Siripunvaraporn 2005, Gribenko and Zhdanov, 2007; Zhdanov, 2009], while the others use approximate forward modeling operators [Zhdanov and Fang, 1996; Golubev et al., 1999; Zhdanov, 2000].
In this paper we introduce a method of rigorous 3D inversion ofMT data, based on the integral equation (IE) method. We use the re-weighted regularized conjugate gradient method (RRCG) for nonlinear MT inversion [Zhdanov, 2002]. The main distinguishing feature of the RRCG algorithm is application of the special stabilization functionals which allow construction of both smooth images of the underground geoelectrical structures and models with sharp geoelectrical boundaries.
Besides the already complicated problem of nonlinear 3D EM inversion, there are computational difficulties associated with the vast inversion domains of regional MT surveys as well as the large amount of data. The so-called receiver "footprint" approach has been applied in airborne EM data inversion [Reid et al.,
2006; Cox and Zhdanov, 2007] in an effort to reduce the computer resources required for 3D inversion of the real data sets. In this paper, we develop the footprint approach for MT data inversion. In the framework of this technique we use data from a receiver to recover conductivity structures only within a certain horizontal distance from this receiver. However, the entire anomalous domain is included in the computation of the predicted fields in the receivers. A footprint approach allows us to reduce dramatically computer memory requirements without loss of accuracy.
Another well known problem in 3D MT data inversion is the removal of static shift. The effect of static shift is due to the presence of small-scale, near-surface inho-mogeneities. It manifests itself as a vertical shifting of the apparent resistivity curve by a frequency-independent factor, without any corresponding change in the phase curve. The amount of static shift varies from site to site, and differs depending on the source polarization. Thus, the interpretation of static-shifted MT data will obviously lead to erroneous results unless static shifts are correctly taken into account. Smith [1996] proposed static divergence correction for static shift. Another approach is to solve simultaneously for both resistivities and static shift parameters during 3D inversion [Sasaki, 2004]. We apply a static shift correction on the first stage of interpretation and then run 3D inversion of the corrected data.
In this paper we present the results of the 3D inversion of principal MT impedances from 139 MT stations collected at 28 frequencies ranging from 0.00006 to 0.3 Hz. As the forward modeling engine, we use the integral equation (IE) method [Hursan and Zhdanov, 2002]. We use the quasi-Born (QB) approximation for the Frechet derivative calculation. This approximation is a reduction of the more general quasi-analytical approximation with variable background (QAVB) and possesses all of its qualities. QAVB is an effective and accurate technique for EM data inversion [Gribenko and Zhdanov, 2007; Zhdanov, 2009]. As a result of QB application, the IE-based method ofMT inversion requires just one forward modeling in every iteration step, which results in a relatively fast but rigorous inversion method.
The major geoelectrical features recovered by the inversion agree well with the existing global geoelectrical models of the lithosphere. The 3D conductivity structures can be attributed to the specifics of the regional geology. Our inversion results show similar conductivity trends compared to other published works.
3D INVERSION OF THE MT DATA IE method for 3D EM modeling
The integral equation (IE) method represents a powerful tool for electromagnetic (EM) numerical modeling and inversion [Weidelt, 1975; Hohmann, 1975; Wannamaker, 1991; Zhdanov, 2002; 2009]. This method is based on the reduction of the system of
Maxwell's equations to a system of integral equations with respect to the electric field within the inhomoge-neity. The IE equation method can be effectively used for forward modeling the MT data.
Fig. 2 represents a sketch of a geoelectrical model used in the IE method. It is well known that the EM field in this model can be presented as a sum of the background (normal) {EA, HA} and anomalous {Ea, Ha} fields:
E = E" + Ea, H = H + H'
(1)
where the background field is a field generated by the given sources in the model with a background distribution of conductivity CTb, and the anomalous field is produced by the anomalous conductivity distribution Aa{r}, r e Vc R3(see Fig. 2). Then, the electric and magnetic fields can be obtained by the following integral expressions:
E(r') = JJjG^r', r)ActE(r)dv + E(r'), H(r') = Jj"jGff(r', r)ActE(r)dv + Hb(r'),
(2)
(3)
E, H (rj)
ab(z)
a(r)
Fig. 2. Sketch of a geoelectrical model used in the IE method.
Z =
7 7
xx xy
7 7
^yx yy_
(4)
The expressions for the components of the full impedance tensor are [Berdichevsky and Zhdanov, 1984]:
where r e Vc R3, and GEand GH are electric and magnetic Green's tensors.
These expressions become equations when r' e V. A solution of each of these equations exists and is unique, as a solution of Fredholm equations of the second kind.
The process of solving the forward electromagnetic problem according to equations (2) and (3) consists of two parts. First, it is necessary to find the electric and magnetic fields inside the domain V (where Act ^ 0), which requires the solution of an integral equation (domain equation) (2) for r' e V. Second, using the dat
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.