Adaptive de-weighting high-precision joint inversion method driven by electromagnetic gradient data

Abstract

The application discloses an adaptive weight-adjusting high-precision joint inversion method based on electromagnetic gradient data driving and relates to the technical field of electromagnetic signal inversion. The application comprises the following steps: calculating the spatial gradient data of the original electromagnetic field observation data by using a Gaussian weighted neighborhood method; combining a conventional data residual term, a data residual term and a regularization term to construct a joint objective function; embedding a gradient sensitive scale matrix into the data residual term of the joint objective function to establish an inversion equation; iteratively calculating the inversion equation and updating model parameters to obtain a three-dimensional inversion result reflecting the spatial distribution of resistivity. The application extends the classic Gauss-Newton method, optimizes the inversion by using the residual terms of the electromagnetic field data and the spatial gradient, introduces the scale matrix sensitive to the gradient in the residual weighting, so that only the model parameters related to the abnormal area can be updated in each inversion iteration, and the model parameters are mapped to the three-dimensional finite element domain for forward calculation.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of electromagnetic signal inversion, in particular to a self-adaptive weight adjustment high-precision joint inversion method based on electromagnetic gradient data driving. BACKGROUND

[0002] The methods for inverting underground resistivity structure based on geophysical electromagnetic data mainly include three types: traditional smooth inversion method, three-dimensional fine grid inversion method and regularization inversion method. Although the three methods can all realize the inversion from electromagnetic data to underground resistivity distribution, they perform differently in terms of inversion accuracy, calculation efficiency, adaptability to complex geology and the like. The traditional smooth inversion method is based on the least square principle and takes the model roughness (or smoothness) as a regularization constraint, and seeks a most smooth resistivity model while fitting the observation data. This method has strong anti-noise ability and high calculation efficiency, and is effective in many scenarios. However, in order to pursue model smoothness, it generally smooths out sharp geological boundaries, and thus is difficult to accurately depict complex geological structures or abrupt interfaces. The three-dimensional fine grid inversion method restores complex three-dimensional geological structures with high resolution and strong adaptability by discretizing the underground space with fine grids. However, due to the need to invert millions of grid cells, it leads to huge memory consumption and calculation cost, and thus is difficult to extend to large-scale calculation domains. In addition, in order to improve the boundary preservation ability, alternative regularization strategies have been developed, including focused regularization and total variation regularization. The focused regularization restores the clear geoelectric structure by introducing the minimum gradient support functional. The total variation regularization preserves the discontinuity of the conductivity by imposing gradient sparsity constraints, so as to produce a piecewise constant solution. Although these technologies effectively enhance the clarity of the boundary, they usually increase the nonlinearity and calculation complexity, which may lead to slower convergence speed and more sensitive parameter adjustment.

[0003] Although great progress has been made in geophysical electromagnetic inversion technology, the existing methods are difficult to simultaneously realize the consistency of physical principles, high boundary resolution, high calculation and acquisition efficiency without the support of multi-physical field integration or repeated measurement. To solve this limitation, the conventional way currently adopted is to integrate the electromagnetic data residual term into the overall inversion residual term, and adjust the weight of the electromagnetic field data and gradient field data residual term by a single constant, to enhance the model convergence speed and structural clarity by using the high sensitivity of electromagnetic field gradient to discontinuous conductivity. However, if a uniform constant is used to constrain the weight, it means that each non-uniform gradient term is given the same regularization strength, which leads to insufficient identification accuracy of abnormal boundaries in the inversion, low noise tolerance in the abnormal area, and easy occurrence of false oscillation or pseudo-anomaly in the model, low model stability and reliability. SUMMARY

[0004] In view of the above existing problems, the present application is proposed.

[0005] Therefore, the present application provides a self-adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, which solves the problem of improving the calculation efficiency while realizing high-resolution structure imaging.

[0006] To solve the above technical problems, the present application provides the following technical solutions:

[0007] The present application provides a self-adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, which includes calculating the spatial gradient data of the original electromagnetic field observation data using the Gaussian weighted neighborhood method; calculating the spatial gradient difference of adjacent frequencies to form the frequency difference spatial gradient; normalizing and combining the frequency difference spatial gradient of all measurement points to obtain the overall gradient intensity normalization matrix; processing the overall gradient intensity normalization matrix using a contrast weighting function with an adjustable exponent to construct the gradient sensitive scale matrix; combining the conventional data residual term, the data residual term and the regularization term to construct the joint objective function; embedding the gradient sensitive scale matrix into the data residual term of the joint objective function to establish the inversion equation; iteratively calculating the inversion equation and updating the resistivity vector of the resistivity model at each grid in the underground to obtain the three-dimensional inversion result reflecting the spatial distribution of the resistivity.

[0008] As a preferred scheme of the self-adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, the specific steps of calculating the spatial gradient data of the original electromagnetic field observation data using the Gaussian weighted neighborhood method are as follows,

[0009] Obtain the electromagnetic field data at the measurement point and the electromagnetic field data at the neighborhood measurement point;

[0010] Calculate the gradient component of each neighborhood point pair according to the spatial distance between the electromagnetic field data at the measurement point under the frequency condition and the electromagnetic field data at the neighborhood measurement point under the same frequency condition.

[0011] Calculate and normalize the distance between the neighborhood measurement points by an exponential function to obtain the normalized Gaussian weight, and use the normalized Gaussian weight to weight process the gradient component of each neighborhood point pair to generate the spatial gradient data.

[0012] As a preferred scheme of the self-adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, the specific steps of calculating the spatial gradient difference of adjacent frequencies to form the frequency difference spatial gradient are as follows,

[0013] Obtain the spatial gradient at the measurement point under the frequency condition and the spatial gradient under the previous frequency condition at the measurement point, and select a frequency function matching the frequency sampling scheme.

[0014] The difference between the frequency at the measuring point and the spatial gradient at the previous frequency at the measuring point is calculated.

[0015] The difference between the frequency at the measuring point and the frequency function value at the previous frequency at the measuring point is calculated.

[0016] The ratio of the difference of the spatial gradient and the difference of the frequency function value is taken as the frequency difference spatial gradient.

[0017] As a preferred scheme of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, the normalized gradient intensity matrix is obtained by normalizing and combining the frequency difference spatial gradients of all measuring points, and the specific steps are as follows,

[0018] The frequency difference spatial gradients of each measuring point are collected to determine the maximum value of the frequency difference spatial gradients of all measuring points.

[0019] Based on the maximum value, the frequency difference spatial gradients of each measuring point are normalized to obtain the normalized frequency difference spatial gradient vector of each measuring point.

[0020] The normalized frequency difference spatial gradient vectors of all measuring points are combined to form the overall gradient intensity normalization matrix.

[0021] As a preferred scheme of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, the specific steps of constructing the gradient sensitive scale matrix are as follows,

[0022] The overall gradient intensity normalization matrix is taken as the input parameter of the contrast weighted function.

[0023] Each element in the overall gradient intensity normalization matrix is processed using a contrast weighted function with adjustable exponent.

[0024] The processed elements are arranged in the form of a diagonal matrix to form the gradient sensitive scale matrix, and the gradient sensitive scale matrix is used to adjust the contribution weight of individual data in the data residual term.

[0025] As a preferred scheme of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, the specific steps of combining the conventional data residual term, the data residual term and the regularization term to construct the joint objective function are as follows,

[0026] The conventional data residual term is constructed according to the difference between the observed electromagnetic field data and the predicted electromagnetic field data of the forward simulation.

[0027] The conventional data residual term, the data residual term and the regularization term are combined by weighting to construct the joint objective function.

[0028] As a preferred scheme of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving provided by the application, the data residual term comprises,

[0029] a processing result of the difference between the diagonal electromagnetic field data weighting matrix and the observed electromagnetic field data and the predicted electromagnetic field data;

[0030] a processing result of the difference between the diagonal gradient field data weighting matrix and the observed gradient data and the predicted gradient data;

[0031] a regulation result of the gradient sensitive scale matrix on the influence degree of a single data in the gradient data residual.

[0032] As a preferred scheme of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving provided by the application, the inversion equation comprises a Jacobian matrix of the predicted electromagnetic field data and a Jacobian matrix of the predicted spatial gradient data.

[0033] As a preferred scheme of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving provided by the application, the data residual term of the joint objective function is embedded in the gradient sensitive scale matrix, and the inversion equation is established, and the specific steps are as follows,

[0034] According to the Jacobian matrix of the predicted electromagnetic field data, the Jacobian matrix of the predicted spatial gradient data and the gradient sensitive scale matrix, a coefficient matrix of the inversion equation is constructed;

[0035] According to the Jacobian matrix of the predicted electromagnetic field data, the Jacobian matrix of the predicted spatial gradient data, the data residual term and the regularization term, a gradient direction of the inversion equation is constructed;

[0036] The coefficient matrix and the gradient direction are combined to establish the inversion equation.

[0037] As a preferred scheme of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving provided by the application, the inversion equation is iteratively calculated and the model parameters are updated to obtain a three-dimensional inversion result reflecting the spatial distribution of resistivity, and the specific steps are as follows,

[0038] The resistivity model update step value of the current iteration is calculated, and the resistivity model parameters are corrected based on the resistivity model update step value;

[0039] According to the corrected resistivity model parameters, a three-dimensional inversion result reflecting the spatial distribution of resistivity is output.

[0040] The present application has the beneficial effects that: by extending the classical Gauss-Newton method, the residual term of the electromagnetic field data and the spatial gradient is used for inversion optimization, and the scale matrix sensitive to the gradient is introduced in the residual weighting, so that only the model parameters related to the abnormal area can be updated in each inversion iteration, and the forward calculation is performed in the three-dimensional finite element domain. This mechanism improves the boundary resolution while significantly reducing the calculation cost and can completely capture the three-dimensional physical response, overcoming the limitations of the traditional constant constraint gradient residual term method in abnormal area boundary identification accuracy and background area stability due to the neglect of grid non-uniformity and physical property spatial heterogeneity, and providing a practical solution for electromagnetic detection in complex geological environment. BRIEF DESCRIPTION OF DRAWINGS

[0041] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the drawings needed in the embodiment description will be briefly introduced. Obviously, the drawings in the following description are only some embodiments of the present application, and other drawings can be obtained by those skilled in the art without creative labor.

[0042] Figure 1 The flowchart of the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving;

[0043] Figure 2 The electromagnetic gradient data and the electromagnetic field data response comparison and their sensitivity to boundary changes, wherein, Figure 2 (a) represents the true three-dimensional model diagram of the electromagnetic gradient field containing the conductive belt anomaly; Figure 2 (b) represents the three-dimensional model diagram of the electromagnetic field data response of the model in Figure 2 (a) with slight displacement and size change; Figure 2 (c) is the relative difference diagram of the electromagnetic field data response between Figure 2 (a) and Figure 2 (b); Figure 2 (d) is the relative difference diagram of the spatial gradient data response between Figure 2 (a) and Figure 2 (b); Figure 2 (e) is the relative difference diagram of the electromagnetic field and gradient field responses of the two models. DETAILED DESCRIPTION

[0044] In order to make the above-mentioned purposes, features and advantages of the present application more obvious and easy to understand, the specific embodiments of the present application will be described in detail below with reference to the drawings of the specification.

[0045] Many specific details are set forth in the following description in order to provide a full understanding of the invention. However, the invention may also be practiced in other ways different from those described herein, and those skilled in the art can make similar extensions without departing from the spirit of the invention. Therefore, the invention is not limited to the specific embodiments disclosed below.

[0046] Secondly, the term "one embodiment" or "embodiment" as used herein refers to a specific feature, structure, or characteristic that may be included in at least one implementation of the present invention. The phrase "in one embodiment" appearing in different places in this specification does not necessarily refer to the same embodiment, nor is it a single or selective embodiment that is mutually exclusive with other embodiments.

[0047] Reference Figures 1-2 This is one embodiment of the present invention, which provides an adaptive weighted high-precision joint inversion method based on electromagnetic gradient data driving, comprising the following steps:

[0048] S1. Calculate the spatial gradient data of the original electromagnetic field observation data using the Gaussian weighted neighborhood method.

[0049] S1.1: Obtain electromagnetic field data at the measuring point and electromagnetic field data at neighboring measuring points;

[0050] It should be noted that the resistivity model response is obtained through forward modeling or through actual detection, and electromagnetic field data at the measuring point and the electromagnetic field data at neighboring measuring points are obtained.

[0051] S1.2: Based on the electromagnetic field data at the measurement point under the same frequency conditions and the spatial distance between the electromagnetic field data at neighboring measurement points under the same frequency conditions, calculate the gradient components of each pair of neighboring points, defined as: ;

[0052] in, Indicates at the measuring point frequency Electromagnetic fields under certain conditions; Indicates at the measuring point Frequency conditions Electromagnetic field under; This represents the x-coordinate of measurement point j in the three-dimensional coordinate system; This represents the x-coordinate of the measuring point p in the three-dimensional coordinate system.

[0053] S1.3: The distance between neighboring measurement points is calculated and normalized using an exponential function to obtain normalized Gaussian weights. The gradient components of each neighboring point pair are weighted using the normalized Gaussian weights to generate spatial gradient data.

[0054] Spatial gradient data is defined as:

[0055] ;

[0056] wherein, is the spatial gradient at the measurement point at frequency ; represents the electromagnetic field at the measurement point at frequency ; represents the electromagnetic field at the measurement point at frequency ; is the set of adjacent measurement points of the measurement point ; is the normalized Gaussian weight.

[0057] In order to simultaneously retain the high noise suppression ability and high boundary definition of the resistivity model, the normalized Gaussian weight expression is:

[0058] ;

[0059] wherein, is the normalized Gaussian weight; represents the x-direction vector at the measurement point in the three-dimensional coordinate system; represents the x-direction vector at the measurement point in the three-dimensional coordinate system; is a weight factor used to adjust the noise suppression ability and boundary sensitivity when calculating the electromagnetic field gradient, when , the Gaussian function is extremely wide, and the weight difference of different points is smoothed, so the weight tends to be uniformly averaged, and the noise suppression ability is enhanced; when , the Gaussian function is extremely narrow, and only the local points have significant weights, so the local sensitivity is extremely prominent, and the boundary sensitivity and definition are enhanced.

[0060] S2, calculate the spatial gradient difference of adjacent frequencies, form the frequency difference spatial gradient, normalize and combine the frequency difference spatial gradient of all measurement points to obtain the overall gradient intensity normalization matrix, and process the overall gradient intensity normalization matrix using a contrast weighting function with adjustable exponent to construct the gradient sensitivity scale matrix.

[0061] S2.1: set the transmission frequency according to the measurement requirement, and use the Gaussian weighted neighborhood method to calculate the electromagnetic field data, obtain the spatial gradient at the measurement point at the frequency and at the previous frequency at the measurement point, and select a mathematical function matched with the frequency sampling scheme as the frequency function according to the frequency sampling mode and the analysis purpose (such as accuracy, resolution);

[0062] S2.2: calculate the difference value of the spatial gradient at the measurement point at the frequency and at the previous frequency at the measurement point, defined as: ;

[0063] in, For measuring points frequency Spatial gradient below;

[0064] It should be noted that, The spatial gradient of frequency at the measurement point. This represents the spatial gradient at the previous frequency at the measurement point.

[0065] S2.3: Calculate the difference between the frequency at the measuring point and the previous frequency at the measuring point, defined as: ;

[0066] in, It is a function of frequency;

[0067] It should be noted that, The frequency function value at the measurement point. This is the frequency function value of the previous frequency at the measurement point.

[0068] S2.4: The ratio of the difference in spatial gradients to the difference in frequency function values ​​is taken as the frequency difference spatial gradient, expressed as:

[0069] ;

[0070] in, For measuring points frequency The gradient in the difference space below, For measuring points frequency Spatial gradient below; The function is a function of frequency. The appropriate function is selected based on the chosen frequency sampling scheme and the requirements for accuracy and resolution (e.g., when using a logarithmic interval sampling frequency). = ).

[0071] It should be noted that frequency normalization suppresses background response, significantly improving boundary sensitivity.

[0072] S2.5: Collect the frequency differential spatial gradient at each measuring point, and determine the maximum value of the frequency differential spatial gradient at all measuring points. The expression is:

[0073] ;

[0074] in, Indicates the measuring point Located in the vertical direction Spatial gradient vectors at each grid point ; The preset minimum number is used to avoid division by zero error. The total number of grids in the longitudinal direction of the measurement point j.

[0075] S2.6: The frequency difference spatial gradient of each measurement point is normalized based on the maximum value to obtain the normalized frequency difference spatial gradient vector of each measurement point, and the expression is:

[0076]

[0077] It should be noted that the frequency difference spatial gradient of each measurement point is normalized by the ratio of the frequency difference spatial gradient of each measurement point to the maximum value to generate the normalized frequency difference spatial gradient vector.

[0078] S2.7: The normalized frequency difference spatial gradient vectors of all measurement points are combined to form the overall gradient intensity normalization matrix.

[0079] The normalized matrix of the gradient intensity of all measurement points is obtained by combining the frequency difference spatial gradients of all measurement points.

[0080]

[0081] Wherein, r is the overall gradient intensity normalization matrix, respectively represent the normalized frequency difference spatial gradient vector of the measurement point 1, the measurement point 2, the measurement point j, and the measurement point N2; is the total number of measurement points.

[0082] S2.8: The overall gradient intensity normalization matrix is used as the input parameter of the contrast weighting function;

[0083]

[0084] Wherein, is an adjustable index used to control the weight between the boundary sensitivity and the stability of the smooth area;

[0085] S2.9: The contrast weighting function with adjustable index is used to process each element in the overall gradient intensity normalization matrix.

[0086] The contrast weighting function with adjustable index is introduced to calculate each element in the overall gradient intensity normalization matrix, which is defined as:

[0087] S2.10: The processed elements are arranged in the form of a diagonal matrix to form a gradient sensitive scale matrix, and the gradient sensitive scale matrix is used to adjust the contribution weight of a single data in the data residual term.

[0088] ​​​​​​​​​​​​​​​​​​​​​​​​​​​​​It should be noted that the gradient sensitivity scale matrix is used to adjust the contribution weight of individual data in the electromagnetic data residual term, to enhance the sensitivity at the electrical boundary to improve the resolution, and to maintain the smooth constraint in the homogeneous area to suppress the noise.

[0089] S3, combining the conventional data residual term, the data residual term and the regularization term, to construct a joint objective function.

[0090] S3.1: According to the difference between the observed electromagnetic field data and the forward simulation predicted electromagnetic field data, the conventional data residual term is constructed as: ;

[0091] wherein, is a diagonal electromagnetic field data weighting matrix; is the observed electromagnetic field data, is the predicted electromagnetic field data; is the observed electromagnetic gradient data, is the predicted electromagnetic gradient data; is a diagonal gradient field data weighting matrix;

[0092] It should be noted that each element value of the diagonal electromagnetic field data weighting matrix corresponds to the inverse of the standard deviation of the noise at each observed electromagnetic field data position; each element value of the diagonal gradient field data weighting matrix corresponds to the inverse of the standard deviation of the noise at each observed electromagnetic gradient data position;

[0093] The predicted electromagnetic field data and the predicted electromagnetic gradient data are obtained by forward simulation, which refers to the simulation of physical laws in the current iteration resistivity model to generate predicted electromagnetic field data, and further derive predicted electromagnetic gradient data.

[0094] S3.2: Combining the conventional data residual term, the data residual term and the regularization term by weighting, to construct a joint objective function.

[0095] The expression of the joint objective function is:

[0096] ;

[0097] wherein, is a resistivity model vector; is a data residual term; is a regularization term, used to impose prior constraints and enhance the smoothness of the resistivity model; is used to balance the weight of the two terms.

[0098] To solve the problem of weak resistivity model updating and slow convergence speed in traditional EM inversion, the sensitivity of the gradient field to the small perturbation of the boundary position is used, and an additional term of the spatial electromagnetic field gradient is introduced into the data residual term to restore the convergence speed and structural clarity of the resistivity model. The data residual term after introducing the spatial electromagnetic field gradient additional term is defined as:

[0099] ;

[0100] wherein, is a diagonal electromagnetic field data weighting matrix; is the observed electromagnetic field data, is the predicted electromagnetic field data; is the observed electromagnetic gradient data, is the predicted electromagnetic gradient data; is a diagonal gradient field data weighting matrix; is a weight parameter, used to adjust the weight of the overall data residual term to the contribution of the data residual term; is a gradient sensitive scale matrix, used to adjust the weight of the individual data to the contribution of the data residual term.

[0101] S3.3: The processing result of the difference between the diagonal electromagnetic field data weighting matrix and the observed electromagnetic field data and the predicted electromagnetic field data is: ;

[0102] It should be noted that the further refinement of the conventional data residual term is described, and it is explicitly pointed out that the conventional data residual term is the processing result of the difference between the diagonal electromagnetic field data weighting matrix and and .

[0103] S3.4: The processing result of the difference between the diagonal gradient field data constraint matrix and the observed gradient data and the predicted gradient data is: ;

[0104] It should be noted that, the core composition of the conventional data residual term is described, and it is explicitly pointed out that the conventional data residual term is the processing result of the difference between the diagonal gradient field data weighting matrix and and .

[0105] S3.5: The adjustment result of the gradient sensitive scale matrix to the influence degree of individual data in the gradient data residual.

[0106] It should be noted that, is a gradient sensitive scale matrix, used to adjust the weight of the individual data to the contribution of the data residual term in the conventional data residual term.

[0107] S4. The gradient sensitivity scale matrix is embedded in a data residual term of a joint objective function, an inversion equation is established, the inversion equation is iteratively calculated and resistivity model parameters are updated to obtain a three-dimensional inversion result reflecting the spatial distribution of resistivity.

[0108] S4.1: The inversion equation includes a Jacobian matrix of predicted electromagnetic field data and a Jacobian matrix of predicted spatial gradient data.

[0109] It should be noted that, is a Jacobian matrix of predicted electromagnetic field data corresponding to the kth iteration resistivity model; is a Jacobian matrix of predicted spatial gradient data corresponding to the kth iteration resistivity model.

[0110] S4.2: A coefficient matrix of the inversion equation is constructed according to the Jacobian matrix of predicted electromagnetic field data, the Jacobian matrix of predicted spatial gradient data and the gradient sensitivity scale matrix;

[0111] ;

[0112] ;

[0113] wherein, represents a left end item of the kth iteration inversion equation, i.e. the coefficient matrix; is a regularization parameter of the kth iteration; is a resistivity model small constraint weighting matrix; represents a Laplace operator for controlling the smoothness of the resistivity model; and are scalar weight coefficients for respectively controlling the balance of small constraints and smooth constraints.

[0114] It should be noted that in the process of constructing the coefficient matrix of the inversion equation, the data fitting, the gradient data and the prior information of the resistivity model are comprehensively considered to determine the direction and step length of the resistivity model update, and the stability and rationality of the inversion are ensured.

[0115] S4.3: A gradient direction of the inversion equation is constructed according to the Jacobian matrix of predicted electromagnetic field data, the Jacobian matrix of predicted spatial gradient data, the data residual term and the regularization term;

[0116] ;

[0117] ;

[0118] wherein, is a right end item of the kth iteration inversion equation, i.e. the gradient direction; is a current resistivity model; to reference resistivity model;

[0119] It should be noted that in the process of constructing the gradient direction of the inversion equation, the right end term is the negative gradient of the objective function with respect to the resistivity model parameters, which determines the direction of resistivity model update in the iterative inversion process, the data residual part reflects the influence of the difference between the observed data and the predicted data on the resistivity model update direction, and the regularization term adjusts the resistivity model update direction with the prior information of the resistivity model as a constraint. Together, the resistivity model update is closer to the observed data, guiding the inversion process to a direction closer to the true resistivity distribution.

[0120] S4.4: Combine the coefficient matrix with the gradient direction to establish the inversion equation.

[0121] ;

[0122] wherein, represents the left end term of the kth iteration inversion equation, that is, the coefficient matrix; is the difference between the resistivity model calculated in the kth step and the resistivity model obtained in the k-1th step; is the right end term of the kth iteration inversion equation, that is, the gradient direction;

[0123] S4.5: Calculate the resistivity model update step value of the current iteration, and correct the resistivity model parameters based on the resistivity model update step value;

[0124] ;

[0125] It should be noted that, is the difference between the resistivity model calculated in the kth step and the resistivity model obtained in the k-1th step, that is, the resistivity model update step value, which is used to correct the current resistivity model to obtain the corrected resistivity model parameters.

[0126] S4.6: Output the three-dimensional inversion result reflecting the spatial distribution of resistivity according to the corrected resistivity model parameters.

[0127] It should be noted that the corrected resistivity model parameters obtained by the last iteration are taken as the three-dimensional inversion result reflecting the spatial distribution of resistivity.

[0128] It should be noted that with reference to Figure 2 , Figure 2 is a schematic diagram of the sensitivity of electromagnetic gradient data to boundary changes. Figure 2 (a) represents a true three-dimensional model of electromagnetic gradient field containing a conductive anomaly; Figure 2 (b) represents a three-dimensional inversion result reflecting the spatial distribution of resistivity obtained by the method of the present application. Figure 2(a) the three-dimensional model in which the model is slightly displaced and changed in size; Figure 2 (c) the relative difference between the two curves in (b) is subtracted from one of them, then divided by the other, and then the logarithm of 10 is taken; Figure 2 (a) and Figure 2 (b) the relative difference between the electromagnetic field data responses of the two models; Figure 2 (d) the relative difference between the two curves in (c) is subtracted from one of them, then divided by the other, and then the logarithm of 10 is taken; Figure 2 (a) and Figure 2 (b) the relative difference between the spatial gradient data responses of the two models; Figure 2 (e) reflects the relative difference between the two models in electromagnetic field and gradient field responses, and the ordinate is the logarithmic value of the relative difference (the pink curve is Figure 2 (c) the relative difference between the two curves is subtracted from one of them, then divided by the other, and then the logarithm of 10 is taken; the blue curve is Figure 2 (d) the relative difference between the two curves is subtracted from one of them, then divided by the other, and then the logarithm of 10 is taken); the x-coordinate reflects the coordinate value of the measurement point on the measurement line relative to the midpoint of the transmitting source, and the coordinate origin X=0 is the midpoint of the transmitting source.

[0129] By comparing Figure 2 (a) and Figure 2 (b) the difference between the electromagnetic field data responses and the spatial gradient data responses of the two models, it is shown that the spatial gradient data response has a sensitivity far exceeding the original electromagnetic field to the change of the anomaly body boundary, and a slight boundary displacement can cause a relative change of more than one order of magnitude in the spatial gradient data response, so that the target body profile can be more clearly delineated.

[0130] The embodiment also provides a computer device suitable for the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving, which comprises a memory and a processor; the memory is used to store computer executable instructions, and the processor is used to execute the computer executable instructions to realize the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving as proposed in the above embodiment.

[0131] The computer device can be a terminal, and the computer device includes a processor, a memory, a communication interface, a display screen and an input device connected by a system bus. The processor of the computer device is used to provide computing and control capabilities. The memory of the computer device includes a non-volatile storage medium and an internal memory. The non-volatile storage medium stores an operating system and a computer program. The internal memory provides an environment for the operating system and the computer program in the non-volatile storage medium to run. The communication interface of the computer device is used to communicate with external terminals in a wired or wireless manner. The wireless manner can be achieved by WIFI, an operator network, NFC (near field communication) or other technologies. The display screen of the computer device can be a liquid crystal display screen or an electronic ink display screen. The input device of the computer device can be a touch layer overlaid on the display screen, or a key, trackball or touchpad arranged on the shell of the computer device, or an external keyboard, touchpad or mouse, etc.

[0132] The embodiment also provides a storage medium having a computer program stored thereon, the program being executed by a processor to implement the adaptive weighting high-precision joint inversion method based on electromagnetic gradient data driving proposed in the above embodiment; and the storage medium can be implemented by any type of volatile or non-volatile storage device or a combination thereof, such as a static random access memory (SRAM), an electrically erasable programmable read-only memory (EEPROM), an erasable programmable read-only memory (EPROM), a programmable read-only memory (PROM), a read-only memory (ROM), a magnetic memory, a flash memory, a magnetic disk or an optical disk.

[0133] In summary, the present application is achieved by: extending the classical Gauss-Newton method, using the residual term of the electromagnetic field data and the spatial gradient for inversion optimization, and introducing the gradient sensitivity scale matrix in the residual weighting, so that only the model parameters related to the abnormal area can be updated in each inversion iteration, and the mapping to the three-dimensional finite element domain for forward calculation. This mechanism improves the boundary resolution while significantly reducing the computational cost and can completely capture the three-dimensional physical response, overcoming the limitations of the traditional constant constraint gradient residual term method in terms of abnormal area boundary identification accuracy and background area stability due to the neglect of grid non-uniformity and physical property spatial heterogeneity, and providing a practical solution for electromagnetic detection in complex geological environments.

[0134] It should be noted that the above examples are only used to illustrate the technical solutions of the present application but not to limit the present application. Although the present application has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present application can be modified or equivalently replaced without departing from the spirit and scope of the present application, and all modifications and equivalents should be included in the scope of the claims of the present application.

Claims

1. A high-precision joint inversion method based on adaptive weighting driven by electromagnetic gradient data, characterized in that: include, Spatial gradient data of the original electromagnetic field observation data were calculated using the Gaussian weighted neighborhood method. The spatial gradient difference between adjacent frequencies is calculated to form a frequency difference spatial gradient. The frequency difference spatial gradients of all measurement points are normalized and combined to obtain the overall gradient intensity normalization matrix. The overall gradient intensity normalization matrix is ​​then processed using a contrast weighting function with an adjustable exponent to construct a gradient sensitive scale matrix. By combining the regular data residual term, the data residual term, and the regularization term, a joint objective function is constructed. The gradient-sensitive scale matrix is ​​embedded into the data residual term of the joint objective function to establish the inversion equation. The inversion equation is iteratively calculated and the resistivity vector of the resistivity model at each grid in the subsurface is updated to obtain the three-dimensional inversion result that reflects the spatial distribution of resistivity. The specific steps for constructing the gradient-sensitive scaling matrix are as follows: The overall gradient intensity normalization matrix is ​​used as the input parameter of the contrast weighting function; Each element in the overall gradient intensity normalization matrix is ​​processed using a contrast weighting function with an adjustable exponent; The processed elements are arranged into a diagonal matrix to form a gradient-sensitive scaling matrix. This matrix is ​​then used to adjust the contribution weights of individual data points in the residual terms.

2. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 1, characterized in that: The specific steps for calculating the spatial gradient data of the original electromagnetic field observation data using the Gaussian weighted neighborhood method are as follows. Acquire electromagnetic field data at the measuring point and electromagnetic field data at neighboring measuring points; Based on the electromagnetic field data under the frequency conditions at the measuring point and the spatial distance between the electromagnetic field data under the same frequency conditions at neighboring measuring points, the gradient components of each pair of neighboring points are calculated. The distance between neighboring measurement points is calculated using an exponential function and normalized to obtain normalized Gaussian weights. These normalized Gaussian weights are then used to weight the gradient components of each pair of neighboring points to generate spatial gradient data.

3. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 1, characterized in that: The specific steps for calculating the spatial gradient difference between adjacent frequencies to form a frequency difference spatial gradient are as follows. Obtain the spatial gradient between the frequency at the measurement point and the previous frequency at the measurement point, and select a frequency function that matches the frequency sampling scheme. Calculate the difference between the frequency at the measuring point and the spatial gradient at the previous frequency at the measuring point; Calculate the difference between the frequency at the measuring point and the previous frequency at the measuring point; The ratio of the difference in spatial gradients to the difference in frequency function values ​​is taken as the frequency difference spatial gradient.

4. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 1, characterized in that: The normalized combination of the frequency difference spatial gradients of all measurement points yields the normalized gradient intensity matrix. The specific steps are as follows: Collect the frequency differential spatial gradient at each measuring point and determine the maximum value of the frequency differential spatial gradient at all measuring points; The frequency difference space gradient of each measuring point is normalized based on the maximum value to obtain the normalized frequency difference space gradient vector of each measuring point. Combine the normalized frequency difference space gradient vectors of all measurement points to form the overall gradient intensity normalized matrix.

5. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 1, characterized in that: The steps for constructing a joint objective function by combining the conventional data residual term, the data residual term, and the regularization term are as follows: Based on the difference between the observed electromagnetic field data and the electromagnetic field data predicted by forward modeling, a conventional data residual term is constructed. A joint objective function is constructed by weighting and combining the regular data residuals, the data residuals, and the regularization terms.

6. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 5, characterized in that: The data residuals include, The processing results of the differences between the weighted matrix of diagonal electromagnetic field data and the observed and predicted electromagnetic field data; The result of processing the difference between the weighted matrix of the diagonal gradient field data and the observed gradient data and the predicted gradient data; The result of adjusting the influence of the gradient-sensitive scaling matrix on the degree of influence of individual data points in the gradient data residual.

7. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 1, characterized in that: The inversion equation includes the Jacobian matrix for predicting electromagnetic field data and the Jacobian matrix for predicting spatial gradient data.

8. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 1, characterized in that: The specific steps for embedding the gradient-sensitive scaling matrix into the data residual term of the joint objective function to establish the inversion equation are as follows: The coefficient matrix of the inversion equation is constructed based on the Jacobian matrix of the predicted electromagnetic field data, the Jacobian matrix of the predicted spatial gradient data, and the gradient-sensitive scale matrix. Based on the Jacobian matrix of the predicted electromagnetic field data, the Jacobian matrix of the predicted spatial gradient data, the data residual terms, and the regularization terms, the gradient direction of the inversion equation is constructed. By combining the coefficient matrix with the gradient direction, an inversion equation is established.

9. The adaptive weighted high-precision joint inversion method based on electromagnetic gradient data-driven method as described in claim 1, characterized in that: The iterative calculation of the inversion equation and the updating of the resistivity model parameters yield a three-dimensional inversion result reflecting the spatial distribution of resistivity. The specific steps are as follows: Calculate the resistivity model update step value for the current iteration, and correct the resistivity model parameters based on the resistivity model update step value; Based on the corrected resistivity model parameters, the output is a three-dimensional inversion result reflecting the spatial distribution of resistivity.

Images