CN113917560A - Three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method - Google Patents
Three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method Download PDFInfo
- Publication number
- CN113917560A CN113917560A CN202111086356.8A CN202111086356A CN113917560A CN 113917560 A CN113917560 A CN 113917560A CN 202111086356 A CN202111086356 A CN 202111086356A CN 113917560 A CN113917560 A CN 113917560A
- Authority
- CN
- China
- Prior art keywords
- inversion
- dimensional
- gravity
- grid
- imaging
- 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.)
- Granted
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V11/00—Prospecting or detecting by methods combining techniques covered by two or more of main groups G01V1/00 - G01V9/00
- G01V11/007—Prospecting or detecting by methods combining techniques covered by two or more of main groups G01V1/00 - G01V9/00 using the seismo-electric effect
-
- Y—GENERAL 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
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention discloses a three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method, which comprises the following steps: s1, acquiring seismic and electromagnetic data of the area to be detected; s2, performing magnetotelluric two-three dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter co-imaging and heavy seismic interface joint inversion according to the seismic and electromagnetic data; s3, performing the combined collaborative imaging of the gravity electromagnetic and seismoelectric physical parameters according to the results of magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion and gravity electromagnetic multi-parameter collaborative imaging; and S4, constructing an imaging result for geological comprehensive interpretation by combining the gravity magnetic seismographic physical property parameter combined collaborative imaging result and the gravity-seismographic interface combined inversion result, and realizing three-dimensional gravity magnetic seismographic multi-parameter collaborative inversion.
Description
Technical Field
The invention belongs to the technical field of geological structure imaging, and particularly relates to a three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method.
Background
Seismic exploration has been the main means of oil and gas exploration and development, mainly due to: the method comprises the following steps that firstly, underground complex geological structure imaging can be obtained by utilizing a seismic imaging technology, and especially, underground structure imaging with higher resolution can be provided by Reverse Time Migration (RTM) and Full Waveform Inversion (FWI); and secondly, extracting physical property parameter information of the reservoir, such as porosity, fluid saturation, permeability and the like, by using a seismic reservoir inversion technology, and further predicting the type of reservoir fluid and evaluating the oil-gas content of the reservoir. However, under complex geological conditions, such as igneous rock areas, carbonate reef development areas, reverse masked fracture zones and the like, the illumination energy of the underlying stratum is weak due to the severe seismic wave scattering effect and shielding effect, and the multi-solution of the underground structure is inferred only by using seismic data. In terms of reservoir description and monitoring, in the case of low velocity difference, such as stratum water and oil, reservoirs with different gas saturation, the reservoir inversion using seismic data has certain ambiguity, because the influence of the change of the oil (gas) saturation on the seismic wave velocity is not obvious, and therefore, the reservoir fluid type and the fluid saturation cannot be well predicted by using the seismic data alone.
Magnetotelluric (MT) methods are particularly important for solving the problem of oil and gas exploration under complex geological conditions (e.g., igneous rock underburden imaging, imaging under salt, etc.), mainly because: firstly, the attenuation of the electromagnetic signal is very weak when the electromagnetic signal passes through a high-resistance body (such as igneous rock), and secondly, the electromagnetic signal is sensitive to a low-resistance sedimentary stratum. Thus, the electromagnetic signals collected contain resistivity information of the underlying sedimentary formation. However, the MT method is not well able to image relatively thin hydrocarbon reservoirs because it is mainly based on low frequency induced electromagnetic fields. The basic equations governing the propagation of seismic and electromagnetic waves are the wave equation and the diffusion equation, respectively, and although both waves propagate with attenuation and scattering properties, the attenuation and scattering degrees are different, i.e., the attenuation and scattering of seismic waves are relatively weak.
Disclosure of Invention
Aiming at the defects in the prior art, the three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method provided by the invention solves the problems of single method multiple solution and low resolution limit in the existing inversion method.
In order to achieve the purpose of the invention, the invention adopts the technical scheme that: a three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method comprises the following steps:
s1, acquiring seismic and electromagnetic data of the area to be detected;
s2, performing magnetotelluric two-three dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter co-imaging and heavy seismic interface joint inversion according to the seismic and electromagnetic data;
s3, performing the combined collaborative imaging of the gravity electromagnetic and seismoelectric physical parameters according to the results of magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion and gravity electromagnetic multi-parameter collaborative imaging;
and S4, constructing an imaging result for geological comprehensive interpretation by combining the gravity magnetic seismographic physical property parameter combined collaborative imaging result and the gravity-seismographic interface combined inversion result, and realizing three-dimensional gravity magnetic seismographic multi-parameter collaborative inversion.
Further, in step S2, the performing magnetotelluric two-dimensional and three-dimensional inversion includes the following sub-steps:
a1, performing magnetotelluric two-dimensional inversion and dip inversion of the spherical column coordinates according to the seismic and electromagnetic data;
a2, performing structured self-adaptive finite element two-dimensional magnetotelluric inversion and modeling according to the seismic and electromagnetic data;
a3, taking the result of the two-dimensional and three-dimensional magnetotelluric inversion of the spherical column coordinates and the result of the two-dimensional and three-dimensional magnetotelluric inversion and modeling of the structured self-adaptive finite elements as the result of the two-dimensional and three-dimensional magnetotelluric inversion;
wherein, step A1 includes the following substeps:
a11, constructing an integrated form of Maxwell equation satisfied by the electric field intensity vector and the magnetic field intensity vector based on the acquired electromagnetic data:
∮H·dl＝∫∫J·dS＝∫∫σE·ds
∮E·dl＝∫∫iωμ0·HdS
wherein E and H are respectively an electric field intensity vector and a magnetic field intensity vector, J is a current density, dl is a contour of a closed integral, dS is an area included by the contour, sigma is a medium conductivity, omega is a circular frequency, and mu0Is a vacuum magnetic conductivity;
a12, dividing the spherical coordinate system into a plurality of grid units by a plurality of parallel spherical surfaces at different intervals in the r, theta and phi directions in the spherical coordinate system, and realizing the discretization of the spherical coordinate system according to a staggered grid subdivision mode;
a13, discretizing and expanding the Maxwell equation integral form in a discretized spherical coordinate system, further determining a staggered grid finite difference equation, and taking the staggered grid finite difference equation as the results of magnetotelluric two-dimensional inversion and dip inversion of the spherical coordinate system;
step a2 includes the following substeps:
a21, for earthquake and electromagnetic data, sequentially processing the earthquake and electromagnetic data by adopting a geoelectromagnetic method self-adaptive vector finite element method, a linear equation set solving method based on an auxiliary space precondition algorithm, a parallel method based on a grid partition technology and a geoelectromagnetic three-dimensional inversion method based on self-adaptive finite element forward modeling;
a22, processing the result obtained by the magnetotelluric three-dimensional inversion method based on the forward modeling of the self-adaptive finite element as the result of the magnetotelluric inversion and modeling of two three-dimensional earth of the structural self-adaptive finite element;
in step a21, the magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling specifically includes the following sub-steps:
t1, obtaining data of the region to be inverted according to the grid partitioning technology-based parallel method, and determining the region to be inverted according to the data;
t2, discretizing the region to be inverted, and taking the discretized grid as an inversion grid;
wherein the inversion grid does not change during the inversion process;
t3, setting a forward net for each frequency in each iteration of the inversion process;
wherein the initial value of the forward grid is the same as the initial value of the inversion grid;
t4, in the iteration of the inversion process, self-adaptive grid encryption is carried out on the forward grid for a plurality of times by using the posterior error estimated value, and a partial derivative is calculated on the final forward grid;
and T5, mapping the partial derivative calculated on the final forward grid back to the corresponding inversion grid, and realizing magnetotelluric three-dimensional inversion.
Further, in step S2, performing three-dimensional quantitative inversion of apparent density and apparent magnetization includes the following sub-steps:
b1, acquiring a gravity magnetic anomaly data body of the whole space of the region to be detected by utilizing a regularization downward extension technology;
b2, changing point constraint or line constraint into volume constraint based on a drilling constraint technology, and searching out a grid interval in which a gravity magnetic anomaly data volume exists by adopting a probability imaging technology to realize dimension reduction of a solution space;
b3, performing three-dimensional quantitative inversion of apparent density and apparent magnetization intensity by adopting a multi-method inversion system for the inversion target function of the solution space after dimensionality reduction;
the method in the multi-method inversion system comprises a golden section inversion method based on correlation search, a three-dimensional gravity inversion method based on rapid near-end objective function optimization and a rock physical relationship constraint inversion method based on machine learning and multivariate statistics.
Further, the rock physics relationship constraint inversion method based on machine learning and multivariate statistics in the step B3 includes the following sub-steps:
r1, establishing a cross variation function reflecting the spatial statistical rock physical relationship by taking the logging data and the seismic frame as constraint conditions;
r2, clustering the velocity and density distribution obtained by independent seismic inversion and gravity inversion respectively, and redistributing the velocity value and the density value to the clustering result according to the statistical result of the rock sample data and the logging data to obtain a new velocity and density model;
and R3, taking the cross variation function as a prior function of an inversion target function, and performing gravity inversion of rock physical relationship constraint based on the new speed and density.
Further, in the step S2, performing the gravity magnetoelectric multi-parameter collaborative imaging including resistivity imaging and multi-parameter imaging;
the resistivity imaging method specifically comprises the following steps: performing joint inversion on the three-dimensional controllable electromagnetic data and the seismic data by adopting a cross gradient operator to obtain resistivity imaging;
the multi-parameter imaging method specifically comprises the following steps:
y1, performing numerical simulation of the seismic wave field and the electromagnetic wave field by adopting a staggered grid finite difference method based on a seismic wave equation and a Maxwell equation;
y2, constructing a cross gradient-based electromagnetic and seismic corresponding joint inversion objective function based on a data fitting term, a model regularization term and a cross gradient constraint of electromagnetic data and seismic data obtained by numerical simulation;
y3, solving the constructed joint inversion target function by adopting a linear conjugate gradient algorithm;
and Y4, performing electromagnetic and seismic joint inversion based on a joint inversion objective function solving process, and realizing multi-parameter imaging.
Further, in the step Y4, when performing joint inversion of the electromagnetic and seismic, the method adopted sequentially includes: and adjusting regularization factors by adopting a self-adaptive regularization method through data fitting based on two adjacent iterative processes, carrying out self-adaptive correction based on K-means clustering and regression analysis, carrying out cross gradient weighting, carrying out frequency division multi-scale inversion and multi-grid division, and setting physical constraint range inversion.
Further, in step S2, the method for performing the joint inversion of the heavy seismic interface specifically includes:
and performing underground density interface inversion by combining a random inversion method and a cross gradient inversion method with gravity and seismic data.
Further, in step S3, the method for performing the gravity magnetic and electric shock physical parameter joint collaborative imaging specifically includes:
based on the comprehensive response characteristics of the gravity, magnetism, electricity and earthquake in the gravity, magnetism, electricity and earthquake physical property parameter joint collaborative imaging result and the gravity and earthquake interface joint inversion result, the gravity, magnetism, electricity and earthquake modeling-inversion is carried out, the physical property parameters of different geophysical fields are inverted, and the physical property parameters of the different geophysical fields of the gravity, magnetism, electricity and earthquake are utilized to carry out mutual constraint inversion collaborative imaging.
The invention has the beneficial effects that:
(1) the invention innovatively constructs heavy-shock, electric-shock and heavy-magnetic-electric-shock collaborative inversion and imaging technologies based on multi-parameter modeling and constraint of density, speed, resistivity, magnetic susceptibility and the like;
(2) the method is applied to the exploration of oil-gas resources and mineral resources, achieves outstanding application effects, is based on the rock physics internal relation, introduces a new algorithm and a new simulation technology, forms a heavy-seismoelectric interface joint inversion technology and a heavy magnetoelectroseismic multi-parameter collaborative simulation imaging technology, and realizes high-precision imaging of various physical parameters (density, resistivity and speed);
(3) by adopting a machine learning method, a multi-parameter joint collaborative simulation imaging technology of density, magnetization intensity, resistivity and the like is developed, the spanning from a single method to a multi-method multi-dimensional collaborative interpretation technology is realized, and the method is an innovative technology for extracting and enhancing the heavy magnetic and electric weak information combination of a geological target body.
Drawings
FIG. 1 is a flow chart of a three-dimensional gravity magnetoseis multi-parameter collaborative inversion method provided by the invention.
Detailed Description
The following description of the embodiments of the present invention is provided to facilitate the understanding of the present invention by those skilled in the art, but it should be understood that the present invention is not limited to the scope of the embodiments, and it will be apparent to those skilled in the art that various changes may be made without departing from the spirit and scope of the invention as defined and defined in the appended claims, and all matters produced by the invention using the inventive concept are protected.
As shown in fig. 1, a three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method includes the following steps:
s1, acquiring seismic and electromagnetic data of the area to be detected;
s2, performing magnetotelluric two-three dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter co-imaging and heavy seismic interface joint inversion according to the seismic and electromagnetic data;
s3, performing the combined collaborative imaging of the gravity electromagnetic and seismoelectric physical parameters according to the results of magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion and gravity electromagnetic multi-parameter collaborative imaging;
and S4, constructing an imaging result for geological comprehensive interpretation by combining the gravity magnetic seismographic physical property parameter combined collaborative imaging result and the gravity-seismographic interface combined inversion result, and realizing three-dimensional gravity magnetic seismographic multi-parameter collaborative inversion.
In step S2 of the present embodiment, performing two-dimensional and three-dimensional magnetotelluric inversion includes the following sub-steps:
a1, performing magnetotelluric two-dimensional inversion and dip inversion of the spherical column coordinates according to the seismic and electromagnetic data;
a2, performing structured self-adaptive finite element two-dimensional magnetotelluric inversion and modeling according to the seismic and electromagnetic data;
a3, taking the result of the two-dimensional and three-dimensional magnetotelluric inversion of the spherical column coordinates and the result of the two-dimensional and three-dimensional magnetotelluric inversion and modeling of the structured self-adaptive finite elements as the result of the two-dimensional and three-dimensional magnetotelluric inversion;
wherein, step A1 includes the following substeps:
a11, constructing an integrated form of Maxwell equation satisfied by the electric field intensity vector and the magnetic field intensity vector based on the acquired electromagnetic data:
∮H·dl＝∫∫J·dS＝∫∫σE·ds
∮E·dl＝∫∫iωμ0·HdS
wherein E and H are respectively an electric field intensity vector and a magnetic field intensity vector, J is a current density, dl is a contour of a closed integral, dS is an area included by the contour, sigma is a medium conductivity, omega is a circular frequency, and mu0Is a vacuum magnetic conductivity;
a12, dividing the spherical coordinate system into a plurality of grid units by a plurality of parallel spherical surfaces at different intervals in the r, theta and phi directions in the spherical coordinate system, and realizing the discretization of the spherical coordinate system according to a staggered grid subdivision mode;
a13, discretizing and expanding the Maxwell equation integral form in a discretized spherical coordinate system, further determining a staggered grid finite difference equation, and taking the staggered grid finite difference equation as the results of magnetotelluric two-dimensional inversion and dip inversion of the spherical coordinate system;
in particular, in the frequency range of magnetotelluric studies (about 10)5～10-5Hz), because the frequency is relatively reduced, the quasi-stable electromagnetic field is satisfied, the displacement current can be ignored, and the magnetic permeability of the underground medium is taken asVacuum magnetic conductivity, and constructing a corresponding Maxwell equation integral form; meanwhile, because the discretization of the difference form of the Max equation set in the spherical coordinate system is complex, the linear equation is obtained by discretizing the integral form of the Max equation set, and therefore mesh generation is needed before discretization.
The step a2 includes the following sub-steps:
a21, for earthquake and electromagnetic data, sequentially processing the earthquake and electromagnetic data by adopting a geoelectromagnetic method self-adaptive vector finite element method, a linear equation set solving method based on an auxiliary space precondition algorithm, a parallel method based on a grid partition technology and a geoelectromagnetic three-dimensional inversion method based on self-adaptive finite element forward modeling;
a22, processing the result obtained by the magnetotelluric three-dimensional inversion method based on the forward modeling of the self-adaptive finite element as the result of the magnetotelluric inversion and modeling of two three-dimensional earth of the structural self-adaptive finite element;
specifically, in step a21, in the magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling, a model established by three-dimensional forward modeling is adopted, the model established by forward modeling and forward modeling is firstly used for adaptively encrypting the grid by a magnetotelluric adaptive vector finite element method, the adaptive grid encryption can automatically encrypt an area with a large error according to posterior error estimation, and the error drop speed is significantly faster than that of the global grid encryption compared with that of the global grid encryption; solving the inverted target function by a linear equation system solving method of an auxiliary space preconditions algorithm; the parallel method of the grid partitioning technology distributes computing tasks to a plurality of processes according to grid partitioning, and each process uses MPI to communicate in the computing process to keep synchronization. The inversion is realized by using a forward algorithm, and the forward and inverse grid separation is carried out, so that the problem of inversion accuracy is solved.
Based on this, in the step a21, the magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling specifically includes the following sub-steps:
t1, obtaining data of the region to be inverted according to the grid partitioning technology-based parallel method, and determining the region to be inverted according to the data;
t2, discretizing the region to be inverted, and taking the discretized grid as an inversion grid;
wherein the inversion grid does not change during the inversion process;
t3, setting a forward net for each frequency in each iteration of the inversion process;
wherein the initial value of the forward grid is the same as the initial value of the inversion grid;
t4, in the iteration of the inversion process, self-adaptive grid encryption is carried out on the forward grid for a plurality of times by using the posterior error estimated value, and a partial derivative is calculated on the final forward grid;
and T5, mapping the partial derivative calculated on the final forward grid back to the corresponding inversion grid, and realizing magnetotelluric three-dimensional inversion.
In step S2 of the present embodiment, the performing of the three-dimensional quantitative inversion of apparent density and apparent magnetization includes the following sub-steps:
b1, acquiring a gravity magnetic anomaly data body of the whole space of the region to be detected by utilizing a regularization downward extension technology;
in this step, the downward continuation of the bit field is a typical ill-posed problem, and due to errors and some interference information in the actual observed value, a high-frequency oscillation effect occurs when the actual observed value is descended to the vicinity of the field source, and useful information is often covered. The regularization method is one of the most stable algorithms, and can solve the instability in the delay;
b2, changing point constraint or line constraint into volume constraint based on a drilling constraint technology, and searching out a grid interval in which a gravity magnetic anomaly data volume exists by adopting a probability imaging technology to realize dimension reduction of a solution space;
in the step, according to known prior information of well drilling and the like, plane control constraints composed of the plane control points can be established, corresponding physical property parameter constraint arrays are established by the control points, constraint conditions are formed, and constraint modeling and constraint inversion of the physical property inversion model are achieved. The geologic body density distribution has certain continuity and has correlation with the gravity anomaly in spatial distribution, and well constraint and constraint expansion can be performed by combining drilling data with the results of gravity anomaly and regularization normal extension. And on the basis of considering a certain known drilling constraint, judging whether the known drilling constraint is drilled in a geological field source body, if so, combining the abnormal distribution of the residual gravity and the regularized extension result to give a quantitative interval constraint range, expanding the constraint, expanding one point of the well constraint into a whole, taking a grid unit passed by the well as a fixed value, and further iteratively modifying other expanded areas along with the inversion iteration process to realize the space dimension reduction.
B3, performing three-dimensional quantitative inversion of apparent density and apparent magnetization intensity by adopting a multi-method inversion system for the inversion target function of the solution space after dimensionality reduction;
the method in the multi-method inversion system comprises a golden section inversion method based on correlation search, a three-dimensional gravity inversion method based on rapid near-end objective function optimization and a rock physical relationship constraint inversion method based on machine learning and multivariate statistics; the inversion method overcomes the difficulty of inversion multi-solution and improves the reliability of the inversion of apparent density and apparent magnetization.
The golden section method is a classical algorithm in optimization calculation and has the characteristics of simple structure and high global optimization precision, and the basic idea of the method is to gradually narrow the search range based on a symmetry principle and an equal-ratio contraction principle.
The mathematical theory basis of the three-dimensional gravity inversion method based on the rapid near-end objective function optimization is Lipschitz continuity in mathematical analysis.
The rock physical relationship constraint inversion method based on machine learning and multivariate statistics forms an algorithm flow by utilizing the machine learning and the multivariate statistics, and the flow can better reflect the spatial statistics rock physical relationship; based on this, the rock physical relationship constraint inversion method based on machine learning and multivariate statistics in the step R3
The method comprises the following steps:
r1, establishing a cross variation function reflecting the spatial statistical rock physical relationship by taking the logging data and the seismic frame as constraint conditions;
r2, clustering the velocity and density distribution obtained by independent seismic inversion and gravity inversion respectively, and redistributing the velocity value and the density value to the clustering result according to the statistical result of the rock sample data and the logging data to obtain a new velocity and density model;
and R3, taking the cross variation function as a prior function of an inversion target function, and performing gravity inversion of rock physical relationship constraint based on the new speed and density.
Compared with the traditional inversion method, the inversion method has the advantages that the model obtained by constraint inversion can better fit the gravity field and the rock physical relationship, and geological interpretation is facilitated.
The heavy magnetoelectric multiparameter collaborative imaging in step S2 of the present embodiment includes resistivity imaging and multiparameter imaging;
the resistivity imaging method specifically comprises the following steps: performing joint inversion on the three-dimensional controllable electromagnetic data and the seismic data by adopting a cross gradient operator to obtain resistivity imaging; the three-dimensional controllable electromagnetic data and the seismic data contain mutually complementary information, and more reliable underground structure imaging can be obtained by jointly inverting the controllable electromagnetic data and the seismic data.
The multi-parameter imaging method specifically comprises the following steps:
y1, performing numerical simulation of the seismic wave field and the electromagnetic wave field by adopting a staggered grid finite difference method based on a seismic wave equation and a Maxwell equation;
y2, constructing a cross gradient-based electromagnetic and seismic corresponding joint inversion objective function based on a data fitting term, a model regularization term and a cross gradient constraint of electromagnetic data and seismic data obtained by numerical simulation;
y3, solving the constructed joint inversion target function by adopting a linear conjugate gradient algorithm;
in the step, a Conjugate Gradient method (Conjugate Gradient) is a method between the steepest descent method and the newton method, only first derivative information is used, but the defect that the steepest descent method is slow in convergence is overcome, and the defect that the newton method needs to store, calculate and invert a Hesse matrix is avoided. Among various optimization algorithms, the conjugate gradient method is a very important one. Its advantages are less memory needed, high step convergence and stability, and no need of external parameters.
And Y4, performing electromagnetic and seismic joint inversion based on a joint inversion objective function solving process, and realizing multi-parameter imaging.
In the step Y4, when performing the joint inversion of the electromagnetism and the earthquake, the method adopted sequentially comprises: and adjusting regularization factors by adopting a self-adaptive regularization method through data fitting based on two adjacent iterative processes, carrying out self-adaptive correction based on K-means clustering and regression analysis, carrying out cross gradient weighting, carrying out frequency division multi-scale inversion and multi-grid division, and setting physical constraint range inversion.
Specifically, the adaptive regularization method is a method for adaptively adjusting the value of a regularization factor according to a certain criterion in the inversion process; the boundary of the geologic body or the stratum is determined based on a K-value clustering algorithm, and the clustering algorithm is a classic clustering method, has the characteristics of simplicity, rapidness, automation and the like, and is widely applied to the fields of image processing, pattern recognition and the like. And after clustering the resistivity and speed models by using a K-means clustering technology, adjusting the resistivity or speed value by using a regression analysis technology. Therefore, the complete adaptive correction method is as follows:
(1) dividing the resistivity and speed model into different regions by utilizing a K-means clustering technology and extracting an interested abnormal body region;
(2) in the resistivity and speed model, corresponding to an abnormal body region of the same geologic body, obtaining the relationship between the resistivity and the speed in the corresponding abnormal body region by utilizing a regression analysis technology;
(3) for the same geologic body in the resistivity and velocity model, the relationship between the resistivity and the velocity obtained by the logarithmic regression analysis is utilized, one physical property with high reliability is utilized to correct the other physical property, and the relative resolving power of electromagnetism (resistivity) and earthquake (velocity) for a certain geologic body is determined according to a specific practical situation.
When using the frequency division multi-scale method, controllable source electromagnetic data of all frequencies is first divided into three frequency bands from low frequency to high frequency: the method comprises the steps of low-frequency band, medium-frequency band and high-frequency band, then, performing controllable source electromagnetic and earthquake joint inversion step by utilizing electromagnetic data of the three frequency bands, and enabling electromagnetic data of each scale to have step-by-step inclusion relation in the joint inversion.
The embodiment provides a specific method for frequency division multi-scale inversion: the method comprises the steps of firstly carrying out joint inversion by using low-frequency-band electromagnetic data and seismic data to obtain a resistivity and velocity model, then using the resistivity and velocity model obtained by joint inversion of the low-frequency-band electromagnetic data and the seismic data as an initial model for joint inversion of the low-frequency-band electromagnetic data and the medium-frequency-band electromagnetic data and the seismic data, using the low-frequency-band data and the medium-frequency-band data instead of the medium-frequency data as the electromagnetic data in the step, and finally using the resistivity and velocity model obtained by joint inversion of the low-frequency-band electromagnetic data and the medium-frequency-band seismic data as the initial model for joint inversion of the low-frequency-band electromagnetic data, the medium-frequency-band electromagnetic data and the seismic data until a termination condition is met and an iteration process is stopped.
The basic idea of the multiple-grid inversion scheme is to divide a model area to be inverted into a series of grids with different scales (thicknesses), a thicker subdivision grid is used for inverting low-frequency controllable source electromagnetic data, the used subdivision grid is finer along with the addition of higher-frequency electromagnetic data, the key of the multiple-grid inversion scheme is the conversion from the coarse grid to the fine grid, model parameters of the coarse grid are directly transmitted to the fine grid contained in the coarse grid after inversion is carried out in the coarse grid, and the coarse grid is used as an initial model for inversion on the fine grid.
In cross-gradient weighting: in the joint inversion process, different weighting factors are given to corresponding grids according to the distance from the section where the seismic survey line is located, so that the resistivity and the speed within a certain distance range from the section where the seismic survey line is located have strong coupling, and the imaging precision of the resistivity data volume in the area without the seismic survey line within the range is improved to a certain extent.
In step S2 of this embodiment, the method for performing the joint inversion of the heavy seismic interface specifically includes:
and performing underground density interface inversion by combining a random inversion method and a cross gradient inversion method with gravity and seismic data. The joint inversion method improves the accuracy of interface depth inversion, and particularly enhances the deep structure identification capability.
In step S3 of this embodiment, the method for performing the gravity magnetic and electric shock physical parameter joint collaborative imaging specifically includes:
based on the combined collaborative imaging result of the gravity, magnetism, electricity and earthquake physical property parameters and the comprehensive response characteristics of the gravity, magnetism, electricity and earthquake in the combined inversion result of the gravity-earthquake interface, the modeling-inversion of the gravity, magnetism, electricity and earthquake is carried out, the physical property parameters of different geophysical fields are inverted, and the mutual constraint inversion collaborative imaging of the physical property parameters of the different geophysical fields of the gravity, magnetism, electricity and earthquake is utilized.
Claims (8)
1. A three-dimensional gravity magnetic-electric shock multi-parameter collaborative inversion method is characterized by comprising the following steps:
s1, acquiring seismic and electromagnetic data of the area to be detected;
s2, performing magnetotelluric two-three dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter co-imaging and heavy seismic interface joint inversion according to the seismic and electromagnetic data;
s3, performing the combined collaborative imaging of the gravity electromagnetic and seismoelectric physical parameters according to the results of magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion and gravity electromagnetic multi-parameter collaborative imaging;
and S4, constructing an imaging result for geological comprehensive interpretation by combining the gravity magnetic seismographic physical property parameter combined collaborative imaging result and the gravity-seismographic interface combined inversion result, and realizing three-dimensional gravity magnetic seismographic multi-parameter collaborative inversion.
2. The three-dimensional gravity magnetoseis multiparameter collaborative inversion method according to claim 1, wherein in the step S2, performing magnetotelluric two-dimensional inversion comprises the following sub-steps:
a1, performing magnetotelluric two-dimensional inversion and dip inversion of the spherical column coordinates according to the seismic and electromagnetic data;
a2, performing structured self-adaptive finite element two-dimensional magnetotelluric inversion and modeling according to the seismic and electromagnetic data;
a3, taking the result of the two-dimensional and three-dimensional magnetotelluric inversion of the spherical column coordinates and the result of the two-dimensional and three-dimensional magnetotelluric inversion and modeling of the structured self-adaptive finite elements as the result of the two-dimensional and three-dimensional magnetotelluric inversion;
wherein, step A1 includes the following substeps:
a11, constructing an integrated form of Maxwell equation satisfied by the electric field intensity vector and the magnetic field intensity vector based on the acquired electromagnetic data:
∮H·dl＝∫∫J·dS＝∫∫σE·ds
∮E·dl＝∫∫iωμ0·HdS
wherein E and H are respectively an electric field intensity vector and a magnetic field intensity vector, J is a current density, dl is a contour of a closed integral, dS is an area included by the contour, sigma is a medium conductivity, omega is a circular frequency, and mu0Is a vacuum magnetic conductivity;
a12, dividing the spherical coordinate system into a plurality of grid units by a plurality of parallel spherical surfaces at different intervals in the r, theta and phi directions in the spherical coordinate system, and realizing the discretization of the spherical coordinate system according to a staggered grid subdivision mode;
a13, discretizing and expanding the Maxwell equation integral form in a discretized spherical coordinate system, further determining a staggered grid finite difference equation, and taking the staggered grid finite difference equation as the results of magnetotelluric two-dimensional inversion and dip inversion of the spherical coordinate system;
step a2 includes the following substeps:
a21, for earthquake and electromagnetic data, sequentially processing the earthquake and electromagnetic data by adopting a geoelectromagnetic method self-adaptive vector finite element method, a linear equation set solving method based on an auxiliary space precondition algorithm, a parallel method based on a grid partition technology and a geoelectromagnetic three-dimensional inversion method based on self-adaptive finite element forward modeling;
a22, processing the result obtained by the magnetotelluric three-dimensional inversion method based on the forward modeling of the self-adaptive finite element as the result of the magnetotelluric inversion and modeling of two three-dimensional earth of the structural self-adaptive finite element;
in step a21, the magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling specifically includes the following sub-steps:
t1, obtaining data of the region to be inverted according to the grid partitioning technology-based parallel method, and determining the region to be inverted according to the data;
t2, discretizing the region to be inverted, and taking the discretized grid as an inversion grid;
wherein the inversion grid does not change during the inversion process;
t3, setting a forward net for each frequency in each iteration of the inversion process;
wherein the initial value of the forward grid is the same as the initial value of the inversion grid;
t4, in the iteration of the inversion process, self-adaptive grid encryption is carried out on the forward grid for a plurality of times by using the posterior error estimated value, and a partial derivative is calculated on the final forward grid;
and T5, mapping the partial derivative calculated on the final forward grid back to the corresponding inversion grid, and realizing magnetotelluric three-dimensional inversion.
3. The three-dimensional gravity magnetoseis multiparameter cooperative inversion method according to claim 2, wherein in step S2, performing three-dimensional quantitative inversion of apparent density and apparent magnetization comprises the following sub-steps:
b1, acquiring a gravity magnetic anomaly data body of the whole space of the region to be detected by utilizing a regularization downward extension technology;
b2, changing point constraint or line constraint into volume constraint based on a drilling constraint technology, and searching out a grid interval in which a gravity magnetic anomaly data volume exists by adopting a probability imaging technology to realize dimension reduction of a solution space;
b3, performing three-dimensional quantitative inversion of apparent density and apparent magnetization intensity by adopting a multi-method inversion system for the inversion target function of the solution space after dimensionality reduction;
the method in the multi-method inversion system comprises a golden section inversion method based on correlation search, a three-dimensional gravity inversion method based on rapid near-end objective function optimization and a rock physical relationship constraint inversion method based on machine learning and multivariate statistics.
4. The three-dimensional gravity magnetoseis multiparameter cooperative inversion method according to claim 3, wherein the rock physical relationship constraint inversion method based on machine learning and multivariate geostatistics in the step B3 comprises the following sub-steps:
r1, establishing a cross variation function reflecting the spatial statistical rock physical relationship by taking the logging data and the seismic frame as constraint conditions;
r2, clustering the velocity and density distribution obtained by independent seismic inversion and gravity inversion respectively, and redistributing the velocity value and the density value to the clustering result according to the statistical result of the rock sample data and the logging data to obtain a new velocity and density model;
and R3, taking the cross variation function as a prior function of an inversion target function, and performing gravity inversion of rock physical relationship constraint based on the new speed and density.
5. The three-dimensional gravity magnetoseis multiparameter cooperative inversion method according to claim 1, wherein in step S2, performing gravity magnetoseis multiparameter cooperative imaging includes resistivity imaging and multiparameter imaging;
the resistivity imaging method specifically comprises the following steps: performing joint inversion on the three-dimensional controllable electromagnetic data and the seismic data by adopting a cross gradient operator to obtain resistivity imaging;
the multi-parameter imaging method specifically comprises the following steps:
y1, performing numerical simulation of the seismic wave field and the electromagnetic wave field by adopting a staggered grid finite difference method based on a seismic wave equation and a Maxwell equation;
y2, constructing a cross gradient-based electromagnetic and seismic corresponding joint inversion objective function based on a data fitting term, a model regularization term and a cross gradient constraint of electromagnetic data and seismic data obtained by numerical simulation;
y3, solving the constructed joint inversion target function by adopting a linear conjugate gradient algorithm;
and Y4, performing electromagnetic and seismic joint inversion based on a joint inversion objective function solving process, and realizing multi-parameter imaging.
6. The three-dimensional gravity electromagnetic positive multi-parameter collaborative inversion method according to claim 5, wherein in the step Y4, when performing the joint inversion of electromagnetism and earthquake, the adopted methods are as follows in sequence: and adjusting regularization factors by adopting a self-adaptive regularization method through data fitting based on two adjacent iterative processes, carrying out self-adaptive correction based on K-means clustering and regression analysis, carrying out cross gradient weighting, carrying out frequency division multi-scale inversion and multi-grid division, and setting physical constraint range inversion.
7. The three-dimensional gravity magnetic-electric-shock multi-parameter collaborative inversion method according to claim 1, wherein in the step S2, the method for performing the joint inversion of the gravity-shock interface specifically comprises:
and performing underground density interface inversion by combining a random inversion method and a cross gradient inversion method with gravity and seismic data.
8. The three-dimensional gravity magnetic and electric shock multi-parameter collaborative inversion method according to claim 1, wherein in the step S3, the method for performing gravity magnetic and electric shock physical parameter joint collaborative imaging specifically comprises:
based on the comprehensive response characteristics of the gravity, magnetism, electricity and earthquake in the gravity, magnetism, electricity and earthquake physical property parameter joint collaborative imaging result and the gravity and earthquake interface joint inversion result, the gravity, magnetism, electricity and earthquake modeling-inversion is carried out, the physical property parameters of different geophysical fields are inverted, and the physical property parameters of the different geophysical fields of the gravity, magnetism, electricity and earthquake are utilized to carry out mutual constraint inversion collaborative imaging.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111086356.8A CN113917560B (en) | 2021-09-16 | 2021-09-16 | Three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111086356.8A CN113917560B (en) | 2021-09-16 | 2021-09-16 | Three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113917560A true CN113917560A (en) | 2022-01-11 |
| CN113917560B CN113917560B (en) | 2023-04-21 |
Family
ID=79235058
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202111086356.8A Active CN113917560B (en) | 2021-09-16 | 2021-09-16 | Three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113917560B (en) |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115859137A (en) * | 2023-02-20 | 2023-03-28 | 中国地质大学（北京） | Magnetic anomaly vector clustering inversion method based on unsupervised machine learning |
| CN116908928A (en) * | 2023-05-15 | 2023-10-20 | 重庆大学 | Stratum adaptive encryption-based magnetotelluric inversion method |
| CN117148457A (en) * | 2023-08-29 | 2023-12-01 | 长安大学 | Magnetic layer magnetization modulus calculation method and system |
| CN117312766A (en) * | 2023-08-30 | 2023-12-29 | 昆明理工大学 | Electromagnetic detection data joint inversion and NSCT mine mining working face water-rich area identification method |
| CN118011514A (en) * | 2024-04-10 | 2024-05-10 | 成都理工大学 | Prediction method and system applied to basin foundation interface fluctuation and density |
Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110264421A1 (en) * | 2010-04-21 | 2011-10-27 | Charlie Jing | Method For Geophysical Imaging |
| US20130179130A1 (en) * | 2012-01-06 | 2013-07-11 | Technoimaging, Llc | Method of simultaneous imaging of different physical properties using joint inversion of multiple datasets |
| CN104216006A (en) * | 2013-06-04 | 2014-12-17 | 中国石油化工股份有限公司 | Method for increasing imaging quality by using gravity and magnetic electric shock synchronous united parameter model |
| CN104635261A (en) * | 2013-11-13 | 2015-05-20 | 中国石油化工股份有限公司 | Gravity, magnetism, electrical method and earthquake combined split modeling method for mountain fronts |
| CN105005097A (en) * | 2015-07-17 | 2015-10-28 | 中国石油化工股份有限公司 | Method for comprehensive recognition of igneous rocks by employing gravity, magnetism, electromagnetism, and earthquake data |
| CN106443822A (en) * | 2016-08-16 | 2017-02-22 | 中国石油化工股份有限公司 | Geological integrated identification method and device based on gravity-magnetic-electric-seismic three-dimensional joint inversion |
| CN107678072A (en) * | 2017-09-22 | 2018-02-09 | 中国石油化工股份有限公司胜利油田分公司勘探开发研究院西部分院 | Based on magnetic force, earthquake, the united igneous reservoirs Forecasting Methodology of drilling well |
| CN107765337A (en) * | 2016-08-19 | 2018-03-06 | 中国石油化工股份有限公司 | Electrical method and synchronized seismic joint inversion method and system |
-
2021
- 2021-09-16 CN CN202111086356.8A patent/CN113917560B/en active Active
Patent Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110264421A1 (en) * | 2010-04-21 | 2011-10-27 | Charlie Jing | Method For Geophysical Imaging |
| US20130179130A1 (en) * | 2012-01-06 | 2013-07-11 | Technoimaging, Llc | Method of simultaneous imaging of different physical properties using joint inversion of multiple datasets |
| CN104216006A (en) * | 2013-06-04 | 2014-12-17 | 中国石油化工股份有限公司 | Method for increasing imaging quality by using gravity and magnetic electric shock synchronous united parameter model |
| CN104635261A (en) * | 2013-11-13 | 2015-05-20 | 中国石油化工股份有限公司 | Gravity, magnetism, electrical method and earthquake combined split modeling method for mountain fronts |
| CN105005097A (en) * | 2015-07-17 | 2015-10-28 | 中国石油化工股份有限公司 | Method for comprehensive recognition of igneous rocks by employing gravity, magnetism, electromagnetism, and earthquake data |
| CN106443822A (en) * | 2016-08-16 | 2017-02-22 | 中国石油化工股份有限公司 | Geological integrated identification method and device based on gravity-magnetic-electric-seismic three-dimensional joint inversion |
| CN107765337A (en) * | 2016-08-19 | 2018-03-06 | 中国石油化工股份有限公司 | Electrical method and synchronized seismic joint inversion method and system |
| CN107678072A (en) * | 2017-09-22 | 2018-02-09 | 中国石油化工股份有限公司胜利油田分公司勘探开发研究院西部分院 | Based on magnetic force, earthquake, the united igneous reservoirs Forecasting Methodology of drilling well |
Non-Patent Citations (1)
| Title |
|---|
| YU PENG等: ""Determination of physical property parameter structure in Middle Guizhou uplift by joint inversion of magnetotelluric, seismic, gravity and magnetic data"" * |
Cited By (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115859137A (en) * | 2023-02-20 | 2023-03-28 | 中国地质大学（北京） | Magnetic anomaly vector clustering inversion method based on unsupervised machine learning |
| CN115859137B (en) * | 2023-02-20 | 2023-04-28 | 中国地质大学（北京） | Magnetic anomaly vector clustering inversion method based on unsupervised machine learning |
| CN116908928A (en) * | 2023-05-15 | 2023-10-20 | 重庆大学 | Stratum adaptive encryption-based magnetotelluric inversion method |
| CN116908928B (en) * | 2023-05-15 | 2024-03-26 | 重庆大学 | Stratum adaptive encryption-based magnetotelluric inversion method |
| CN117148457A (en) * | 2023-08-29 | 2023-12-01 | 长安大学 | Magnetic layer magnetization modulus calculation method and system |
| CN117312766A (en) * | 2023-08-30 | 2023-12-29 | 昆明理工大学 | Electromagnetic detection data joint inversion and NSCT mine mining working face water-rich area identification method |
| CN118011514A (en) * | 2024-04-10 | 2024-05-10 | 成都理工大学 | Prediction method and system applied to basin foundation interface fluctuation and density |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113917560B (en) | 2023-04-21 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| USRE49507E1 (en) | Faulted geological structures having unconformities | |
| CN113917560B (en) | Three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method | |
| Oldenburg et al. | Three dimensional inversion of multisource time domain electromagnetic data | |
| US11042676B2 (en) | Representing structural uncertainty in a mesh representing a geological environment | |
| US11175434B2 (en) | Geologic stratigraphy via implicit and jump functions | |
| US8612194B2 (en) | Updating a subterranean model using at least electromagnetic data | |
| Wang et al. | On a new method of estimating shear wave velocity from conventional well logs | |
| US11249208B2 (en) | Geologic structural model generation | |
| WO2021168230A1 (en) | Physics-constrained deep learning joint inversion | |
| US20150066460A1 (en) | Stratigraphic function | |
| US20140222403A1 (en) | Geologic model via implicit function | |
| Soleimani et al. | Integrated petrophysical modeling for a strongly heterogeneous and fractured reservoir, Sarvak Formation, SW Iran | |
| Zhang et al. | 2D joint inversion of geophysical data using petrophysical clustering and facies deformation | |
| Bhark et al. | A multiscale workflow for history matching in structured and unstructured grid geometries | |
| Yin et al. | Enhancement of dynamic reservoir interpretation by correlating multiple 4D seismic monitors to well behavior | |
| WO2023130074A1 (en) | Geologic modeling framework | |
| Liu et al. | Performance investigations of auxiliary‐space Maxwell solver preconditioned iterative algorithm for controlled‐source electromagnetic induction problems with electrical anisotropy | |
| CN115880455A (en) | Three-dimensional intelligent interpolation method based on deep learning | |
| Zou et al. | Log-constrained inversion based on a conjugate gradient-particle swarm optimization algorithm | |
| Marchant et al. | 3D inversion of electromagnetic logging-while-drilling data | |
| Gao et al. | Research on reservoir inversion method of Paleogene in few well area | |
| Shen et al. | Predication of super thin reservoir based on geostatistical inversion | |
| Singh et al. | Hybrid Strategy for Porosity Distribution Mapping on a Heterogeneous Reservoir Using Artificial Intelligence Method | |
| Hu et al. | Analyzing Fractures Using Time-Lapse Electric Potential Data | |
| WO2024064657A1 (en) | Geologic modeling framework |
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 |