CN113311494A - Satellite gravity field inversion method - Google Patents

Satellite gravity field inversion method Download PDF

Info

Publication number
CN113311494A
CN113311494A CN202110580140.0A CN202110580140A CN113311494A CN 113311494 A CN113311494 A CN 113311494A CN 202110580140 A CN202110580140 A CN 202110580140A CN 113311494 A CN113311494 A CN 113311494A
Authority
CN
China
Prior art keywords
gravity
matrix
inversion
roughness
grid
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
Application number
CN202110580140.0A
Other languages
Chinese (zh)
Other versions
CN113311494B (en
Inventor
鲁光银
张东兴
朱自强
张亮
毛雅静
郭晓旺
曹书锦
许刚
马致远
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Central South University
Original Assignee
Central South University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Central South University filed Critical Central South University
Priority to CN202110580140.0A priority Critical patent/CN113311494B/en
Publication of CN113311494A publication Critical patent/CN113311494A/en
Application granted granted Critical
Publication of CN113311494B publication Critical patent/CN113311494B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

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

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 satellite gravity field inversion method which includes the steps of obtaining gravity, gravity vectors, gravity gradient tensor data and corresponding kernel matrixes through equal parameter transformation, further generating corresponding depth weighting matrixes, constructing a multi-component joint inversion target function, and finally obtaining an inversion result. Compared with the traditional spherical shell unit, the dggid used by the method can improve the reuse capacity of the kernel function and reduce the memory occupation of the kernel matrix. Meanwhile, the method utilizes the characteristic that the dggard grids have consistent adjacency, and the central lines on three opposite sides of the dggard grids are used as the roughness matrix construction direction, so that the purposes of optimizing the roughness matrix construction and simplifying the calculation of the smoothness coefficient of the roughness matrix are achieved. In addition, in the inversion optimization solving process, the depth weighting matrix irrelevant to the component attenuation characteristic is used, the joint inversion setting of the components of the multiple force/gravity vector/gravity gradient tensor is optimized, the convergence speed is accelerated, and the inversion accuracy is improved.

Description

Satellite gravity field inversion method
Technical Field
The invention belongs to the technical field of geophysical exploration, and particularly relates to a satellite gravity field inversion method based on isoparametric transformation under global discrete grid spherical coordinates.
Background
In a spherical coordinate system, the physical property inversion is along the latitude
Figure BDA0003085742200000011
The three directions of longitude lambda and radial r are respectively divided into unit blocks at equal intervals, and the units are called Tesseroid unit bodies of longitude and latitude grids. The longitude and latitude grid has the advantages that the grid coding process is simple, and the transformation with the three-dimensional space coordinate is easy to realize. However, with the increase of the latitude, the longitude and latitude grids can generate a degradation phenomenon from the equator to the two poles in the spherical coordinate system, quadrangles near the poles can degrade into triangles, and a large amount of data redundancy is generated, so that the tesseroiid grids are difficult to realize global physical inversion.
In the prior art, the problem that Tesseroid units are equally spaced in longitude and latitude and are unequally spaced exists in latitude grid division, and the difference of corresponding kernel functions is obvious due to the influence of the potential field volume effect. In addition, the prior art constructs a roughness matrix W in a certain direction of the corresponding inversionmIn time, the smoothness coefficient needs to be recalculated to balance W according to different measuring region positionsr,
Figure BDA0003085742200000013
And WλThe difference caused by unequal length due to equal longitude and latitude intervals of the Tesseroid units. And the coefficient is influenced by the longitude and latitude of the observation grid, so that the inversion density model is easy to cause
Figure BDA0003085742200000012
The directional and lambda direction resolutions are poor.
In the prior art, in the gravity inversion, kappa is generally between 1.5 and 2.0; in the gravity gradient inversion, k is generally between 2.0 and 3.0; in the susceptibility inversion, κ ═ 3. This renders κ difficult to determine in the joint inversion of the multiple component data of the heavy magnetic bit field.
When the traditional tesseroid unit is inverted, all direction weight coefficients of the roughness matrix need to be debugged again in the longitude and latitude ranges of different measuring areas, and the implementation speed of the traditional inversion method in the longitude and latitude ranges of different measuring areas is reduced.
The longitude and latitude coordinate axes of the traditional tesseloid unit are warped, so that the process of calculating the smoothness coefficient of the roughness matrix is complex, the constraint capability of the fitting item of the inversion model is low, and the convergence speed is low, so that the resolution of the inversion result in the horizontal and depth directions cannot be ensured.
Based on the above description, a new satellite gravity field inversion method is needed to solve the problems that in the prior art, it is difficult to realize global physical property inversion by using tesseoid grids, the memory occupancy of a forward and backward evolution kernel matrix is large, and an inversion density model is in existence
Figure BDA0003085742200000021
The direction and lambda direction resolution ratio is poor, the inversion speed is slow, the horizontal and depth resolution ratio of the inversion result cannot be ensured, the multi-component joint inversion is difficult to effectively realize, and the inversion accuracy is low.
Disclosure of Invention
The invention aims to provide a satellite gravity field inversion method, which can reduce the memory occupation of a forward and backward kernel matrix and greatly improve the inversion speed and the inversion accuracy.
In order to solve the technical problems, the specific technical scheme of the satellite gravity field inversion method provided by the invention is as follows:
a satellite gravity field inversion method based on isoparametric transformation global discrete grid spherical coordinates comprises the following steps:
s1, determining the number and spatial position of the integral points of the hexagonal prism, constructing a regular hexagonal prism function, and calculatingThe weight coefficient of each integral point is calculated by using isoparametric transformation to obtain the abnormal response of each dggrid grid at each observation point, and the gravity and gravity gradient forward and backward error check matrix G is obtainedk；
S2, taking the normal and radial directions of three pairs of edges of the spherical shell hexagonal prism as the roughness matrix construction directions, optimizing the construction method, wherein the roughness matrix expression is as follows:
Figure BDA0003085742200000031
wherein p is norm order, and p is more than or equal to 0;
s3, constructing a depth weighting function Z:
Figure BDA0003085742200000032
s4, constructing a target function P of the multi-component joint inversion of gravity, gravity vector and gravity gradient tensora(m)：
Pa(m)＝φd+αφm→min
Wherein phi isdIs the norm of the difference between the observed data and the predicted data and the covariance matrix sigma of the observed data, i.e. the data functional, phimA stable functional or a model functional;
and S5, solving the optimal solution of the objective function to obtain an inversion result.
Preferably, in step S1, the gravity and gravity gradient forward/backward error check matrix GkThe following were used:
d＝Gm
wherein d is the length NdThe vector of (a); m is a physical property parameter vector, namely the density value of the physical property grid; g is Nd×NmPositive and negative evolution kernel matrix of size, NmAnd NdThe number of the subdivision grids and the number of the observation points are respectively.
Preferably, in step S2, the roughness matrix expression derivation process is:
the method for optimizing the construction of the dggrid grid by utilizing the characteristic that the dggrid grid has consistent adjacency and taking the central lines and the radial directions of three opposite sides of the dggrid grid as the construction directions of the roughness matrix is optimized, wherein the expression of the roughness matrix is as follows:
Figure BDA0003085742200000033
wherein, alSmoothness coefficient of the corresponding edge-to-edge roughness matrix, arSmoothness coefficient of radial roughness matrix, m is subsurface spatial density distribution, nlDispersing the normal vectors of the corresponding opposite sides along the normal directions of the three opposite sides, and introducing a reference geological model mrefThe following can be obtained:
Figure BDA0003085742200000041
wherein p is norm order, and p is more than or equal to 0.
Preferably, in step S4, the data functional φdComprises the following steps:
Figure BDA0003085742200000042
said model functional phimComprises the following steps:
Figure BDA0003085742200000043
wherein, 1 is not more than p<∞，
Figure BDA0003085742200000044
Is a k-order finite difference operator.
Preferably, a multi-component data fit is added to the data fit difference part, and the data fit term can be rewritten as:
Figure BDA0003085742200000045
in the formula, k is the number of components of the joint inversion of the components of the multiple force gradient tensor; scalkThe scaling factors of the inverted kernel matrix corresponding to the different components.
Preferably, D is taken as an identity matrix, a discrete gradient operator, or a discrete laplacian operator, corresponding to zeroth, first, and second order regularization terms, respectively.
The invention has the beneficial effects that:
the satellite gravity field inversion method provided by the invention fully utilizes the isoparametric transformation to convert a complex shape in a system coordinate system into a regular shape in local coordinates, and converts the integral interval of the spherical shell hexagonal prism into the integral interval of the regular hexagonal prism through the isoparametric transformation. The method has the following advantages:
1. in the prior art, the degradation phenomenon exists at two poles of a spherical coordinate system in a longitude and latitude grid, so that the tesserioid grid is difficult to realize the global physical inversion. However, the dggid unit used in this scheme does not have this problem.
2. The method comprises the steps of firstly obtaining gravity, gravity vector, gravity gradient tensor data and corresponding kernel matrix through isoparametric transformation. Secondly, as the spherical shell regular hexagonal prism grids are in symmetrical relation along each diagonal line and each edge line, compared with the traditional spherical shell unit forward modeling, the dggid used by the method can improve the repeated use capability of the kernel function and greatly reduce the memory occupation of the forward and backward modeling kernel matrix.
3. As the problem of unequal intervals due to longitude and latitude intervals in latitude grid division in the prior art is influenced by the volume effect of a position field, the difference of corresponding kernel functions is obvious. In addition, due to the fact that a smoothness coefficient needs to be introduced when a roughness matrix in a certain direction is constructed and inverted, and the coefficient is influenced by the longitude and latitude of an observation grid, extreme differences of an inverted density model are prone to occurring in the x direction and the y direction respectively. However, the distances between the spherical shell hexagonal prisms used in the scheme are equal, and the constructed roughness matrix has no obvious difference, so that the effect of inverting the density model is not influenced.
4. Compared with the traditional tesseroid unit inversion, the method does not need to debug the weight coefficients of the roughness matrix in all directions again in the longitude and latitude ranges of different measuring areas, so that the implementation speed of the inversion method provided by the invention in the measuring areas in different longitude and latitude ranges is accelerated.
5. According to the scheme, the grid of the dggrid has the characteristic of consistent adjacency, and the roughness matrix is constructed in the normal direction and the radial direction of three pairs of edges along the dggrid grid, so that the problem of buckling and warping of longitude and latitude coordinate axes of a traditional tesseroid unit is solved, the calculation of the smoothness coefficient of the roughness matrix is simplified, the constraint capacity of a fitting item of an inversion model is improved, the convergence speed is accelerated, and the resolution of an inversion result in the horizontal direction and the depth direction is ensured.
6. In the inversion optimization solving process, due to the fact that the depth weighting matrix irrelevant to the component attenuation characteristics is used, the joint inversion setting of the components of the multiple force/gravity vector/gravity gradient tensor is optimized, the convergence speed is accelerated, and the inversion accuracy is improved.
Drawings
FIG. 1 is a flow chart of the present invention;
FIG. 2 is a schematic diagram of the dggid grid forward model of the present invention;
FIG. 3 is a schematic diagram of a model of the inversion result of the dggrid grid of the present invention;
FIG. 4 is a graph fitting the dggid grid inversion data (gravity gradient tensor components) of the present invention;
FIG. 5 is a graph fitting the dggid grid inversion data of the present invention (gravity vector components);
FIG. 6 is an iterative convergence curve of the dggid grid inversion of the present invention.
Detailed Description
In order to better understand the method and steps of the present invention, the following describes the method for inverting the satellite gravitational field in detail with reference to the accompanying drawings.
As shown in fig. 1, the satellite gravity field inversion method of the present invention is based on isoparametric transformation global discrete network, and includes the following specific steps:
determining the thickness H of the spherical shell in the research area, and equally dividing the thickness H along the radius direction of the spherical shellIs nzAnd (2) constructing a first-level grid on the spherical surface (and further constructing an nth-level grid on the spherical surface) according to the definition of the global discrete grids, wherein 12 pentagonal global discrete grids are formed, 30 hexagonal global discrete grids are formed, and the spherical surface is used for discretizing the spherical shells of the research area layer by layer along the radial direction, namely, each layer of spherical shells is discretized into 12 global discrete grids of pentagonal prisms, and 30 global discrete grids of hexagonal prisms, which are referred to as global discrete grid spherical shell units.
And analyzing the relation between each observation point and each dggid grid coordinate distance r, and storing the relation in cor (i, j, r), wherein i corresponds to the observation point serial number and j is the physical grid serial number.
And determining the number and spatial position of the integral points of the hexagonal prism, further constructing a regular hexagonal prism function, calculating the weight coefficient of each integral point, and determining the integral formulas of each component of the gravity, the gravity vector and the gravity gradient tensor.
Optionally, the heavy magnetic satellite data or the earth and moon harmonic model data are selected to perform inversion of actual data; or selecting any one or more global discrete grid spherical shell units, utilizing a random number generation function to endow random residual density values to the selected global discrete grid spherical shell units so as to simulate a complex geological structure in the earth crust, utilizing an isoparametric transformation to calculate the gravity vector and gravity gradient tensor abnormal response of all the dggrid grids at each observation point, and carrying out theoretical research on an inversion method. The latter method is selected here, as shown in fig. 2, to generate forward data.
Optionally, different levels of white gaussian noise may be added to the data to evaluate the noise immunity of the inversion method.
Optionally, finally, the abnormal response of each dggrid grid at each observation point is obtained by utilizing equal parameter transformation calculation, and the gravity and gravity gradient forward and backward kernel matrix G is obtainedkAnd is computed only once for all the kernel functions that satisfy cer (i, j, r).
According to the principle of universal gravitation, a certain component kernel matrix G in multi-component joint inversionkKernel function
Figure BDA0003085742200000071
Is a function of the coordinate distance r between the observation point and the dggard grid, so that the observation point and the dggard grid with the same spatial relationship have the same kernel function, thereby improving the quality of the observation point and the dggard grid
Figure BDA0003085742200000072
The reuse capacity of the core matrix G is reducedkThe memory usage of.
The grid of the dggrid has the characteristic of consistent adjacency, the middle lines and the radial directions of three opposite sides of the grid are used as the construction directions of the roughness matrix, the construction method is optimized, and the expression of the roughness matrix is as follows:
Figure BDA0003085742200000081
wherein, alSmoothness coefficient of the corresponding edge-to-edge roughness matrix, arSmoothness coefficient of radial roughness matrix, m is subsurface spatial density distribution, nlDispersing the normal vectors of the corresponding opposite sides along the normal directions of the three opposite sides, and introducing a reference geological model mrefThe following can be obtained:
Figure BDA0003085742200000082
wherein, the norm order p can be made to be 2, and the smoothness coefficient a of the corresponding edge roughness matrixlSmoothness coefficient of radial roughness matrix a 1r＝10-6Normal vector n of the corresponding opposite sidelObtained for real-time calculation. Under the condition of no prior geological information, a reference geological model mref＝0。
Due to the different attenuation characteristics of gravity, gravity vector and gravity gradient tensor components, it is difficult to determine the attenuation coefficient k in the joint inversion of the multi-component data of the potential field. For this purpose, a sensitivity matrix is used
Figure BDA0003085742200000083
The main diagonal elements of (2) enable observation point attenuation characteristics above the physical property grid to be effectively reflected, thereby constructing a depth weighting matrix. For small scale geophysical models:
Figure BDA0003085742200000084
for large scale geophysical models, it is difficult to compute or store GTG, for this purpose, the above formula is rewritten as:
Figure BDA0003085742200000085
due to the fact that
Figure BDA0003085742200000086
Is a function of the distance r, then Z is calculated11And Z22Only exist when
Figure BDA0003085742200000087
And
Figure BDA0003085742200000088
the slight difference is:
Figure BDA0003085742200000089
Figure BDA00030857422000000810
…
Figure BDA0003085742200000091
wherein the coefficient ajFor each layer of kernel functions, Z can therefore be calculated quickly using catch-up.
Constructing a gravity, gravity vector and gravity gradient tensor multicomponent unionInverted objective function Pa(m)：
Pa(m)＝φd+αφm→min
Wherein phi isdIs the norm of the covariance matrix sigma of the observed data and the difference between the observed data and the predicted data, namely a data functional; phi is amData functiona phi for stabilizing the model objective function or model functiona-a priori constraints, i.e. model fitting terms, usually penalized by L1 or L2d：
Figure BDA0003085742200000092
Sum model functional phim：
Figure BDA0003085742200000093
Wherein, p may be made 2,
Figure BDA0003085742200000094
is a second order finite difference operator.
Because the abnormal amplitudes of the components such as the first derivative of the gravity, the gravity gradient tensor and the like are not consistent, only the composition of inversion data is considered in the traditional multi-component joint inversion, and the self characteristics of the components are not considered. For the above case, the multi-component data fitting is added to the data fitting difference part, and the data fitting term can be rewritten as:
Figure BDA0003085742200000095
wherein k is the number of components of the joint inversion of the components of the multi-force gradient tensor; scalkAnd scaling factors of the inversion kernel matrix corresponding to different components so as to balance and consider the problem of inconsistent abnormal amplitudes of the components.
As shown in fig. 3, finally, the optimal solution of the objective function is solved by using an optimization algorithm, so as to obtain an inversion result, and it can be seen from the inversion result that when a physical grid with a smaller density is filtered, the inversion density distribution and the density model are consistent, which indicates the correctness of the method provided by the present invention.
As shown in fig. 4, the residual error between the forward result and the data calculated by the inversion is small, which also indicates the correctness of the method provided by the present invention.
As shown in FIG. 5 and FIG. 6, it is the objective function P of the joint inversion of gravity, gravity vector and gravity gradient tensor multicomponenta(m) a convergence curve in the optimization iteration solving process shows that the inversion target function is stable in convergence after 6 times of inversion iteration, which shows that the method provided by the invention has the capability of fast convergence.
It will be understood by those skilled in the art that various changes may be made and equivalents may be substituted for elements thereof without departing from the spirit and scope of the invention. In addition, many modifications may be made to adapt a particular situation to the teachings of the invention without departing from the essential scope thereof. Therefore, it is intended that the invention not be limited to the particular embodiment disclosed, but that the invention will include all embodiments falling within the scope of the appended claims.

Claims (6)

1. A satellite gravity field inversion method is characterized by comprising the following steps:
s1, determining the number and spatial position of the integral points of the hexagonal prism, constructing a regular hexagonal prism function, calculating the weight coefficient of each integral point, obtaining the abnormal response of each dggrid grid at each observation point by using isoparametric transformation calculation, and obtaining the gravity and gravity gradient forward and reverse error checking matrix Gk；
S2, taking the normal and radial directions of three pairs of edges of the spherical shell hexagonal prism as the construction directions of the roughness matrix, and further optimizing the construction method, wherein the roughness matrix expression is as follows:
Figure FDA0003085742190000011
wherein p is norm order, and p is more than or equal to 0;
s3, constructing a depth weighting function Z:
Figure FDA0003085742190000012
s4, constructing a target function P of the multi-component joint inversion of gravity, gravity vector and gravity gradient tensora(m)：
Pa(m)＝φd+αφm→min
Wherein phi isdIs the norm of the difference between the observed data and the predicted data and the covariance matrix sigma of the observed data, i.e. the data functional, phimA stable functional or a model functional;
and S5, solving the optimal solution of the objective function to obtain an inversion result.
2. The satellite gravitational field inversion method of claim 1, wherein in step S1, the gravity and gravity gradient forward and backward kernel matrix GkThe following were used:
d＝Gm
wherein d is the length NdThe vector of (a); m is a physical property parameter vector, namely the density value of the physical property grid; g is Nd×NmPositive and negative evolution kernel matrix of size, NmAnd NdThe number of the subdivision grids and the number of the observation points are respectively.
3. The satellite gravity field inversion method of claim 1, wherein in step S2, the roughness matrix expression derivation process is:
the grid of the dggrid has the characteristic of consistent adjacency, the middle line and the radial direction of three opposite sides of the grid are used as the construction direction of the roughness matrix, and the construction method is optimized, wherein the expression of the roughness matrix is as follows:
Figure FDA0003085742190000021
wherein, alSmoothness coefficient of the corresponding edge-to-edge roughness matrix, arSmoothness coefficient of radial roughness matrix, m is subsurface spatial density distribution, nlDispersing the normal vectors of the corresponding opposite sides along the normal directions of the three opposite sides, and introducing a reference geological model mrefThe following can be obtained:
Figure FDA0003085742190000022
wherein p is norm order, and p is more than or equal to 0.
4. The satellite gravitational field inversion method of claim 1, wherein in step S4, said data functional ΦdComprises the following steps:
Figure FDA0003085742190000023
said model functional phimComprises the following steps:
Figure FDA0003085742190000024
wherein, 1 is not more than p<∞，
Figure FDA0003085742190000025
Is a k-order finite difference operator.
5. The satellite gravity field inversion method of claim 4, wherein a multi-component data fit is added to a data fit difference part, and the data fit term can be rewritten as:
Figure FDA0003085742190000026
in the formula, k is the number of components of the joint inversion of the components of the multiple force gradient tensor; scalkThe scaling factors of the inverted kernel matrix corresponding to the different components.
6. The satellite gravitational field inversion method of claim 4, wherein said D is taken as an identity matrix, a discrete gradient operator, or a discrete Laplace operator, corresponding to zero, first, and second order regularization terms, respectively.
CN202110580140.0A 2021-05-26 2021-05-26 Satellite gravity field inversion method Active CN113311494B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110580140.0A CN113311494B (en) 2021-05-26 2021-05-26 Satellite gravity field inversion method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110580140.0A CN113311494B (en) 2021-05-26 2021-05-26 Satellite gravity field inversion method

Publications (2)

Publication Number Publication Date
CN113311494A true CN113311494A (en) 2021-08-27
CN113311494B CN113311494B (en) 2022-04-26

Family

ID=77375088

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110580140.0A Active CN113311494B (en) 2021-05-26 2021-05-26 Satellite gravity field inversion method

Country Status (1)

Country Link
CN (1) CN113311494B (en)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114966878A (en) * 2022-04-12 2022-08-30 成都理工大学 Three-dimensional gravity field inversion method based on mixed norm and cross-correlation constraint
CN115220119A (en) * 2022-06-21 2022-10-21 广州海洋地质调查局 Gravity inversion method suitable for large-scale data
CN115236759A (en) * 2022-02-28 2022-10-25 中国人民解放军战略支援部队信息工程大学 Hexagonal grid subdivision method for determining earth gravity field
CN116205115A (en) * 2023-05-05 2023-06-02 南京航空航天大学 Structural form inversion precision improving method based on inverse element method and virtual-real combination technology

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20040201585A1 (en) * 2003-03-31 2004-10-14 Council Of Scientific And Industrial Research Generation of three dimensional fractal subsurface structure by Voronoi Tessellation and computation of gravity response of such fractal structure
CN108663718A (en) * 2018-06-14 2018-10-16 中国人民解放军61540部队 The least squares collocation model and modeling method of Satellite gravity field tensor diagonal line three-component inverting earth gravitational field
CN110045432A (en) * 2018-12-05 2019-07-23 中南大学 Gravitational field forward modeling method and 3-d inversion method under spherical coordinate system based on 3D-GLQ
CN110286416A (en) * 2019-05-13 2019-09-27 吉林大学 A kind of fast two-dimensional inversion of Density method based on physical property function
CN111400654A (en) * 2020-03-13 2020-07-10 中南大学 Gravity field rapid forward modeling method and inversion method based on Toplitz nuclear matrix
CN111999778A (en) * 2020-07-13 2020-11-27 国家海洋信息中心 Antarctic continental Mohuo surface depth inversion method based on satellite gravity gradient data

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20040201585A1 (en) * 2003-03-31 2004-10-14 Council Of Scientific And Industrial Research Generation of three dimensional fractal subsurface structure by Voronoi Tessellation and computation of gravity response of such fractal structure
CN108663718A (en) * 2018-06-14 2018-10-16 中国人民解放军61540部队 The least squares collocation model and modeling method of Satellite gravity field tensor diagonal line three-component inverting earth gravitational field
CN110045432A (en) * 2018-12-05 2019-07-23 中南大学 Gravitational field forward modeling method and 3-d inversion method under spherical coordinate system based on 3D-GLQ
CN110286416A (en) * 2019-05-13 2019-09-27 吉林大学 A kind of fast two-dimensional inversion of Density method based on physical property function
CN111400654A (en) * 2020-03-13 2020-07-10 中南大学 Gravity field rapid forward modeling method and inversion method based on Toplitz nuclear matrix
CN111999778A (en) * 2020-07-13 2020-11-27 国家海洋信息中心 Antarctic continental Mohuo surface depth inversion method based on satellite gravity gradient data

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
李红蕾等: "基于代数重构算法的重力梯度张量反演", 《地球物理学进展》 *
梁生贤: "互相关系数自约束的重力三维反演与高效求解", 《吉林大学学报(地球科学版)》 *

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115236759A (en) * 2022-02-28 2022-10-25 中国人民解放军战略支援部队信息工程大学 Hexagonal grid subdivision method for determining earth gravity field
CN115236759B (en) * 2022-02-28 2023-09-05 中国人民解放军战略支援部队信息工程大学 Hexagonal grid subdivision method for determining earth gravity field
CN114966878A (en) * 2022-04-12 2022-08-30 成都理工大学 Three-dimensional gravity field inversion method based on mixed norm and cross-correlation constraint
CN115220119A (en) * 2022-06-21 2022-10-21 广州海洋地质调查局 Gravity inversion method suitable for large-scale data
CN115220119B (en) * 2022-06-21 2023-02-24 广州海洋地质调查局 Gravity inversion method suitable for large-scale data
CN116205115A (en) * 2023-05-05 2023-06-02 南京航空航天大学 Structural form inversion precision improving method based on inverse element method and virtual-real combination technology
CN116205115B (en) * 2023-05-05 2023-08-25 南京航空航天大学 Structural form inversion precision improving method based on inverse element method and virtual-real combination technology

Also Published As

Publication number Publication date
CN113311494B (en) 2022-04-26

Similar Documents

Publication Publication Date Title
CN113311494B (en) Satellite gravity field inversion method
Gibou et al. A review of level-set methods and some recent applications
CN110045432B (en) Gravity field forward modeling method and three-dimensional inversion method under spherical coordinate system based on 3D-GLQ
CN111158059B (en) Gravity inversion method based on cubic B spline function
Aricò et al. MAST-2D diffusive model for flood prediction on domains with triangular Delaunay unstructured meshes
US10685482B2 (en) System and method for 3D restoration of complex subsurface models
CN113671570B (en) Seismic surface wave travel time and gravity anomaly joint inversion method and system
CN113239567A (en) Gravity field model assisted inverse distance weighting quasi-geoid grid interpolation method
Čunderlík et al. Numerical solution of the linearized fixed gravimetric boundary-value problem
CN114943178A (en) Three-dimensional geological model modeling method and device and computer equipment
CN106646645A (en) Novel gravity forward acceleration method
CN107561592A (en) A kind of density is the gravitational field forward modeling method of polynomial arbitrary polyhedron
CN109636912A (en) Tetrahedron subdivision finite element interpolation method applied to three-dimensional sonar image reconstruction
Galley et al. Geophysical inversion for 3D contact surface geometry
CN111679336B (en) Calculation method and device for bump gravity abnormal value of aviation gravity measurement point
CN113109883B (en) Method for forward modeling satellite gravity field based on isoparametric transformation global discrete grid spherical coordinates
Long et al. Three-dimensional forward modelling of gravity data using mesh-free methods with radial basis functions and unstructured nodes
Bai et al. 3-D non-linear travel-time tomography: imaging high contrast velocity anomalies
Yin et al. A fast 3D gravity forward algorithm based on circular convolution
Martyshko et al. Computationally Effective Gravity
CN114966878B (en) Three-dimensional gravity field inversion method based on mixed norm and cross-correlation constraint
Zhebel et al. Solving the 3D acoustic wave equation with higher-order mass-lumped tetrahedral finite elements
Gavete et al. A numerical comparison of two different approximations of the error in a meshless method
Mostafa Finite cube elements method for calculating gravity anomaly and structural index of solid and fractal bodies with defined boundaries
Burak et al. Determining density of regular grid for creating DTM using bicubic spline interpolation

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