CN113917560B - Three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method - Google Patents
Three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method Download PDFInfo
- Publication number
- CN113917560B CN113917560B CN202111086356.8A CN202111086356A CN113917560B CN 113917560 B CN113917560 B CN 113917560B CN 202111086356 A CN202111086356 A CN 202111086356A CN 113917560 B CN113917560 B CN 113917560B
- Authority
- CN
- China
- Prior art keywords
- inversion
- dimensional
- heavy
- 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.)
- Active
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 heavy electromagnetic vibration multi-parameter collaborative inversion method, which comprises the following steps: s1, acquiring earthquake and electromagnetic data of a region to be detected; s2, performing magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter collaborative imaging and heavy seismology interface joint inversion according to the earthquake and electromagnetic data respectively; s3, carrying out combined collaborative imaging on the physical parameters of the heavy electromagnetic vibration according to the results of the two-dimensional inversion of the magnetotelluric, the three-dimensional quantitative inversion of the apparent density and the apparent magnetization and the collaborative imaging of the heavy electromagnetic vibration multiple parameters; s4, constructing an imaging result for geological comprehensive interpretation by combining the combined collaborative imaging result of the physical parameters of the heavy magnetic and electric shock and the combined inversion result of the heavy shock interface, and realizing three-dimensional heavy magnetic and electric shock multi-parameter collaborative inversion.
Description
Technical Field
The invention belongs to the technical field of geologic structure imaging, and particularly relates to a three-dimensional heavy electromagnetic vibration multi-parameter collaborative inversion method.
Background
Seismic exploration has been the main means of oil and gas exploration and development, mainly due to: (1) the imaging of the underground complex geological structure can be obtained by utilizing the seismic imaging technology, and particularly, the imaging of the underground structure with higher resolution can be provided by Reverse Time Migration (RTM) and Full Waveform Inversion (FWI); (2) the physical property parameter information of the reservoir, such as porosity, fluid saturation, permeability and the like, can be extracted by utilizing the seismic reservoir inversion technology, so that the reservoir fluid type can be further predicted, and the oil-gas property of the reservoir can be evaluated. However, under complex geological conditions, such as igneous rock regions, carbonate rock development regions, reverse-masked fracture zones, etc., severe seismic wave scattering effects and shielding result in weaker illumination energy of the underlying formation, and only seismic data is used to infer that the subsurface structure is multi-solvable. In terms of reservoir description and monitoring, reservoirs of different gas saturation, such as formation water and oil, are used with some ambiguity in reservoir inversion using seismic data because the change in the oil (gas) saturation has no significant effect on the seismic wave velocity and therefore reservoir fluid type and fluid saturation cannot be well predicted using seismic data alone.
The Magnetotelluric (MT) method is particularly important for solving hydrocarbon exploration problems under complex geological conditions (e.g., igneous subsurface formation imaging, sub-salt imaging, etc.), mainly due to: (1) the electromagnetic signal decays weakly as it passes through the high resistance (e.g., igneous rock), (2) the electromagnetic signal is relatively sensitive to low resistance sedimentary formations. Thus, the acquired electromagnetic signals contain resistivity information of the underlying sedimentary formations. However, since MT methods are based primarily on low frequency induced electromagnetic fields, relatively thin hydrocarbon reservoirs cannot be imaged well. The basic equations governing the propagation of seismic waves and electromagnetic waves are wave equation and diffusion equation, respectively, and although both waves propagate with attenuation and scattering properties, the degree of attenuation and scattering is different, i.e. the attenuation and scattering of seismic waves is relatively weak.
Disclosure of Invention
Aiming at the defects in the prior art, the three-dimensional heavy electromagnetic vibration multi-parameter collaborative inversion method provided by the invention solves the problem of low limitation of single-method multi-resolution and resolution in the existing inversion method.
In order to achieve the aim of the invention, the invention adopts the following technical scheme: a three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method comprises the following steps:
s1, acquiring earthquake and electromagnetic data of a region to be detected;
s2, performing magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter collaborative imaging and heavy seismology interface joint inversion according to the earthquake and electromagnetic data respectively;
s3, carrying out combined collaborative imaging on the physical parameters of the heavy electromagnetic vibration according to the results of the two-dimensional inversion of the magnetotelluric, the three-dimensional quantitative inversion of the apparent density and the apparent magnetization and the collaborative imaging of the heavy electromagnetic vibration multiple parameters;
s4, constructing an imaging result for geological comprehensive interpretation by combining the combined collaborative imaging result of the physical parameters of the heavy magnetic and electric shock and the combined inversion result of the heavy shock interface, and realizing three-dimensional heavy magnetic and electric shock multi-parameter collaborative inversion.
Further, in the step S2, performing magnetotelluric two-dimensional inversion includes the following sub-steps:
a1, performing geomagnetic two-dimensional inversion and inclination inversion of spherical column coordinates according to earthquake and electromagnetic data;
a2, carrying out structural self-adaptive finite element two-dimensional three-dimensional magnetotelluric inversion and modeling according to the earthquake and electromagnetic data;
a3, taking a result of the magnetotelluric two-dimensional inversion of the spherical column coordinates and a result of the magnetotelluric two-dimensional inversion and modeling of the structured self-adaptive finite element as a result of the magnetotelluric two-dimensional inversion;
wherein, step A1 comprises the following sub-steps:
a11, constructing a Maxwell equation integral form which is satisfied by an electric field intensity vector and a 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 electric field intensity vector and magnetic field intensity vector, J is current density, dl is enclosing line of closed integration, dS is area contained in the enclosing line, sigma is medium conductivity, omega is circular frequency, mu 0 Is vacuum magnetic permeability;
a12, in the spherical coordinate system, dividing the spherical coordinate system into a plurality of grid units by a plurality of parallel spherical surfaces at different intervals along the r, theta and phi directions, so as to realize discretization of the spherical coordinate system according to an interlaced grid subdivision form;
a13, performing discretization expansion on the integral form of the Maxwell equation under a discretized spherical coordinate system, and further determining a staggered grid finite difference equation, wherein the staggered grid finite difference equation is used as the result of magnetotelluric two-dimensional inversion and the titler inversion of the spherical coordinate system;
step A2 comprises the following sub-steps:
a21, for seismic and electromagnetic data, sequentially adopting a magnetotelluric 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 magnetotelluric three-dimensional inversion method based on a self-adaptive finite element forward algorithm to process the seismic and electromagnetic data;
a22, processing a result obtained by a magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling, and taking the result as a structural adaptive finite element two-dimensional magnetotelluric inversion and modeling result;
in the step a21, the magnetotelluric three-dimensional inversion method based on the adaptive finite element forward modeling specifically includes the following sub-steps:
t1, obtaining data of a region to be inverted according to a parallel method based on a grid partition technology, and determining the region to be inverted according to the data;
t2, discretizing a 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 grid for each frequency in each iteration of the inversion process;
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, performing adaptive grid encryption on the forward grid for a plurality of times by using the posterior error estimation value, and calculating a partial derivative on the final forward grid;
and T5, mapping the partial derivative calculated on the final forward grid back to the corresponding inversion grid, so as to realize magnetotelluric three-dimensional inversion.
Further, in the step S2, performing three-dimensional quantitative inversion of apparent density and apparent magnetization includes the following sub-steps:
b1, acquiring a heavy magnetic abnormal data body of the whole space of the region to be detected by using a regularized downward delay technology;
b2, changing point constraint or line constraint into volume constraint based on a well drilling constraint technology, and searching a grid interval in which a heavy magnetic abnormal data body exists by adopting a probability imaging technology to realize dimension reduction of a solution space;
b3, performing three-dimensional quantitative inversion on the apparent density and apparent magnetization intensity by adopting a multi-method inversion system for the inversion objective function of the solution space after dimension 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 petrophysical relation constraint inversion method based on machine learning and multivariate statistics.
Further, the petrophysical relationship constraint inversion method based on machine learning and multivariate statistics in the step B3 comprises the following sub-steps:
r1, taking logging data and an earthquake frame as constraint conditions, and establishing a cross variation function reflecting a space statistics rock physical relationship;
r2, clustering the speed and density distribution obtained by independent seismic inversion and gravity inversion respectively, and reassigning the speed value and the density value to a clustering result according to the statistical result of the rock sample data and the logging data to obtain a new speed and density model;
and R3, taking the cross variation function as a priori function of an inversion objective function, and carrying out gravity inversion of the petrophysical relationship constraint based on the new speed and the new density.
Further, in the step S2, the performing of the magneto-electric multi-parameter collaborative imaging includes resistivity imaging and multi-parameter imaging;
the resistivity imaging method specifically comprises the following steps: performing joint inversion on three-dimensional controllable electromagnetic data and seismic data by adopting an operator based on a cross gradient to obtain resistivity imaging;
the multi-parameter imaging method specifically comprises the following steps:
y1, carrying out numerical simulation on a seismic wave field and an electromagnetic wave field by adopting a staggered grid finite difference method based on a seismic wave equation and a Maxwell equation;
y2, constructing a joint inversion objective function corresponding to electromagnetism and earthquake based on a cross gradient based on a data fitting term, a model regularization term and a cross gradient constraint of the electromagnetic data and the earthquake data obtained through numerical simulation;
y3, solving the constructed joint inversion objective function by adopting a linear conjugate gradient algorithm;
and Y4, carrying out joint inversion of electromagnetism and earthquake based on a joint inversion objective function solving process, and realizing multi-parameter imaging.
Further, in the step Y4, when performing electromagnetic and seismic joint inversion, the method adopted sequentially includes: the data fitting based on the adjacent two iterative processes is adopted, the regularization factor is adjusted by adopting an adaptive regularization method, the adaptive correction is carried out based on K-means clustering and regression analysis, the cross gradient weighting is carried out, the frequency division multi-scale inversion and the multiple grid division are carried out, and the inversion of the physical constraint range is set.
Further, in the step S2, the method for performing the joint inversion of the heavy-vibration interface specifically includes:
and carrying out underground density interface inversion by adopting a random inversion method and a cross gradient inversion method and combining gravity and seismic data.
Further, in the step S3, the method for performing the combined collaborative imaging of the physical parameters of the heavy electromagnetic vibration specifically includes:
and carrying out modeling-inversion of the heavy, magnetic, electric and earthquake based on the comprehensive response characteristics of the heavy, magnetic, electric and earthquake in the heavy electromagnetic and earthquake physical parameter joint collaborative imaging result and the heavy earthquake interface joint inversion result, inverting different geophysical field physical parameters, and carrying out inversion collaborative imaging by utilizing the mutual constraint of the heavy, magnetic, electric and earthquake different geophysical field physical parameters.
The beneficial effects of the invention are as follows:
(1) The invention innovatively constructs the heavy earthquake, electric earthquake and heavy electromagnetic earthquake collaborative inversion and imaging technology based on multi-parameter modeling and constraint of density, speed, resistivity, magnetic susceptibility and the like;
(2) The method is applied to oil gas resource and mineral resource exploration, achieves outstanding application effects, introduces new algorithm and simulation technology based on the internal relation of petrophysical, forms a heavy-vibration electric interface joint inversion technology and a heavy-magnetic electric vibration 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 such as density, magnetization intensity, resistivity and the like is developed, the span from a single method to a multi-method multi-dimensional collaborative interpretation technology is realized, and the method is an innovative technology for combined extraction and enhancement of heavy magnetic and electric weak information of a geological target body.
Drawings
FIG. 1 is a flow chart of a three-dimensional heavy electromagnetic vibration 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 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 all the inventions which make use of the inventive concept are protected by the spirit and scope of the present invention as defined and defined in the appended claims to those skilled in the art.
As shown in FIG. 1, the three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method comprises the following steps:
s1, acquiring earthquake and electromagnetic data of a region to be detected;
s2, performing magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter collaborative imaging and heavy seismology interface joint inversion according to the earthquake and electromagnetic data respectively;
s3, carrying out combined collaborative imaging on the physical parameters of the heavy electromagnetic vibration according to the results of the two-dimensional inversion of the magnetotelluric, the three-dimensional quantitative inversion of the apparent density and the apparent magnetization and the collaborative imaging of the heavy electromagnetic vibration multiple parameters;
s4, constructing an imaging result for geological comprehensive interpretation by combining the combined collaborative imaging result of the physical parameters of the heavy magnetic and electric shock and the combined inversion result of the heavy shock interface, and realizing three-dimensional heavy magnetic and electric shock multi-parameter collaborative inversion.
In step S2 of the present embodiment, performing magnetotelluric two-dimensional inversion includes the following sub-steps:
a1, performing geomagnetic two-dimensional inversion and inclination inversion of spherical column coordinates according to earthquake and electromagnetic data;
a2, carrying out structural self-adaptive finite element two-dimensional three-dimensional magnetotelluric inversion and modeling according to the earthquake and electromagnetic data;
a3, taking a result of the magnetotelluric two-dimensional inversion of the spherical column coordinates and a result of the magnetotelluric two-dimensional inversion and modeling of the structured self-adaptive finite element as a result of the magnetotelluric two-dimensional inversion;
wherein, step A1 comprises the following sub-steps:
a11, constructing a Maxwell equation integral form which is satisfied by an electric field intensity vector and a 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 electric field intensity vector and magnetic field intensity vector, J is current density, dl is enclosing line of closed integration, dS is area contained in the enclosing line, sigma is medium conductivity, omega is circular frequency, mu 0 Is vacuum magnetic permeability;
a12, in the spherical coordinate system, dividing the spherical coordinate system into a plurality of grid units by a plurality of parallel spherical surfaces at different intervals along the r, theta and phi directions, so as to realize discretization of the spherical coordinate system according to an interlaced grid subdivision form;
a13, performing discretization expansion on the integral form of the Maxwell equation under a discretized spherical coordinate system, and further determining a staggered grid finite difference equation, wherein the staggered grid finite difference equation is used as the result of magnetotelluric two-dimensional inversion and the titler inversion of the spherical coordinate system;
in particular, in the frequency range of magnetotelluric studies (about 10 5 ～10 -5 Hz), because the frequency is relatively reduced, the quasi-stable electromagnetic field is satisfied, the displacement current can be ignored, the magnetic permeability of the underground medium is taken as the vacuum magnetic permeability, and a corresponding Maxwell equation integral form is constructed; meanwhile, as the discretization of the difference form of the Max equation set under the spherical coordinate system is complex, the integral form of the Max equation set is discretized to obtain a linear equation, and therefore mesh subdivision is needed before discretization.
The step A2 comprises the following sub-steps:
a21, for seismic and electromagnetic data, sequentially adopting a magnetotelluric 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 magnetotelluric three-dimensional inversion method based on a self-adaptive finite element forward algorithm to process the seismic and electromagnetic data;
a22, processing a result obtained by a magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling, and taking the result as a structural adaptive finite element two-dimensional magnetotelluric inversion and modeling result;
specifically, in step a21, in the magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling, a three-dimensional forward modeling model is adopted, the model built by forward modeling carries out adaptive encryption on grids through a magnetotelluric method adaptive vector finite element method, the adaptive grid encryption can automatically encrypt a region with larger error according to posterior error estimation, and compared with global grid encryption, the error reduction speed is obviously faster than that of global grid encryption; solving an inverted objective function by a linear equation set solving method of an auxiliary space precondition algorithm; the parallel approach of the grid partition technique distributes computing tasks to multiple processes according to grid partitions, with the individual processes communicating using MPI to maintain synchronization during the computing process. The inversion is realized by using a forward algorithm, and the forward and inversion grid separation is performed, so that the inversion precision problem is solved.
Based on this, in the above step a21, the magnetotelluric three-dimensional inversion method based on the adaptive finite element forward modeling specifically includes the following sub-steps:
t1, obtaining data of a region to be inverted according to a parallel method based on a grid partition technology, and determining the region to be inverted according to the data;
t2, discretizing a 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 grid for each frequency in each iteration of the inversion process;
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, performing adaptive grid encryption on the forward grid for a plurality of times by using the posterior error estimation value, and calculating a partial derivative on the final forward grid;
and T5, mapping the partial derivative calculated on the final forward grid back to the corresponding inversion grid, so as to realize magnetotelluric three-dimensional inversion.
In step S2 of the present embodiment, performing three-dimensional quantitative inversion of apparent density and apparent magnetization includes the following sub-steps:
b1, acquiring a heavy magnetic abnormal data body of the whole space of the region to be detected by using a regularized downward delay technology;
in this step, bit field down-extension is a typical uncertainty problem, and due to errors in the actual observations and some disturbance information, high frequency oscillation effects occur when down-extending to the vicinity of the field source, often covering the useful information. The regularization method is one of the most stable algorithms, and can solve instability in downlinks;
b2, changing point constraint or line constraint into volume constraint based on a well drilling constraint technology, and searching a grid interval in which a heavy magnetic abnormal data body exists by adopting a probability imaging technology to realize dimension reduction of a solution space;
in the step, according to the known prior information such as drilling and the like, plane control constraints consisting of the plane control points can be established, corresponding physical parameter constraint arrays are established by the control points, constraint conditions are formed, and constraint modeling and constraint inversion of a physical inversion model are realized. The density distribution of the geologic body has certain continuity and has correlation with the gravity anomaly in the space distribution, and well constraint and constraint expansion can be carried out by utilizing drilling data and combining the gravity anomaly and regularized forward-extending result. Considering that whether the known drilling constraint meets a geological field source body is judged on the basis of knowing a certain drilling constraint, if the known drilling constraint meets the geological field source body, a quantitative interval constraint range is given by combining the abnormal distribution of residual gravity and regularized delay results, the constraint is expanded, one point of the well constraint is expanded into a body, grid cells penetrated by a well are fixed values, and other expanded areas are further modified in an iterative manner along with an inversion iterative process, so that the dimensional reduction of a solution space is realized.
B3, performing three-dimensional quantitative inversion on the apparent density and apparent magnetization intensity by adopting a multi-method inversion system for the inversion objective function of the solution space after dimension reduction;
the method in the multi-method inversion system comprises a golden section inversion method based on related search, a three-dimensional gravity inversion method based on rapid near-end objective function optimization and a petrophysical relation constraint inversion method based on machine learning and multivariate statistics; the inversion method overcomes the difficulty of inversion of multiple solutions and improves the reliability of inversion of apparent density and apparent magnetization intensity.
The golden section method is a classical algorithm in optimization calculation and has the characteristics of simple structure and high global optimizing precision, and the basic idea of the method is to gradually narrow the searching range based on a symmetry principle and an isopiestic 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 relation constraint inversion method based on machine learning and multivariate statistics utilizes the machine learning and multivariate statistics to form an algorithm flow, and the flow can better reflect the space statistical rock physical relation; based on this, the petrophysical relation constraint inversion method based on machine learning and multivariate statistics in the above step R3
Comprises the following sub-steps:
r1, taking logging data and an earthquake frame as constraint conditions, and establishing a cross variation function reflecting a space statistics rock physical relationship;
r2, clustering the speed and density distribution obtained by independent seismic inversion and gravity inversion respectively, and reassigning the speed value and the density value to a clustering result according to the statistical result of the rock sample data and the logging data to obtain a new speed and density model;
and R3, taking the cross variation function as a priori function of an inversion objective function, and carrying out gravity inversion of the petrophysical relationship constraint based on the new speed and the new density.
Compared with the traditional inversion method, the inversion method has the advantages that the model obtained by constraint inversion can be better fit with the gravity field and the petrophysical relationship, and geological interpretation is facilitated.
The heavy magnetism and electricity multi-parameter collaborative imaging in step S2 of the present embodiment includes resistivity imaging and multi-parameter imaging;
the resistivity imaging method specifically comprises the following steps: performing joint inversion on three-dimensional controllable electromagnetic data and seismic data by adopting an operator based on a cross gradient to obtain resistivity imaging; the three-dimensional controllable electromagnetic data and the seismic data contain complementary information, and the joint inversion of the controllable electromagnetic data and the seismic data can obtain more reliable underground structure imaging.
The multi-parameter imaging method specifically comprises the following steps:
y1, carrying out numerical simulation on a seismic wave field and an electromagnetic wave field by adopting a staggered grid finite difference method based on a seismic wave equation and a Maxwell equation;
y2, constructing a joint inversion objective function corresponding to electromagnetism and earthquake based on a cross gradient based on a data fitting term, a model regularization term and a cross gradient constraint of the electromagnetic data and the earthquake data obtained through numerical simulation;
y3, solving the constructed joint inversion objective function by adopting a linear conjugate gradient algorithm;
in this step, the conjugate gradient method (Conjugate Gradient) is one method between the steepest descent method and the newton method, which only needs to use the first derivative information, but overcomes the defect of slow convergence of the steepest descent method, and avoids the defect that the newton method needs to store and calculate the Hesse matrix and perform inversion, and the conjugate gradient method is one of the most useful methods for solving large linear equation sets and one of the most effective algorithms for solving large nonlinear optimization. Among various optimization algorithms, the conjugate gradient method is a very important one. Its advantages are small memory, high step convergence and stability, and no need of external parameters.
And Y4, carrying out joint inversion of electromagnetism and earthquake based on a joint inversion objective function solving process, and realizing multi-parameter imaging.
In the step Y4, when the joint inversion of the electromagnetic and the earthquake is performed, the method adopted sequentially comprises the following steps: the data fitting based on the adjacent two iterative processes is adopted, the regularization factor is adjusted by adopting an adaptive regularization method, the adaptive correction is carried out based on K-means clustering and regression analysis, the cross gradient weighting is carried out, the frequency division multi-scale inversion and the multiple grid division are carried out, and the inversion of the physical constraint range is set.
Specifically, the self-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 stratum is determined based on a K value clustering algorithm, and the clustering algorithm is a classical 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. After the resistivity and the speed model are clustered by using a K-means clustering technology, a regression analysis technology is adopted to adjust the resistivity or the speed value. Therefore, the complete adaptive correction method is as follows:
(1) Dividing the resistivity and speed model into different areas by using a K-means clustering technology and extracting an interested abnormal body area;
(2) In the resistivity and speed model, an abnormal body region corresponding to the same geologic body is obtained by using a regression analysis technology, and the relationship between the resistivity and the speed in the abnormal body region is obtained correspondingly;
(3) For the same geologic body in the resistivity and velocity model, the relationship between the resistivity and the velocity obtained by the logistic 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) to a certain geologic body is determined according to specific practical situations.
When the frequency division multi-scale method is utilized, controllable source electromagnetic data of all frequencies are firstly divided into three frequency bands from low frequency to high frequency: and then carrying out controlled source electromagnetic and seismic joint inversion step by utilizing electromagnetic data of the three frequency bands, wherein the electromagnetic data of each scale has a step-by-step inclusion relation in the joint inversion.
The specific implementation of the frequency division multi-scale inversion method is given in the embodiment: the method comprises the steps of performing joint inversion by using low-frequency electromagnetic data and seismic data to obtain a resistivity and velocity model, then performing joint inversion on the low-frequency electromagnetic data and the seismic data to obtain a resistivity and velocity model as initial models of joint inversion on the low-frequency electromagnetic data, the medium-frequency electromagnetic data and the seismic data, wherein the electromagnetic data are low-frequency electromagnetic data, medium-frequency electromagnetic data instead of medium-frequency electromagnetic data, and finally performing joint inversion on the low-frequency electromagnetic data, the medium-frequency electromagnetic data and the seismic data to obtain a resistivity and velocity model as initial models of joint inversion on the low-frequency electromagnetic data, the medium-frequency electromagnetic data, the high-frequency electromagnetic data (namely full-frequency electromagnetic data) and the seismic data until a termination condition is met, and the iterative process is stopped.
The basic idea of the multi-grid inversion scheme is that a model area to be inverted is divided into a series of grids with different scales (thicknesses), coarse divided grids are used for inverting low-frequency controllable source electromagnetic data, the used divided grids are finer and finer along with the addition of higher-frequency electromagnetic data, the key of the multi-grid inversion scheme is that the conversion from the coarse grids to the fine grids is carried out, after inversion is carried out in the coarse grids, model parameters of the coarse grids are directly transferred to the fine grids contained in the coarse grids to serve as an initial model for inversion on the fine grids, and the method has the characteristics of small calculated amount, easiness in implementation and the like.
When cross gradient weighting is performed: in the joint inversion process, different weight factors are given to corresponding grids according to the distance from the section of the seismic line, so that the resistivity and the speed within a certain distance range from the section of the seismic line have stronger coupling, and the imaging precision of the resistivity data volume of the area without the seismic line within the range is improved to a certain extent.
In step S2 of the present embodiment, the method for performing the joint inversion of the heavy-seismic interface specifically includes:
and carrying out underground density interface inversion by adopting a random inversion method and a cross gradient inversion method and combining gravity and seismic data. The joint inversion method improves the accuracy of interface depth inversion, and particularly enhances the identification capability of deep layer structures.
In step S3 of the present embodiment, the method for performing the combined collaborative imaging of the physical parameters of the heavy electromagnetic vibration specifically includes:
and carrying out modeling-inversion of the heavy, magnetic, electric and earthquake based on the comprehensive response characteristics of the heavy, magnetic, electric and earthquake in the heavy electromagnetic and earthquake physical parameter joint collaborative imaging result and the heavy earthquake interface joint inversion result, inverting different geophysical field physical parameters, and carrying out inversion collaborative imaging by utilizing the mutual constraint of the heavy, magnetic, electric and earthquake different geophysical field physical parameters.
Claims (4)
1. The three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method is characterized by comprising the following steps of:
s1, acquiring earthquake and electromagnetic data of a region to be detected;
s2, performing magnetotelluric two-dimensional inversion, apparent density and apparent magnetization three-dimensional quantitative inversion, heavy magnetoelectric multi-parameter collaborative imaging and heavy seismology interface joint inversion according to the earthquake and electromagnetic data respectively;
s3, carrying out combined collaborative imaging on the physical parameters of the heavy electromagnetic vibration according to the results of the two-dimensional inversion of the magnetotelluric, the three-dimensional quantitative inversion of the apparent density and the apparent magnetization and the collaborative imaging of the heavy electromagnetic vibration multiple parameters;
s4, combining physical parameters of the heavy electromagnetic vibration, a collaborative imaging result and a heavy vibration interface joint inversion result, constructing an imaging result for geological comprehensive interpretation, and realizing three-dimensional heavy electromagnetic vibration multi-parameter collaborative inversion;
in the step S2, performing magnetotelluric two-dimensional inversion comprises the following sub-steps:
a1, performing geomagnetic two-dimensional inversion and inclination inversion of spherical column coordinates according to earthquake and electromagnetic data;
a2, carrying out structural self-adaptive finite element two-dimensional three-dimensional magnetotelluric inversion and modeling according to the earthquake and electromagnetic data;
a3, taking a result of the magnetotelluric two-dimensional inversion of the spherical column coordinates and a result of the magnetotelluric two-dimensional inversion and modeling of the structured self-adaptive finite element as a result of the magnetotelluric two-dimensional inversion;
wherein, step A1 comprises the following sub-steps:
a11, constructing a Maxwell equation integral form which is satisfied by an electric field intensity vector and a 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 electric field intensity vector and magnetic field intensity vector, J is current density, dl is enclosing line of closed integration, dS is area contained in the enclosing line, sigma is medium conductivity, omega is circular frequency, mu 0 Is vacuum magnetic permeability;
a12, in the spherical coordinate system, dividing the spherical coordinate system into a plurality of grid units by a plurality of parallel spherical surfaces at different intervals along the r, theta and phi directions, so as to realize discretization of the spherical coordinate system according to an interlaced grid subdivision form;
a13, performing discretization expansion on the integral form of the Maxwell equation under a discretized spherical coordinate system, and further determining a staggered grid finite difference equation, wherein the staggered grid finite difference equation is used as the result of magnetotelluric two-dimensional inversion and the titler inversion of the spherical coordinate system;
step A2 comprises the following sub-steps:
a21, for seismic and electromagnetic data, sequentially adopting a magnetotelluric 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 magnetotelluric three-dimensional inversion method based on a self-adaptive finite element forward algorithm to process the seismic and electromagnetic data;
a22, processing a result obtained by a magnetotelluric three-dimensional inversion method based on adaptive finite element forward modeling, and taking the result as a structural adaptive finite element two-dimensional magnetotelluric inversion and modeling result;
in the step a21, the magnetotelluric three-dimensional inversion method based on the adaptive finite element forward modeling specifically includes the following sub-steps:
t1, obtaining data of a region to be inverted according to a parallel method based on a grid partition technology, and determining the region to be inverted according to the data;
t2, discretizing a 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 grid for each frequency in each iteration of the inversion process;
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, performing adaptive grid encryption on the forward grid for a plurality of times by using the posterior error estimation value, and calculating a partial derivative on the final forward grid;
t5, mapping the partial derivative calculated on the final forward grid back to the corresponding inversion grid to realize magnetotelluric three-dimensional inversion;
in the step S2, the three-dimensional quantitative inversion of the apparent density and the apparent magnetization intensity includes the following sub-steps:
b1, acquiring a heavy magnetic abnormal data body of the whole space of the region to be detected by using a regularized downward delay technology;
b2, changing point constraint or line constraint into volume constraint based on a well drilling constraint technology, and searching a grid interval in which a heavy magnetic abnormal data body exists by adopting a probability imaging technology to realize dimension reduction of a solution space;
b3, performing three-dimensional quantitative inversion on the apparent density and apparent magnetization intensity by adopting a multi-method inversion system for the inversion objective function of the solution space after dimension reduction;
the method in the multi-method inversion system comprises a golden section inversion method based on related search, a three-dimensional gravity inversion method based on rapid near-end objective function optimization and a petrophysical relation constraint inversion method based on machine learning and multivariate statistics;
in the step S2, performing the multi-parameter cooperative imaging of the heavy magnetism and the electricity comprises resistivity imaging and multi-parameter imaging;
the resistivity imaging method specifically comprises the following steps: performing joint inversion on three-dimensional controllable electromagnetic data and seismic data by adopting an operator based on a cross gradient to obtain resistivity imaging;
the multi-parameter imaging method specifically comprises the following steps:
y1, carrying out numerical simulation on a seismic wave field and an electromagnetic wave field by adopting a staggered grid finite difference method based on a seismic wave equation and a Maxwell equation;
y2, constructing a joint inversion objective function corresponding to electromagnetism and earthquake based on a cross gradient based on a data fitting term, a model regularization term and a cross gradient constraint of the electromagnetic data and the earthquake data obtained through numerical simulation;
y3, solving the constructed joint inversion objective function by adopting a linear conjugate gradient algorithm;
y4, carrying out joint inversion of electromagnetism and earthquake based on a joint inversion objective function solving process, and realizing multi-parameter imaging; in the step S2, the method for performing the joint inversion of the heavy-seismic interface specifically includes:
and carrying out underground density interface inversion by adopting a random inversion method and a cross gradient inversion method and combining gravity and seismic data.
2. The three-dimensional heavy electromagnetic vibration multi-parameter collaborative inversion method according to claim 1, wherein the petrophysical relationship constraint inversion method based on machine learning and multivariate statistics in the step B3 comprises the following sub-steps:
r1, taking logging data and an earthquake frame as constraint conditions, and establishing a cross variation function reflecting a space statistics rock physical relationship;
r2, clustering the speed and density distribution obtained by independent seismic inversion and gravity inversion respectively, and reassigning the speed value and the density value to a clustering result according to the statistical result of the rock sample data and the logging data to obtain a new speed and density model;
and R3, taking the cross variation function as a priori function of an inversion objective function, and carrying out gravity inversion of the petrophysical relationship constraint based on the new speed and the new density.
3. The three-dimensional heavy electromagnetic vibration multi-parameter collaborative inversion method according to claim 1, wherein in the step Y4, when performing joint inversion of electromagnetic and earthquake, the method adopted sequentially comprises: the data fitting based on the adjacent two iterative processes is adopted, the regularization factor is adjusted by adopting an adaptive regularization method, the adaptive correction is carried out based on K-means clustering and regression analysis, the cross gradient weighting is carried out, the frequency division multi-scale inversion and the multiple grid division are carried out, and the inversion of the physical constraint range is set.
4. The three-dimensional heavy electromagnetic vibration multi-parameter collaborative inversion method according to claim 1, wherein in the step S3, the method for performing the heavy electromagnetic vibration physical parameter joint collaborative imaging specifically comprises:
and carrying out modeling-inversion of the heavy, magnetic, electric and earthquake based on the comprehensive response characteristics of the heavy, magnetic, electric and earthquake in the heavy electromagnetic and earthquake physical parameter joint collaborative imaging result and the heavy earthquake interface joint inversion result, inverting different geophysical field physical parameters, and carrying out inversion collaborative imaging by utilizing the mutual constraint of the heavy, magnetic, electric and earthquake different geophysical field physical parameters.
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 CN113917560A (en) | 2022-01-11 |
| CN113917560B true 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) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115859137B (en) * | 2023-02-20 | 2023-04-28 | 中国地质大学（北京） | Magnetic anomaly vector clustering inversion method based on unsupervised machine learning |
| CN116908928B (en) * | 2023-05-15 | 2024-03-26 | 重庆大学 | Stratum adaptive encryption-based magnetotelluric inversion method |
Family Cites Families (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US8498845B2 (en) * | 2010-04-21 | 2013-07-30 | Exxonmobil Upstream Research Company | Method for geophysical imaging |
| US10242126B2 (en) * | 2012-01-06 | 2019-03-26 | Technoimaging, Llc | Method of simultaneous imaging of different physical properties using joint inversion of multiple datasets |
| CN104216006B (en) * | 2013-06-04 | 2017-05-17 | 中国石油化工股份有限公司 | Method for increasing imaging quality by using gravity and magnetic electric shock synchronous united parameter model |
| CN104635261B (en) * | 2013-11-13 | 2018-04-13 | 中国石油化工股份有限公司 | Mountain front weight magnetoelectricity shake joint split modeling method |
| CN105005097B (en) * | 2015-07-17 | 2017-07-07 | 中国石油化工股份有限公司 | Igneous rock method is comprehensively recognized using gravity, magnetic force, electromagnetism, seismic 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 |
| CN107765337B (en) * | 2016-08-19 | 2020-02-21 | 中国石油化工股份有限公司 | Electric method and earthquake synchronous joint inversion method and system |
| CN107678072B (en) * | 2017-09-22 | 2019-08-20 | 中国石油化工股份有限公司胜利油田分公司勘探开发研究院西部分院 | Based on magnetic force, earthquake, the united igneous reservoirs prediction technique of drilling well |
-
2021
- 2021-09-16 CN CN202111086356.8A patent/CN113917560B/en active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN113917560A (en) | 2022-01-11 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| USRE49507E1 (en) | Faulted geological structures having unconformities | |
| Oldenburg et al. | Three dimensional inversion of multisource time domain electromagnetic data | |
| Gueting et al. | Imaging and characterization of facies heterogeneity in an alluvial aquifer using GPR full-waveform inversion and cone penetration tests | |
| US11042676B2 (en) | Representing structural uncertainty in a mesh representing a geological environment | |
| US20140222403A1 (en) | Geologic model via implicit function | |
| US20150066460A1 (en) | Stratigraphic function | |
| US20090043554A1 (en) | Updating a Subterranean Model Using at Least Electromagnetic Data | |
| US20180113235A1 (en) | Geologic Stratigraphy Via Implicit and Jump Functions | |
| CN113917560B (en) | Three-dimensional heavy magnetic electric shock multi-parameter collaborative inversion method | |
| US11249208B2 (en) | Geologic structural model generation | |
| Soleimani et al. | Integrated petrophysical modeling for a strongly heterogeneous and fractured reservoir, Sarvak Formation, SW Iran | |
| Erdoğan et al. | The conductivity structure of the Gediz Graben geothermal area extracted from 2D and 3D magnetotelluric inversion: Synthetic and field data applications | |
| Zhang et al. | 2D joint inversion of geophysical data using petrophysical clustering and facies deformation | |
| Wang et al. | 2D joint inversion of CSAMT and magnetic data based on cross-gradient theory | |
| Cai et al. | Effective 3D-transient electromagnetic inversion using finite-element method with a parallel direct solver | |
| Xiao et al. | Fast 2.5 D and 3D inversion of transient electromagnetic surveys using the octree-based finite-element method | |
| WO2023130074A1 (en) | Geologic modeling framework | |
| Wang et al. | 3-D crosswell electromagnetic inversion based on general measures | |
| CN115880455A (en) | Three-dimensional intelligent interpolation method based on deep learning | |
| Liu et al. | Performance investigations of auxiliary‐space Maxwell solver preconditioned iterative algorithm for controlled‐source electromagnetic induction problems with electrical anisotropy | |
| Marchant et al. | 3D inversion of electromagnetic logging-while-drilling data | |
| Ma et al. | Ant colony optimization inversion using the L1 norm in advanced tunnel detection | |
| Xu et al. | Bayesian prediction of potential depressions in the Erlian Basin based on integrated geophysical parameters | |
| Gao et al. | Research on reservoir inversion method of Paleogene in few well area | |
| Zhang et al. | Gravity inversion with a binary structure constraint imposed by magnetotelluric data and its application to a Carlin-type gold deposit |
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 |