CN111856596B - Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method - Google Patents

Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method Download PDF

Info

Publication number
CN111856596B
CN111856596B CN202010774885.6A CN202010774885A CN111856596B CN 111856596 B CN111856596 B CN 111856596B CN 202010774885 A CN202010774885 A CN 202010774885A CN 111856596 B CN111856596 B CN 111856596B
Authority
CN
China
Prior art keywords
inversion
model
resistivity
parameter
matrix
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010774885.6A
Other languages
Chinese (zh)
Other versions
CN111856596A (en
Inventor
罗鸣
李予国
裴建新
卢杰
李卓轩
乔荷
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Ocean University of China
Original Assignee
Ocean University of China
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Ocean University of China filed Critical Ocean University of China
Priority to CN202010774885.6A priority Critical patent/CN111856596B/en
Publication of CN111856596A publication Critical patent/CN111856596A/en
Application granted granted Critical
Publication of CN111856596B publication Critical patent/CN111856596B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/38Processing data, e.g. for analysis, for interpretation, for correction
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30Assessment of water resources

Landscapes

  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

The invention provides a layered medium resistivity anisotropy ocean controllable source electromagnetic fast inversion method, which comprises the following steps: reading and converting inversion data; setting inversion execution parameters; setting parameters of an inversion initial model; constructing an anisotropic resistivity inversion target function; obtaining a Jacobian matrix and a Hessian matrix of the electromagnetic field relative to the anisotropic resistivity; adaptively calculating a regularization factor based on the inversion parameter characteristics; calculating the updating amount of the model; calculating a target function fitting difference of the inversion iteration model; judging whether inversion requirements are met; and outputting the final inversion model. The method can simultaneously invert the anisotropic resistivity of the seabed medium and the thickness of the seabed stratum, so that relatively complex seabed stratum distribution can be obtained by using an inversion model with a small number of layers, and the method has the characteristics of high inversion speed and high efficiency, and is suitable for inversion of complex geoelectric models. The inversion method provides a feasible technical means for rapidly detecting data quality and acquiring inversion interpretation results in the field.

Description

Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method
Technical Field
The invention relates to the technical field of marine geophysical inversion, in particular to a rapid inversion method of marine controllable source electromagnetism with anisotropic laminar medium resistivity.
Background
Marine Controlled Source Electromagnetic (CSEM) is a marine geophysical prospecting method for the exploration of subsea hydrocarbon and mineral resources. Research shows that about 30% of oil and gas resources in the world exist in the electrically anisotropic stratum, the electrical anisotropy of the seabed medium is an important influence factor for obtaining the correct electrical distribution of the seabed medium, and if the resistivity anisotropy is neglected in explaining marine electromagnetic data, a reasonable seabed geoelectric model cannot be obtained. However, the inversion interpretation methods widely applied at present are still traditional geophysical inversion methods, which often assume that the seabed medium is electrically isotropic when inverting marine electromagnetic data acquired by marine resource exploration, and such inversion methods may provide inversion results with huge errors for later data interpretation, which may also result in obtaining incorrect interpretation results. For this reason, the properties of resistivity anisotropy of the seafloor medium must be considered when interpreting marine electromagnetic data.
At present, the geophysical inversion interpretation method still takes isotropic inversion as a main method, and the resistivity anisotropic inversion algorithm is not reported much. Analyzing sensitivity characteristics of resistivity vertical anisotropy by Ramanan jaona (2011), and explaining laminar resistivity vertical anisotropy marine CSEM data by using an Occam inversion method; rogover (2016) utilizes a Gauss-Newton method to achieve a one-dimensional vertical anisotropy marine CSEM inversion algorithm, and can better reconstruct the submarine electrical structure of multi-frequency data analyzed by multiple emission sources. However, the above methods only invert the resistivity information of the seabed medium, but do not invert the thickness of the seabed stratum, and the seabed stratum is set to be divided into a plurality of thin layers for inversion when the inversion model is set, so that the inversion speed is low, and the efficiency is low.
Disclosure of Invention
The invention aims to provide a method for quickly inverting marine controllable source electromagnetism with anisotropic resistivity of a layered medium, which can improve the accuracy and efficiency of marine electromagnetic inversion.
In order to achieve the above purpose, the invention provides the following technical scheme:
1. a method for quickly inverting marine controllable source electromagnetism with anisotropic resistivity of a layered medium is characterized by mainly comprising the following steps:
s1, reading and converting ocean controllable source electromagnetic field data participating in inversion, wherein the data comprises real part and imaginary part data of an electromagnetic field, amplitude and phase data, and polarization ellipse long axis and short axis parameters;
s2, inversion execution parameters are set, wherein the parameters comprise inversion maximum iteration times, target fitting difference, maximum iteration step length, step length proportionality coefficient, penalty function type and regularization attenuation coefficient;
s3, setting inversion initial model parameters, wherein the parameters comprise a background layer resistivity parameter, a thickness parameter and an observation system parameter;
s4, constructing an anisotropic resistivity inversion target function;
s5, solving a Jacobian matrix and a Hessian matrix of the electromagnetic field relative to the anisotropic resistivity;
s6, adaptively calculating a regularization factor based on the inversion parameter characteristics;
s7, calculating the model updating amount;
s8, calculating a target function fitting difference of the inversion iterative model;
s9, judging whether inversion requirements are met, if so, turning to S10, and if not, turning to S5;
and S10, outputting the final inversion model.
2. The method for rapidly inverting the marine controllable source electromagnetism of the layered dielectric resistivity as claimed in claim 1, wherein the objective function of constructing the anisotropic resistivity inversion is:
Figure BDA0002618041550000031
wherein φ is the objective function of the inversion algorithm; m is a model inversion parameter vector including a seafloor anisotropic resistivity parameter m ρ And the thickness m of the formation h I.e. m = m ρ +m h
Figure BDA0002618041550000033
Is the gradient of the model parameter vector; i | · | | is a standard deviation operator; d is an observation data vector used for inversion; w d Weighting the matrix for the data; w m Weighting the model with a matrix; f (m) represents a forward response operator of the model m; mu.s ρ And mu h Respectively a seabed anisotropic resistivity parameter m in an inversion model ρ And the thickness m of the formation h The regularization factor of (1).
3. The method for the marine controllable source electromagnetic fast inversion of the laminar medium resistivity anisotropy according to claim 1, wherein the adaptive regularization factor is determined according to the following formula:
Figure BDA0002618041550000032
wherein i represents the ith iteration of inversion; mu.s i Is a regularization factor; max | is the element with the maximum matrix absolute value; a is mj Is a matrix product [ (W) d J) T (W d J)]The element (b); m is the matrix product [ (W) d J) T (W d J)]Dimension of (d); χ is the attenuation coefficient; lambda is the transverse resistivity rho of the submarine stratum of the inversion model h Resistivity p in the vertical direction v And a weighting factor for the formation thickness h, determined using the following equation
Figure BDA0002618041550000041
Wherein, gamma is a proportionality coefficient, and m is an inversion model parameter.
4. The method for rapidly inverting the marine controllable source electromagnetism of laminar medium resistivity anisotropy according to claim 1, wherein the method for calculating the model parameter update amount based on the inversion parameters comprises the following steps:
Figure BDA0002618041550000042
wherein, Δ m is the model parameter update amount of the next iteration, i is the ith inversion iteration, H i Is Hessian matrix, g i Is the gradient of the objective function.
5. The inversion method of sea bottom anisotropic resistivity and emitter position according to claim 4, wherein the method of obtaining the minimum value of the objective function comprises:
Figure BDA0002618041550000043
where i is the ith iteration of inversion, J i Is a jacobian matrix.
Compared with the prior art, the method provided by some embodiments of the invention has the beneficial effects that:
the invention mainly aims at the problems of seabed medium anisotropy widely existing in the detection of a marine controllable source and the problems of low inversion speed and low efficiency of the anisotropy of the electromagnetic resistivity of the marine controllable source, and provides a layered medium resistivity anisotropic marine controllable source electromagnetic fast inversion method. The efficiency of the inversion algorithm has a direct relation with the complexity of the inversion model, in order to obtain a more accurate inversion result, the traditional inversion algorithm needs to set an inversion model consisting of dozens of or even hundreds of thin layers for inversion, and the inversion efficiency is seriously influenced due to the large size of the inversion model. However, the invention only needs to set the seabed stratum into a simple model with several layers, continuously modifies the thickness and anisotropic resistivity of the stratum in the inversion iteration process, fits the stratum resistivity distribution and layer thickness of the real model, and finally approaches the real model continuously. Compared with the traditional inversion algorithm, the algorithm provided by the invention has the characteristics of high inversion speed and high efficiency, and is suitable for the inversion of a complex earth-electric model. Meanwhile, the inversion method provides a feasible technical means for rapidly detecting data quality and acquiring inversion interpretation results in the field.
Drawings
In order to more clearly illustrate the technical solutions in the embodiments of the present invention, the drawings needed to be used in the embodiments or the prior art descriptions will be briefly described below, and it is obvious that the drawings in the following description are only some embodiments of the present invention, and it is obvious for those skilled in the art to obtain other drawings based on these drawings without inventive exercise.
FIG. 1 is a block flow diagram of the process of the present invention;
FIG. 2 is a structural diagram of a one-dimensional resistivity isotropic model;
FIG. 3a is a diagram illustrating sensitivity of a horizontal electric field Ey of a one-dimensional resistivity isotropic model with respect to lateral conductivity of a seabed medium;
FIG. 3b is a diagram showing sensitivity of a horizontal electric field Ey of a one-dimensional resistivity isotropic model with respect to vertical conductivity of a seabed medium;
FIG. 3c is a graph showing sensitivity of a one-dimensional resistivity isotropic model horizontal electric field Ey to isotropic conductivity of the seafloor medium;
FIG. 3d is a diagram showing the relative error of the horizontal electric field Ey of the one-dimensional resistivity isotropic model with respect to the sensitivity of the isotropic conductivity of the seabed medium;
FIG. 4 is a structural diagram of a one-dimensional resistivity anisotropy model;
FIG. 5a is a variation curve of resistivity anisotropy inversion data fitting difference and inversion iteration number;
FIG. 5b is a schematic diagram of the resistivity anisotropy inversion results;
FIG. 6a is a schematic diagram of a data fit between Ey component observation data participating in resistivity anisotropy inversion and inversion model forward response;
FIG. 6b is a schematic diagram of a data fit between Bx component observation data participating in resistivity anisotropy inversion and inversion model forward response;
FIG. 6c is a schematic diagram of data fitting of Ez component observed data participating in resistivity anisotropy inversion and inversion model forward response;
Detailed Description
In order to make the technical problems, technical solutions and advantageous effects to be solved by the present invention more clearly apparent, the present invention is further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the invention and are not intended to limit the invention.
Example 1
The invention provides a method for quickly inverting the anisotropic ocean controllable source electromagnetic wave of the layered medium resistivity, which can simultaneously invert the anisotropic resistivity of the ocean bottom medium and the thickness of the ocean bottom stratum, so that relatively complex ocean bottom stratum distribution can be obtained by using an inversion model with fewer layers. Meanwhile, the inversion method provides a feasible technical means for detecting data quality in the field and rapidly acquiring inversion interpretation results.
Referring to fig. 1, a method for rapidly inverting the marine controllable source electromagnetic wave with anisotropic resistivity of a layered medium mainly comprises the following steps:
s1, reading and converting ocean controllable source electromagnetic field data participating in inversion, wherein the data comprises real part data and imaginary part data of an electromagnetic field, amplitude data and phase data, and parameters of a long axis and a short axis of a polarization ellipse.
According to user requirements, converting the marine controllable source electromagnetic data file participating in inversion into a data file with a specific format, converting the input electromagnetic field data into specific parameters (real part and imaginary part of an electromagnetic field, amplitude and phase, polarization ellipse major axis and minor axis parameters and the like) according to the user requirements, completing reading of inversion data and starting inversion.
S2, inversion execution parameters are set, wherein the parameters comprise inversion maximum iteration times, target fitting difference, maximum iteration step length, step length proportion coefficient, penalty function type and regularization attenuation coefficient.
The inversion process is a complex and huge operation system, and a large number of parameters need to be applied in the process of realizing inversion iteration until a final inversion model is obtained. Different inversion execution parameters can obtain different inversion effects, and a user can set the inversion execution parameters according to the characteristics of the data participating in inversion and the specific requirements needed to be used. The inversion parameters mainly involved include: the method comprises the following steps of inversion maximum iteration times, target fitting difference, maximum iteration step length, step length proportionality coefficient, penalty function type, data components participating in inversion, regularized attenuation coefficient and the like.
S21, inversion maximum iteration times: maximum value of inversion iteration number. The inversion program can be exited after the inversion process reaches the maximum iteration number by setting the parameters, continuous cycle operation of the inversion program is avoided, and the program calculation efficiency is improved.
S22, target fitting difference: and the target value of the fitting degree of the inversion model and the real model. The final goal of inversion is to make the inversion model approach the real model continuously, and the poor fitting of data is to determine the fitting standard of the inversion model and the real model. Setting a target fitting difference by a user according to the data quality, and when the fitting difference of the inversion model reaches or is lower than the target fitting difference, exiting the inversion; if the fitting difference of the inversion model does not reach the target fitting difference, the inversion procedure will continue.
S23, maximum iteration step length: and setting the maximum value of the variable quantity of the model parameters in the inversion process by a user. In the inversion, because the magnitude span of the ocean controllable source electromagnetic field (or other parameters) is large, the inversion operation is usually performed after logarithms are taken. If the iterative model variation is too large, the result is changed very greatly, for example, the variation of the model parameter (after taking the logarithm) is 1.0, i.e., the true value is increased by an order of magnitude, and such a variation is usually unreasonable, so the parameter is designed to limit the variation of the model parameter.
S24, step length proportion coefficient: and setting a proportionality coefficient of the model parameter transformation quantity in the inversion process by a user. Inversion involves many inversion parameters, and in the inversion iteration process, the variation degrees of different inversion parameters are different, so that the variation of the model parameter can be multiplied by a ratio coefficient (i.e. step-size ratio coefficient) to ensure the steady convergence of the inversion model in order to make the overall model variation in a reasonable range.
S25, penalty function type: the manner of constraints between the model parameters. In the inversion process, all model parameters participate in inversion at the same time, certain relation exists among different inversion parameters, restriction modes among the model parameters can be realized by selecting different penalty functions, and the model can be smoother or the curvature of the model can be increased.
S26, regularized attenuation coefficient: the attenuation coefficient of the regularization factor is adaptively selected. In the inversion iteration process, in order to ensure the ratio relationship between different parameters, the regularization factor needs to be changed correspondingly with the change of the inversion times, and the parameter is the attenuation coefficient of the regularization factor with the iteration times.
And S3, setting parameters of an inversion initial model, wherein the parameters comprise a background layer resistivity parameter, a thickness parameter and an observation system parameter.
The initial inversion model is an initial value of the inversion model parameters, and usually needs to be set according to an observation system of input data and known water depth data and other data.
And S4, constructing an anisotropic resistivity inversion target function.
In view of the difference between different resistivity anisotropy parameters of the seabed medium and the stratum thickness, the inversion process of the ocean controllable source electromagnetic seabed anisotropy resistivity and the stratum thickness parameter is stabilized by utilizing regularization constraint. The objective function used was:
Figure BDA0002618041550000091
in the formula, phi is the target function of the inversion algorithm; m is a model inversion parameter vector including a sea floor anisotropic resistivity parameter m ρ And the thickness m of the formation h I.e. m = m ρ +m h
Figure BDA0002618041550000092
Is the gradient of the model parameter vector; i | · | | is a standard deviation operator; d is an observation data vector used for inversion; w d Weighting the matrix for the data; w m Weighting the model with a matrix; f (m) represents a forward response operator of the model m; mu.s ρ And mu h Respectively a seabed anisotropic resistivity parameter m in an inversion model ρ And the thickness m of the formation h The regularization factor of (1).
Inversion parameter m of seabed resistivity of inversion model ρ Has the following forms
Figure BDA0002618041550000104
In the formula, ρ h And ρ v Respectively the transverse resistivity and the vertical resistivity of the seabed stratum, and M is the number of layers of the seabed stratum.
Thickness m of stratum h Has the following forms
m h =[log 10 h 1 …log 10 h n T
In the formula, n is the layer number of the seabed stratum of the inversion model.
S5, solving a Jacobian matrix and a Hessian matrix of the electromagnetic field relative to the anisotropic resistivity.
S51. Jacobian matrix J i Is a matrix of partial derivatives of the electromagnetic field with respect to inversion parameters, different inversion parameters having different forms:
Figure BDA0002618041550000101
wherein i is the ith inversion iteration; j. the design is a square i A Jacobian matrix that is the forward response F (m); ρ = (ρ) hv ) Is the anisotropic resistivity distribution of the formation, and h is the thickness of the subsea formation.
If the resistivity tensor rho is the transverse resistivity rho of the inversion model h And vertical resistivity ρ v I.e. ρ = f (ρ) hv ) Thus, the Jacobian matrix of the electromagnetic field with respect to the resistivity tensor ρ has the form:
Figure BDA0002618041550000102
s52. Hessian matrix H i For the second derivative of the objective function to the inversion parameters, neglecting the second derivative terms and the asymmetry terms, hessian matrix H i Can be simplified as follows:
Figure BDA0002618041550000103
where i is the ith iteration of inversion, J i A Jacobian matrix that is the forward response F (m); m is a unit of ρ And m h Respectively, a seabed resistivity parameter and a seabed stratum thickness;
Figure BDA0002618041550000105
is the gradient of the model parameter vector; the | is a standard deviation operator; d is an observation data vector used for inversion; w d Weighting the matrix for the data; w m Weighting the model with a matrix; mu.s ρ And mu h Respectively a seabed anisotropic resistivity parameter m in an inversion model ρ And the thickness m of the seabed stratum h The regularization factor of (1).
S6, adaptively calculating a regularization factor based on the inversion parameter characteristics.
The selection mode of the regularization factor mu is the key for solving the problem of reasonability or not of the inversion problem of the electromagnetic seabed anisotropic resistivity of the ocean controllable source and the electromagnetic emission source parameter. In geophysical inversion methods, many regularization factor selection methods have emerged. When the inversion of the electromagnetic seabed anisotropic resistivity and the electromagnetic emission source parameters of the marine controllable source is carried out, different regularization factors are selected according to different inversion parameters. Firstly, selecting characteristic parameters related to a Jacobian matrix in the inversion process as a base number, and then regulating the regularization factors of different inversion parameters by combining the characteristics of the inversion parameters and the relationship among the parameters, thereby realizing the self-adaptive selection of the regularization factors.
Regularization factor mu i The following can be written:
Figure BDA0002618041550000111
wherein i represents the ith iteration of inversion; mu.s i Is a regularization factor; max | is an element for solving the maximum matrix absolute value; a is a mj Is a matrix product [ (W) d J) T (W d J)]An element of (1); m is the matrix product [ (W) d J) T (W d J)]Dimension of (d); χ is the attenuation coefficient; lambda is the inverse model transverse resistivity ρ h Resistivity p in the vertical direction v And a weighting factor for the formation thickness h, determined using the following equation:
Figure BDA0002618041550000112
wherein beta is a proportionality coefficient, m is an inversion model parameter, and the initial value of the weight vector is lambda 1 Is a unit vector. When resistance exists in the stratum of the inverse modelWhen the resistivity anisotropy exists, the regularization factor can be adjusted in a self-adaptive mode according to the resistivity anisotropy rate, so that the weight among elements in the target function is adjusted, and the effect of adjusting the inversion process in a self-adaptive mode is achieved.
And S7, obtaining the model updating amount.
By finding the minimum of the objective function, the model update can be calculated, i.e. the gradient g of the objective function is calculated i =0, can be determined by the following formula:
Figure BDA0002618041550000121
where i is the ith inversion iteration, J i A Jacobian matrix that is the forward response F (m); m is a unit of ρ And m h Respectively, a seabed resistivity parameter and a seabed stratum thickness;
Figure BDA0002618041550000124
is the gradient of the model parameter vector; i | · | | is a standard deviation operator; d is an observation data vector used for inversion; w is a group of d Weighting the matrix for the data; w m A model weighting matrix; mu.s ρ And mu h Respectively a seabed anisotropic resistivity parameter m in an inversion model ρ And the thickness m of the seabed stratum h The regularization factor of (1).
The model parameter update quantity of the next iteration is calculated to be Δ m, which can be expressed as:
Figure BDA0002618041550000122
wherein, Δ m is the model parameter updating amount of the next iteration, i is the ith inversion iteration, H i Is Hessian matrix, g i Is the gradient of the objective function.
And S8, calculating the target function fitting difference of the inversion iteration model.
In the inversion iteration process, the inversion iteration model gradually converges to the real model, and the target function fitting difference psi gradually converges to 1.0. The objective function fitting difference ψ has the following form:
Figure BDA0002618041550000123
wherein, i is inversion iteration times, and d is an observation data vector used for inversion; f (m) represents a forward response operator of the model m; n is the number of inversion data; k is the number of the inversion data; delta. For the preparation of a coating k Is the standard deviation of the kth data.
And S9, judging whether the inversion requirements are met, if so, turning to S9, and if not, turning to S5.
There are two criteria for determining whether to exit inversion: 1) Whether an inversion target fitting difference is satisfied; 2) Whether the inversion iteration number reaches the maximum iteration number or not. If not, turning to S5, and continuing to iteratively solve the model updating amount; if yes, the process goes to step S10.
And S10, outputting a final inversion model.
Example 2
To verify the correctness of the sensitivity matrix, referring to fig. 2, a one-dimensional typical earth model diagram is shown, assuming that the depth of the resistivity isotropic seawater layer is 10km and the resistivity is 0.3 Ω m. The buried depth of the isotropic high-resistance thin layer with the thickness of 100m is 1km, and the resistivity is rho =40 Ω m. Overburden and bedrock resistivities are also isotropic formations, with the resistivities of overburden and bedrock being ρ =2 Ω m and ρ =4 Ω m, respectively. Assuming that the line is coincident with the y-axis, 76 transmission sources are equally spaced within the range of 0m-15000m from the line, all located 50m directly above the seafloor, with a transmission frequency of 0.25Hz, and the receiving station is located at the seafloor with coordinates (x =0m, y =0m, z =10000 m). Taking the one-dimensional resistivity isotropic model shown in fig. 2 as an example, the sensitivity correctness of the resistivity anisotropy analytically calculated by the method is verified.
Referring to fig. 3, sensitivity comparison of the one-dimensional resistivity isotropic model horizontal electric field Ey with respect to the transverse conductivity (fig. 3 a), the vertical conductivity (fig. 3 b) and the isotropic conductivity (fig. 3 c) of the seafloor medium, and relative error (fig. 3 d). Wherein, the three solid lines are respectively the sensitivity of the electric field relative to the conductivity of the covering layer and the conductance of the electric field relative to the high-resistance layerSensitivity of the rate and sensitivity of the electric field with respect to the conductivity of the matrix. It can be seen from the figure that, under the condition of small receiving and transmitting distance, the sensitivity of the covering layer is larger than that of other layers; under the condition of large transmitting-receiving distance, the sensitivity of the electric field on the transverse conductivity of the bedrock is higher than that of the covering layer and the high-resistance layer, and the vertical conductivity of the electric field on the high-resistance layer is higher than that of the other two layers; the vertical conductivity sensitivity is more consistent with the overall trend of isotropic sensitivity, whereas the lateral conductivity is much different from both of them. The electrical resistivity isotropy is a special case of electrical resistivity anisotropy, and the isotropic electrical conductivity sigma can be converted iso As a function of the lateral conductivity σ h And vertical conductivity σ v A function of (a), i.e. sigma iso =diag[σ h σ h σ v ]Where σ is iso =σ h =σ v . Therefore, for a resistivity isotropic earth model, the sensitivity of the electromagnetic field with respect to conductivity is the sum of the electromagnetic field with respect to transverse and vertical conductivities, i.e.
Figure BDA0002618041550000141
Wherein sigma iso =σ h =σ v
The solid line in fig. 3c is the sum of the electric field calculated by this patent with respect to the transverse conductivity and the vertical conductivity, and the cross line is the result calculated by the conductivity isotropy program, and as can be seen from the figure, both errors are very small, both being less than 0.1% (fig. 3 d). It is thus demonstrated that the sensitivity calculated by this patent is accurate.
Example 3
This example is resistivity isotropic half-space model inversion. Because the underground medium in the half-space model has no layering condition, the influence of the stratum thickness on the inversion can be removed. Two half-space models are established, and the seabed half-space is assumed to be a resistivity isotropic medium, and the resistivities are 2 ohm m and 40 ohm m respectively. In real geophysical inversion, it is often not known how many layers of media are on the seafloor, and the number of formation layers is usually set to be greater than the number of real formation layers. Therefore, the seabed of the initial inversion model is set to be 5 layers, and the resistivity and the layer thickness of the seabed medium formation (resistivity isotropic inversion) and the transverse resistivity, the vertical resistivity and the layer thickness (resistivity anisotropic inversion) are respectively inverted.
In the axial observation mode, a horizontal electric dipole source with y-orientation generates a horizontal magnetic field component B x And two electric field components: a horizontal electric field component E y And a vertical electric field component E z The real and imaginary components of the three electromagnetic field components are used as the inversion data of the composite model.
Referring to table 1, the inversion results of resistivity isotropy and resistivity anisotropy under the condition that the seabed half space is 2 Ω m are shown. Referring to table 2, the resistivity isotropy and resistivity anisotropy inversion results are obtained for the case where the seafloor half space is 40 Ω m. As can be seen from tables 1 and 2, the fitting difference of inversion data of resistivity isotropy and resistivity anisotropy finally converges to be near 1, although the number of layers of the seabed medium is set to be more than that of the real seabed stratum in the inversion initial model, the resistivity and the stratum thickness are inverted simultaneously, and the inversion result still conforms to the resistivity distribution condition of the real seabed medium.
TABLE 1 resistivity Isotropic half-space model inversion results (2. Omega. M)
Figure BDA0002618041550000151
TABLE 2 resistivity Isotropic half-space model inversion results (40. Omega. M)
Figure BDA0002618041550000161
Example 4
This example performs conductivity anisotropy inversion for conductivity anisotropy data. Referring to fig. 4, a resistivity anisotropy model is shown. The seawater layer at a depth of 10km was assumed to be isotropic in resistivity, with a resistivity of 0.3 Ω m. The buried depth of the high-resistance thin layer with the vertical anisotropy and the thickness of 100m is 1km, and the transverse resistivity is rho h =20 Ω m, verticalResistivity is rho v =80 Ω m; the resistivity of the overburden and the bedrock also have vertical anisotropy, and the transverse resistivity and the vertical resistivity of the overburden are respectively rho h =1 Ω m and ρ v =4 Ω m, and the transverse resistivity and the vertical resistivity of the bedrock are ρ n =2 Ω m and ρ v =8 Ω m. Assuming that the survey line coincides with the y-axis, the transmitting source is located in the sea water 50m from the sea bottom, 76 receiving stations are arranged at equal intervals in the range of 0m-15000m from the survey line, and the transmitting frequency varies in the range of 0.01-100 Hz. In the inversion, the inversion initial model is the same as the inversion example, the sea water resistivity and the depth are fixed, and only the conductivity and the stratum thickness of the seabed medium are inverted.
Referring to fig. 5a, a variation curve of difference versus inversion iteration number is fitted to the resistivity anisotropy inversion data of the resistivity anisotropy model forward response data. As can be seen, the inversion eventually converges through nearly 80 iterations. Referring to fig. 5b, resistivity isotropic inversion results, wherein a dark dotted line is a lateral resistivity of the real model, a light dotted line is a vertical resistivity of the real model, dark and light banding are implemented as inversion iteration results of conductivity isotropic inversion from 60 th to 80 th times, a solid line with lighter color is an inversion result with less iteration times, and a solid line with darker color is an inversion result with more iteration times. As can be seen from FIG. 5, the transverse resistivity and the vertical resistivity of the overburden and the bedrock are consistent with the real values, the buried depth, the thickness and the vertical resistivity of the high-resistivity thin layer are very close to the real values, and the transverse resistivity of the inverted high-resistivity layer is lower than the real values. As can be seen from the results of 60-80 inversion iterations, the resistivity of each inversion layer hardly changes at the later stage of inversion, but the thickness of each layer is slightly adjusted.
Referring to fig. 6, a schematic diagram is fitted to observed data of the Ey component (fig. 6 a), bx component (fig. 6 b) and Ez component (fig. 6 c) and data of the forward response of the inversion model, which participate in the resistivity anisotropy inversion. In the figure, the "+" line and the "o" line are respectively the comparison between the observed data and the forward response of the model, and the lower data fitting condition (Misfit) is the fitting difference D between the observed data and the forward response of the inverse model at each transceiving distance d Of the ith dataThe fitting difference is defined as follows
Figure BDA0002618041550000171
Wherein, F mod And F obs Respectively are forward modeling response and observation data of the inversion model, and delta is a standard deviation. When the inverse model is identical to the real model, the data fitting difference D d Equal to 1.0. As can be seen from FIG. 6, except for the data with the electromagnetic field amplitude smaller than the background noise, the forward response of the inverse model and the observed data are well fitted, and the data fitting difference D d Are distributed around 1.0. The inversion result shows that the ocean controllable source resistivity anisotropy inversion provided by the method can better reconstruct the ocean bottom anisotropy resistivity distribution.
The above description is intended to be illustrative of the preferred embodiment of the present invention and should not be taken as limiting the invention, but rather, the intention is to cover all modifications, equivalents, and alternatives falling within the spirit and scope of the invention.

Claims (3)

1. A layered medium resistivity anisotropy ocean controllable source electromagnetic fast inversion method is characterized by comprising the following steps:
s1, reading and converting ocean controllable source electromagnetic field data participating in inversion, wherein the data comprises real part data and imaginary part data of an electromagnetic field, amplitude data and phase data, and parameters of a long axis and a short axis of a polarization ellipse;
s2, inversion execution parameters are set, wherein the parameters comprise inversion maximum iteration times, target fitting difference, maximum iteration step length, step length proportionality coefficient, penalty function type and regularization attenuation coefficient;
s3, setting inversion initial model parameters, wherein the parameters comprise a background layer resistivity parameter, a thickness parameter and an observation system parameter;
s4, constructing an anisotropic resistivity inversion target function,
the anisotropic resistivity inversion objective function is:
Figure FDA0004049273990000011
wherein φ is the objective function of the inversion algorithm; m is a model inversion parameter vector including a sea floor anisotropic resistivity parameter m ρ And the thickness m of the formation h I.e. m = m ρ +m h
Figure FDA0004049273990000012
Is the gradient of the model parameter vector; i | · | | is a standard deviation operator; d is the observation data vector used for inversion; w d Weighting the matrix for the data; w is a group of m A model weighting matrix; f (m) represents a forward response operator of the model m; mu.s ρ And mu h Respectively a seabed anisotropic resistivity parameter m in an inversion model ρ And the thickness m of the formation h A regularization factor of (a);
s5, solving a jacobian matrix and a hessian matrix of the electromagnetic field relative to the anisotropic resistivity; the expression of the Jacobian matrix is as follows:
Figure FDA0004049273990000013
wherein i is the ith inversion iteration; ji is the Jacobian matrix of the forward response F (m); ρ = (ρ h, ρ v) is anisotropic resistivity distribution of the formation, h is the thickness of the subsea formation;
the expression of the Hessian matrix is as follows:
Figure FDA0004049273990000014
wherein i is the ith inversion iteration, and Ji is the Jacobian matrix of forward response F (m); the mp and mh are respectively a seabed resistivity parameter and a seabed stratum thickness;
Figure FDA0004049273990000015
is the gradient of the model parameter vector; i | · | | is a standard deviation operator; d is an observation data vector used for inversion; wd is a data weighting matrix; wm is a model weighting matrix; mu rho and mu h are regularization factors of a seabed anisotropic resistivity parameter m rho and a seabed stratum thickness mh in the inversion model respectively;
s6, adaptively calculating regularization factors based on inversion parameter characteristics, wherein the regularization factors of different inversion parameters are adjusted by selecting characteristic parameters related to a Jacobian matrix in the inversion process as a base number and combining the characteristics of the inversion parameters and the relationship among the parameters; the adaptive regularization factor is determined according to the following formula:
Figure FDA0004049273990000021
wherein i represents the ith iteration of inversion; mu.s i Is a regularization factor; max | is the element with the maximum matrix absolute value; a is a mj Is a matrix product [ (W) d J) T (W d J)]The element (b); m is a matrix product [ (W) d J) T (W d J)]Dimension (d); χ is the attenuation coefficient; lambda is the transverse resistivity rho of the submarine stratum of the inversion model h Resistivity p in the vertical direction v And a weighting factor for the formation thickness h, determined using the following formula,
Figure FDA0004049273990000022
wherein gamma is a proportionality coefficient, and m is an inversion model parameter;
s7, calculating a model updating amount, enabling the gradient gi of the target function to be =0, and calculating the model updating amount;
s8, calculating a target function fitting difference of the inversion iteration model;
s9, judging whether inversion requirements are met, if so, turning to S10, and if not, turning to S5; the criteria for judging whether the inversion requirements are met include: 1) Whether an inversion target fitting difference is satisfied; 2) Whether the inversion iteration times reach the maximum iteration times or not;
and S10, outputting the final inversion model.
2. The method for rapidly inverting the marine controllable source electromagnetism of laminar medium resistivity anisotropy according to claim 1, wherein the method for calculating the model parameter update amount based on the inversion parameters comprises the following steps:
Figure FDA0004049273990000023
wherein, Δ m is the model parameter updating amount of the next iteration, i is the ith inversion iteration, H i Is Hessian matrix, g i Is the gradient of the objective function.
3. The method for rapidly inverting the marine controllable source electromagnetism of the layered dielectric resistivity anisotropy according to claim 2, wherein the method for solving the minimum value of the objective function is as follows:
Figure FDA0004049273990000031
where i is the ith iteration of inversion, J i Is a jacobian matrix.
CN202010774885.6A 2020-08-05 2020-08-05 Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method Active CN111856596B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010774885.6A CN111856596B (en) 2020-08-05 2020-08-05 Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010774885.6A CN111856596B (en) 2020-08-05 2020-08-05 Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method

Publications (2)

Publication Number Publication Date
CN111856596A CN111856596A (en) 2020-10-30
CN111856596B true CN111856596B (en) 2023-03-28

Family

ID=72953347

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010774885.6A Active CN111856596B (en) 2020-08-05 2020-08-05 Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method

Country Status (1)

Country Link
CN (1) CN111856596B (en)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112578470B (en) * 2020-11-10 2022-07-29 中国海洋大学 Product function-based ocean controllable source electromagnetic and magnetotelluric joint inversion method
CN112733449B (en) * 2021-01-11 2022-12-02 中国海洋大学 CNN well-seismic joint inversion method, CNN well-seismic joint inversion system, CNN well-seismic joint inversion storage medium, CNN well-seismic joint inversion equipment and CNN well-seismic joint inversion application
CN113190957B (en) * 2021-03-24 2024-03-22 中国海洋大学 Controllable source electromagnetic simulation wave number sequence optimization method based on elimination strategy
CN115407423B (en) * 2021-05-27 2024-02-06 中国自然资源航空物探遥感中心 Three-dimensional inversion method and device for weight and magnetic measurement data
CN118033764B (en) * 2024-04-12 2024-06-07 中国地质大学(武汉) Multi-scale multi-method three-dimensional imaging method and system for resistivity of geologic body

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2008013613A2 (en) * 2006-07-25 2008-01-31 Exxonmobil Upstream Research Company Method for determining physical properties of structures
US20090083006A1 (en) * 2007-09-20 2009-03-26 Randall Mackie Methods and apparatus for three-dimensional inversion of electromagnetic data
CN104375195B (en) * 2013-08-15 2017-03-15 中国石油天然气集团公司 Many source multi-component three-dimensional joint inversion methods of time-frequency electromagnetism
CN106019394B (en) * 2016-04-27 2019-04-05 中国地质科学院矿产资源研究所 Three-dimensional parallel inversion method for nonlinear conjugate gradient of ocean magnetotelluric field
US10871590B2 (en) * 2017-06-16 2020-12-22 Pgs Geophysical As Electromagnetic data inversion
CN111103628B (en) * 2020-01-14 2022-03-11 桂林理工大学 Two-dimensional inversion method and device for electric field data in magnetotelluric (TE) polarization mode
CN111103627B (en) * 2020-01-14 2022-03-11 桂林理工大学 Two-dimensional inversion method and device for electric field data by magnetotelluric (TM) polarization mode

Also Published As

Publication number Publication date
CN111856596A (en) 2020-10-30

Similar Documents

Publication Publication Date Title
CN111856596B (en) Layered medium resistivity anisotropy ocean controllable source electromagnetic rapid inversion method
Colombo et al. Physics-driven deep-learning inversion with application to transient electromagnetics
Key 1D inversion of multicomponent, multifrequency marine CSEM data: Methodology and synthetic studies for resolving thin resistive layers
Grayver et al. 3D inversion and resolution analysis of land-based CSEM data from the Ketzin CO 2 storage formation
Persova et al. Three‐dimensional inversion of airborne data with applications for detecting elongated subvertical bodies overlapped by an inhomogeneous conductive layer with topography
Zhdanov et al. Anisotropic 3D inversion of towed-streamer electromagnetic data: Case study from the Troll West Oil Province
CN111666534B (en) Electrical random anisotropic electromagnetic field decoupling method
Jaysaval et al. Fast multimodel finite-difference controlled-source electromagnetic simulations based on a Schur complement approach
Nguyen et al. Comparing large-scale 3D Gauss–Newton and BFGS CSEM inversions
Chen et al. Integration of principal-component-analysis and streamline information for the history matching of channelized reservoirs
Amaya et al. A low-rank approximation for large-scale 3D controlled-source electromagnetic Gauss-Newton inversion
Li et al. Kirchhoff migration using eikonal-based computation of traveltime source derivatives
Başokur Automated 1D interpretation of resistivity soundings by simultaneous use of the direct and iterative methods
CN111880235B (en) Ocean electromagnetic formation anisotropic resistivity and emission source posture joint inversion method
Matveev et al. Geology realism control in automated history matching
Olalekan et al. Particle swarm optimization method for stochastic inversion of MTEM data
Li et al. 3-D marine CSEM forward modeling with general anisotropy using an adaptive finite-element method
Colombo et al. Machine-learning inversion via adaptive learning and statistical sampling: Application to airborne micro-TEM for seismic sand corrections
CN111856597B (en) Towed marine electromagnetic formation resistivity and receiving station position joint inversion method
CN110764154A (en) Time domain aviation electromagnetic one-dimensional inversion method based on improved particle swarm optimization
Carcione A spectral numerical method for electromagnetic diffusion
Bhuyian et al. Controlled source electromagnetic three‐dimensional grid‐modelling based on a complex resistivity structure of the seafloor: effects of acquisition parameters and geometry of multi‐layered resistors
Wirianto et al. Applying essentially non‐oscillatory interpolation to controlled‐source electromagnetic modelling
Kumar et al. Regularization analysis of three-dimensional magnetotelluric inversion
Luan et al. Ground-wire source TEM 3D full time multinary inversion using adaptive regulation

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant