CN112363236B - Gravity field data equivalent source continuation and data type conversion method based on PDE - Google Patents
Gravity field data equivalent source continuation and data type conversion method based on PDE Download PDFInfo
- Publication number
- CN112363236B CN112363236B CN202011103444.XA CN202011103444A CN112363236B CN 112363236 B CN112363236 B CN 112363236B CN 202011103444 A CN202011103444 A CN 202011103444A CN 112363236 B CN112363236 B CN 112363236B
- Authority
- CN
- China
- Prior art keywords
- gravity
- data
- equivalent source
- pde
- representing
- 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
- G01V7/00—Measuring gravitational fields or waves; Gravimetric prospecting or detecting
- G01V7/02—Details
- G01V7/06—Analysis or interpretation of gravimetric records
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 provides a gravity field data equivalent source continuation and data type conversion method based on PDE, which comprises the following steps: inputting gravity data and establishing a topographic relief curved surface; determining a space range of mesh generation according to the elevation information of the fluctuating observation curved surface and the preset inversion maximum depth, performing continuous structural non-uniform mesh generation, and determining an equivalent source inversion mesh space; according to the gravity background field, carrying out PDE three-dimensional inversion calculation with depth normalization factors, positive value constraint terms and normalization terms on gravity field data based on equivalent source inversion grid space to obtain a multilayer equivalent source model; different observation parameters are given, and the model is utilized to carry out three-dimensional forward modeling calculation on the gravity field of the nonlinear PDE so as to obtain one or more of gravity abnormal plumb data, gravity abnormal field data and gravity gradient tensor data. The invention has the beneficial effects that: the method can perform self-adaptive and high-precision continuation and data type conversion calculation on the underground gravity anomaly data in the complex environment.
Description
Technical Field
The invention relates to the technical field of geophysical surveying, in particular to a gravity field data equivalent source continuation and data type conversion method based on PDE.
Background
In gravity detection, people usually detect the vertical first derivative of a gravity field, but in actual data interpretation, the data often needs to be converted into required parameters and types, such as conversion of different observation heights (continuation) and gradient tensor data, and the like.
The prior art mainly comprises the following steps: (1) the equivalent source method utilizes a virtual field source to simulate the actual measurement abnormity and can be used for space continuation (including curved surface continuation), gradient calculation, component conversion and the like of the bit field data; selecting a single layer of equivalent sources and placing them near the surface is a key feature of the equivalent source method.
(2) A single-layer equivalent source is used for analyzing gravity field data, and the single-layer equivalent source is used for denoising gravity gradient tensor data; in order to reconstruct the potential field with high precision while ensuring the calculation efficiency, the multilayer equivalent source method is a relatively reasonable choice,
(3) the document "plum ends; aging; fipronil; carrying out roof beam green treatment; a velvet sea dragon; a gravity data processing method [ C ]//2017 annual meeting of China Union of geoscience, 0. "based on multilayer equivalent sources provides a method for dividing underground equivalent sources into a shallow layer, a middle layer and a deep layer to realize the conversion of gravity field data,
(4) the document "Chentao, Zhang Guibin, equivalent source technology [ J ] for calculating gravity gradient by using gravity anomaly.
The application of these technologies enables people to perform polarization and data type conversion operation of ground gravity data through equivalent source technology, but there are the following problems: 1) the construction of equivalent source grids is only three layers at most, and the grid division is discontinuous, in the prior art, the number of the set layers of equivalent sources is limited, and the grid division of an underground space is unreasonable, so that the precision and the speed of gravity field continuation and data type conversion are influenced; 2) the depth position of each equivalent layer needs to be estimated independently and then placed independently; 3) the equivalent source is calculated by adopting a linear forward modeling method and a linear model inversion method, and the anti-interference capability and the calculation accuracy are poor.
Disclosure of Invention
In view of the above, the technical problems to be solved by the present invention are: how to perform equivalent source extension and data type conversion on the gravity data of the irregular measurement position d and improve the calculation accuracy;
therefore, the invention provides a gravity data equivalent source continuation and data type conversion method based on PDE, which adopts continuous structural non-uniform mesh subdivision, introduces a depth normalization factor based on a forward and backward theoretical framework of nonlinear Differential PDE (Partial Differential Equations, and can also comprise finite elements, finite volumes and finite difference methods), directly determines the depth and distribution of an equivalent source in the inversion process, and can calculate the gravity field data of the whole three-dimensional free space at one time.
The invention provides a gravity data equivalent source continuation and data type conversion method based on PDE, which comprises the following steps:
s1, acquiring gravity field data d on the undulation observation curved surface0Establishing a topographic relief curved surface according to the topographic height information of the area where the gravity field data is located; the undulation observation curved surface is a recorded known value during observation;
s2, determining a space range of mesh subdivision according to the elevation information of the undulation observation surface and the preset inversion maximum depth, and performing continuous structural non-uniform mesh subdivision on the space range according to the undulation surface of the terrain to determine an equivalent source inversion mesh space;
s3, according to the gravity background field U0Inverting the grid space versus gravity field data d based on the equivalent source0Carrying out PDE three-dimensional inversion calculation with depth normalization factors, positive value constraint terms and normalization terms to obtain a multilayer equivalent source model of the gravity abnormal body;
s4, giving different observation parameters, and performing three-dimensional forward calculation of the gravity field of the nonlinear PDE by using the multilayer equivalent source model to obtain gravity abnormal plumb data and gravity abnormal field data U generated after the gravity abnormal body is extended at different observation positionss' and gravity gradient tensor data.
The technical scheme provided by the invention has the beneficial effects that: the technical scheme provided by the invention can be used for carrying out high-precision polar and data type conversion calculation on the gravity anomaly data generated by the underground gravity anomaly body obtained at the irregular measurement position on the undulating observation curved surface; the invention adopts continuous mesh subdivision, the number of equivalent source layers is usually more than 3, the calculation precision is higher, meanwhile, a depth normalization factor is added, the depth and the range of the equivalent layer do not need to be independently estimated, the gravity anomaly can be subjected to self-adaption, the required data after type conversion and continuation data can be rapidly, efficiently and accurately generated, and the stability and the precision are higher.
Drawings
FIG. 1 is a flowchart of a method for equivalent source extension and data type conversion based on gravity data of a PDE according to an embodiment of the present invention;
fig. 2 is a schematic diagram illustrating an effect of non-uniform mesh generation according to an embodiment of the present invention.
Detailed Description
In order to make the objects, technical solutions and advantages of the present invention more apparent, embodiments of the present invention will be further described with reference to the accompanying drawings.
Referring to fig. 1, an embodiment of the present invention provides a PDE-based gravity data equivalent source extension and data type conversion method, including the following steps:
s1, acquiring gravity field data d on the undulation observation curved surface0Establishing a topographic relief curved surface according to topographic height information of the area where the gravity field data is located; the gravity field data d0The data may be gravity abnormally vertical data or gravity gradient tensor data, and the gravity abnormally vertical data is taken as an example in this embodiment.
S2, determining a space range of mesh generation according to altitude information of the aerial survey undulation observation surface and the set inversion maximum depth, performing continuous structural non-uniform mesh generation on the space range according to the aerial survey undulation observation surface, and further determining an equivalent source inversion mesh space.
Specifically, the space range of the mesh generation includes an upper top surface and a lower bottom surface, wherein the upper top surface is determined according to the maximum height of the undulation observation curved surface, then the maximum depth of the gravity anomaly which can exist is estimated based on the existing detection technology or actual experience, and the inversion maximum depth is set according to the maximum depth, so that the lower bottom surface is determined;
after the spatial range of mesh generation is determined, dividing the spatial range by the lowest point of the topographic relief surface, carrying out uniform mesh generation on the spatial range above the lowest point to obtain a fine mesh, and carrying out non-uniform expanded mesh generation on the spatial range below the lowest point to obtain an expanded mesh; preferably, if the vertical edge of the fine mesh is 1 length unit (the specific value of the length unit may be set according to the size of the observation area space, for example, 100m is 1 length unit), the vertical edge of the extended mesh is increased at a speed 1.2 times that of the vertical edge of the fine mesh, and the maximum speed increase is set to 1.5 times, so that the calculation amount is reduced on the basis of ensuring a certain inversion accuracy. Please refer to fig. 2, which is a schematic diagram of a result of non-uniform mesh generation in the present embodiment, wherein a spatial range from a relief surface to a bottom surface forms an equivalent source inversion mesh space.
S3, according to the gravity background field U0Inversion of grid space gravity based on equivalent sourcesField data d0And carrying out PDE three-dimensional inversion calculation with depth normalization factors, positive value constraint terms and normalization terms to obtain a multilayer equivalent source model of the gravity abnormal body, wherein the multilayer equivalent source model specifically means that an inversion model solving space comprises a plurality of model depth surfaces, and based on the scheme, the number of model depth surfaces obtained by solving is usually more than 3 layers.
Preferably, the objective function of the PDE-based three-dimensional inversion calculation is:
Us=F(U0,m)
m≥0
where phi denotes an optimization objective (i.e., error),a numerical constraint representing the objective function is represented,representing model constraint of a target function, wherein m represents a density matrix of a multilayer equivalent source model to be solved, and considering that the physical property of an object leads the density to be a positive value, carrying out positive value constraint, namely m is more than or equal to 0; f (-) represents the three-dimensional PDE forward calculation of the multilayer equivalent source model, and can adopt a finite volume or finite element method, UsRepresenting gravity abnormal field data obtained by forward operation, and T (-) representing a conversion function from the gravity abnormal field data to the gravity abnormal plumb data; qx、Qy、QzRespectively representing interpolation functions in the north direction, the east direction and the vertical direction, wherein the interpolation functions comprise observation position information, and a Krigin interpolation function and the like can be selected; u shape0Representing the gravity background field intensity vector, U0x、U0y、U0zRespectively showing its north, east and vertical directionsA directional component; beta represents a regularization factor added according to actual requirements, and if not required, beta is 1; m isrefDensity matrix, W, representing a reference modelrRepresents the depth normalization factor:
wherein z represents the distance from the equivalent source to the relief surface, z0Representing the relief surface height, and r represents the depth factor, typically 3.
In this embodiment, the forward computation adopts a finite volume PDE method, and it should be noted that, in order to satisfy the finite volume solution condition, the equivalent source inversion grid space needs to be expanded, please refer to fig. 2, and the grid expansion is performed on the horizontal space below the set inversion maximum depth, above the undulation observation curved surface, and around the equivalent source inversion grid space according to the finite volume method, so as to meet the solution requirement.
In the forward calculation, the gravity field PDE equation is:
wherein U is U ═ U0+UsGamma denotes a gravitational constant, and ρ denotes the density of the multi-layered equivalent source model; based on the equation, gravity abnormal field data U is obtainedsForward calculation formula of (c):
Us=F(U0,m);
and (3) iteratively solving an objective function, namely minimizing the error phi, wherein a new density matrix m (namely the distribution of the density rho) is obtained after each iteration is completed and is used for fitting the measured gravity field data d0And finally obtaining a density matrix m of the multilayer equivalent source model. The forward calculation obtains three-component data, and the three-component data can be directly converted into plumb data or gradient tensor data in the inversion process so as to meet the requirements of different observation data types.
S4, given different observation parameters, using multiple layersThe equivalent source model performs three-dimensional forward calculation of the gravity field of the nonlinear PDE to obtain gravity abnormal plumb data and gravity abnormal field data U generated after the gravity abnormal body is extended at different observation positionssAnd one or more of gravity gradient tensor data.
Performing forward calculation according to the density matrix m determined in the step S3 to obtain extended gravity field data Us′:
Us′=F(U0′,m);
The formula for calculating the gravitational field tensor conversion is as follows:
wherein the content of the first and second substances,represents a gradient operator, U ═ Ux,Uy,Uz]TAnd representing the gravity field data subjected to the gravity field tensor conversion, wherein factors in a rightmost matrix of the equation are different tensors obtained by converting the gravity field data.
While the present invention has been described with reference to the embodiments shown in the drawings, the present invention is not limited to the embodiments, which are illustrative and not restrictive, and it will be apparent to those skilled in the art that various changes and modifications can be made therein without departing from the spirit and scope of the invention as defined in the appended claims.
Claims (7)
1. A gravity data equivalent source continuation and data type conversion method based on PDE is characterized by comprising the following steps:
s1, acquiring gravity field data d on the undulation observation curved surface0Establishing a topographic relief curved surface according to the topographic height information of the area where the gravity field data is located; the undulation observation curved surface is a recorded known value during observation;
s2, determining a space range of mesh subdivision according to the elevation information of the undulation observation surface and the preset inversion maximum depth, and performing continuous structural non-uniform mesh subdivision on the space range according to the undulation surface of the terrain to determine an equivalent source inversion mesh space;
s3, according to the gravity background field U0Inverting the grid space versus gravity field data d based on the equivalent source0Carrying out PDE three-dimensional inversion calculation with depth normalization factors, positive value constraint terms and normalization terms to obtain a multilayer equivalent source model of the gravity abnormal body;
s4, giving different observation parameters, and performing three-dimensional forward calculation of the gravity field of the nonlinear PDE by using the multilayer equivalent source model to obtain gravity abnormal plumb data and gravity abnormal field data U generated after the gravity abnormal body is extended at different observation positionss' and gravity gradient tensor data;
the number of layers of the model depth surface of the multilayer equivalent source model is more than 3;
in step S3, the objective function of the PDE three-dimensional inversion calculation is shown in formula (1):
Us=F(U0,m)
m≥0
phi denotes an optimization objective, i.e. an error;a numerical constraint representing an objective function;model constraints representing the objective function, m representing the value ofA density matrix of a multilayer equivalent source model solved by iterative optimization; m isrefA density matrix representing a reference model; f (-) represents the three-dimensional PDE forward calculation of the multilayer equivalent source model, UsRepresenting gravity abnormal field data obtained by forward operation, and T (-) representing a conversion function from the gravity abnormal field data to the gravity abnormal plumb data; qx、Qy、QzRespectively representing interpolation functions in the north direction, the east direction and the vertical direction, wherein the interpolation functions comprise observation position information; u shape0Representing the gravity background field intensity vector, U0x、U0y、U0zRespectively representing the north, east and vertical components; beta represents a preset weight parameter;
the depth normalization factor is as shown in formula (2):
in the formula (2), z represents the distance from the equivalent source to the relief surface, and z0Representing the relief surface height and r representing the depth factor.
2. The PDE-based gravity data equivalent source extension and data type conversion method according to claim 1, wherein in step S2, the spatial range of the mesh subdivision includes an upper top surface and a lower bottom surface, wherein the upper top surface is a plane determined by a maximum height of the undulation observation surface, and the lower bottom surface is a plane determined by a preset inversion maximum depth.
3. The PDE-based gravity data equivalent source continuation and data type conversion method according to claim 2, wherein in step S2, according to the topographic relief surface, continuous structured non-uniform mesh subdivision is performed on the spatial range to determine an equivalent source inversion mesh space, specifically:
dividing the space range according to the lowest point of the relief curved surface;
carrying out uniform mesh subdivision on the space range above the lowest point to obtain a fine mesh;
carrying out non-uniform mesh subdivision on the space range below the lowest point to obtain an expanded mesh;
and the spatial range from the topographic relief surface to the lower bottom surface, namely the fine grid and the expanded grid form an equivalent source inversion grid space.
4. The PDE-based gravity data equivalent source continuation and data type conversion method according to claim 3, wherein said vertical edge of said extended grid is alpha to said vertical edge of said fine grid1The speed is increased by multiple times and the maximum speed increase is set to alpha2Wherein α is2>α1>1。
5. The method for PDE-based gravity data equivalent source extension and data type conversion according to claim 1, wherein in step S3, the density matrix m and the gravity background field U are generated iteratively during optimization0Performing PDE forward calculation to obtain gravity anomaly field data generated by the gravity anomaly body:
Us=F(U0,m)。
6. the method for PDE-based gravity data equivalent source extension and data type conversion according to claim 1, wherein in step S4, forward calculation is performed according to the density matrix m determined in step S3 to obtain extended gravity anomaly field data Us′:
Us′=F(U0′,m)。
7. The PDE-based gravity data equivalent source extension and data type conversion method according to claim 1, wherein the calculation of the gravity field tensor data is as follows (3):
wherein the content of the first and second substances,represents a gradient operator, U ═ Ux,Uy,Uz]TAnd (3) representing the gravity field data subjected to the gravity field tensor conversion, wherein the factors in the rightmost matrix of the formula (3) are different tensor data obtained by converting the gravity field data.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202011103444.XA CN112363236B (en) | 2020-10-15 | 2020-10-15 | Gravity field data equivalent source continuation and data type conversion method based on PDE |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202011103444.XA CN112363236B (en) | 2020-10-15 | 2020-10-15 | Gravity field data equivalent source continuation and data type conversion method based on PDE |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112363236A CN112363236A (en) | 2021-02-12 |
CN112363236B true CN112363236B (en) | 2022-04-01 |
Family
ID=74506805
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202011103444.XA Active CN112363236B (en) | 2020-10-15 | 2020-10-15 | Gravity field data equivalent source continuation and data type conversion method based on PDE |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112363236B (en) |
Families Citing this family (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112949133B (en) * | 2021-03-08 | 2024-03-22 | 中国地质大学(武汉) | Heavy magnetic joint inversion method based on PDE |
CN113587921B (en) * | 2021-05-31 | 2024-03-22 | 中国人民解放军61540部队 | Gravity gradient field and gravity anomaly field submersible vehicle fusion positioning method and system |
CN113627051B (en) * | 2021-07-23 | 2024-01-30 | 中国地质科学院地球物理地球化学勘查研究所 | Gravity anomaly field separation method, system, storage medium and electronic equipment |
CN114167511B (en) * | 2021-11-26 | 2023-08-11 | 兰州大学 | Bit field data rapid inversion method based on continuous expansion downward continuation |
CN114814967B (en) * | 2022-04-26 | 2024-04-02 | 中国人民解放军61540部队 | High-resolution submarine topography nonlinear method for inversion of local sea area disturbance gravity data |
CN116047617B (en) * | 2023-03-10 | 2023-06-27 | 中国地质科学院地球物理地球化学勘查研究所 | Method and device for identifying geological features between wells |
Family Cites Families (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP3033689B1 (en) * | 2013-04-24 | 2021-09-15 | The University of British Columbia | A penalty method for pde-constrained optimization |
WO2019013906A1 (en) * | 2017-07-13 | 2019-01-17 | Exxonmobil Upstream Research Company | Visco-pseudo-elastic tti fwi/rtm formulation and implementation |
CN107402409A (en) * | 2017-09-26 | 2017-11-28 | 西南石油大学 | A kind of three-dimensional irregular stratum fluctuating interface gravity forward modeling method |
CN109375280B (en) * | 2018-12-10 | 2020-03-31 | 中南大学 | Gravity field rapid high-precision forward modeling method under spherical coordinate system |
CN110221344B (en) * | 2019-06-17 | 2020-08-28 | 中国地质大学(北京) | Seismic full-waveform and gravity joint inversion method for three-dimensional density structure of earth crust |
CN110389391B (en) * | 2019-08-01 | 2020-12-15 | 自然资源部第二海洋研究所 | Heavy magnetic bit field analytic extension method based on spatial domain |
CN111142169A (en) * | 2020-02-25 | 2020-05-12 | 中国地质大学(北京) | Submarine topography inversion method based on gravity gradient data |
-
2020
- 2020-10-15 CN CN202011103444.XA patent/CN112363236B/en active Active
Also Published As
Publication number | Publication date |
---|---|
CN112363236A (en) | 2021-02-12 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN112363236B (en) | Gravity field data equivalent source continuation and data type conversion method based on PDE | |
CN105549106B (en) | A kind of gravity multiple solutions inversion method | |
CN112147709B (en) | Gravity gradient data three-dimensional inversion method based on partial smoothness constraint | |
CN111856599B (en) | Magnetic measurement data equivalent source pole and type conversion method based on PDE | |
US9915742B2 (en) | Method and system for geophysical modeling of subsurface volumes based on label propagation | |
CN111856598B (en) | Magnetic measurement data multilayer equivalent source upper extension and lower extension method | |
WO2009092992A1 (en) | Geophysical data processing systems | |
WO2014099204A1 (en) | Method and system for geophysical modeling of subsurface volumes based on computed vectors | |
WO2014099201A1 (en) | Geophysical modeling of subsurface volumes based on horizon extraction | |
CN109902315B (en) | Method for delineating deep boundary of hidden granite rock mass | |
CN112346139B (en) | Gravity data multilayer equivalent source continuation and data conversion method | |
CN114943178A (en) | Three-dimensional geological model modeling method and device and computer equipment | |
CN112526625B (en) | Computing device for abnormal value of Bragg gravity of aviation gravity measurement point | |
CN115437027A (en) | Method and device for calculating bump gravity anomaly by using geological information variable density forward modeling | |
CN111859251B (en) | Magnetic measurement data equivalent source extension and extension method based on PDE | |
CN109444973B (en) | Gravity forward modeling acceleration method under spherical coordinate system | |
US20160334537A1 (en) | System and method for two dimensional gravity modeling with variable densities | |
CN107942374A (en) | Diffracted wave field extracting method and device | |
CN112666614B (en) | Debris flow source static reserve calculation method based on electrical prospecting and digital elevation model | |
EP3281044B1 (en) | Method for estimating elastic parameters of subsoil | |
CN111208558B (en) | Method and device for establishing ultra-deep low-amplitude three-dimensional geological structure | |
CN112748471B (en) | Gravity-magnetic data continuation and conversion method of unstructured equivalent source | |
WO2012021938A1 (en) | A method of analysing data obtained using a gravity gradiometer | |
CN107730586B (en) | Method and system for modeling stratum | |
CN114966878B (en) | Three-dimensional gravity field inversion method based on mixed norm and cross-correlation constraint |
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 |