Time shift induced polarization method four-dimensional inversion method based on finite memory quasi-Newton method

Abstract

The invention provides a time shift induced polarization method four-dimensional inversion method based on a finite memory quasi-Newton method. The time shift induced polarization method four-dimensional inversion method comprises the following steps: 1) preprocessing observation data; 2) reading a file required by inversion; 3) calculating by adopting a master-slave mode designed by an MPI parallel program; (4) the forward modeling response of the time-shifting polarizability model is calculated in parallel through MPI; 5) calculating an objective function and a gradient thereof in parallel through MPI; 6) performing approximate calculation by using a finite memory quasi-Newton method to obtain an inverse matrix of the Hessian matrix; 7) calculating the updating direction of the time-shifting polarizability model; 8) searching a proper step length alpha in the current iteration; 9) updating the model; and 10) judging whether an iteration termination condition is met or not, if so, outputting a result, otherwise, returning to the step 3) to carry out circulation until an inversion result meets the iteration termination condition. According to the method, the calculation of a Hessian matrix is avoided, the inversion calculation amount and the required memory are reduced, and the problem that the time required by inversion is too long is solved by adopting a multi-thread synchronous calculation mode.

Description

Technical Field

[0001] The invention relates to the technical field of geophysical exploration, and in particular to a four-dimensional inversion method of a time-shifted induced polarization method based on a finite-memory quasi-Newton method. Background Art

[0002] In actual work, the phenomenon of additional electric fields that vary with time will be generated when charging and discharging underground is called the induced polarization effect (IP effect for short). The induced polarization method detects information such as the polarizability of underground media by observing the IP effect of rocks and ores. The induced polarization method has the advantages of low cost, high acquisition efficiency, and convenient work. It has great advantages in solving problems that require long-term monitoring, such as environmental protection, ecological restoration, seawater intrusion, and geological disasters. In addition to the three-dimensional spatial dimension, the monitoring data also includes the time dimension. In order to process the monitoring data of the induced polarization method, foreign scholars have developed a four-dimensional inversion method of the time-shifted induced polarization method. By simultaneously inverting the IP data of the same area collected at different times and using time-shift terms for constraints, it can monitor small changes in underground structures.

[0003] Since the scale of monitoring data is several times that of a single moment, there are three problems in the operation process: the inversion calculation occupies too much memory and is difficult to run on a single computer; the inversion calculation efficiency is too low and takes a lot of time; the observation data is interfered by external noise, and the degree of noise interference at different moments is different. The polarizability results obtained by inversion are different from the real model to varying degrees, making it difficult for the monitoring data inversion results to maintain continuity over time and unable to reflect the slight changes in the underground medium. Therefore, it is of great significance to study methods to reduce the computational complexity and memory requirements of the four-dimensional inversion of the time-shifted induced polarization method, improve its computational efficiency, and suppress external noise interference to promote the development and application of the four-dimensional inversion of the time-shifted induced polarization method. Summary of the invention

[0004] In view of the technical problems existing in the above-mentioned background technology, the present invention proposes a four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method. The method has a reasonable conception and utilizes limited memory quasi-Newton method. The inverse matrix of the Hessian matrix is ​​approximated by the gradient and model of the objective function in recent times, thereby avoiding the calculation of the Hessian matrix and reducing the amount of calculation and the required memory for inversion. At the same time, MPI parallel programming is introduced in the inversion calculation process, and a multi-threaded synchronous calculation method is adopted to solve the problem of too long time required for inversion and improve the inversion efficiency. In addition, the introduction of time-shift terms improves the anti-noise ability of the inversion, improves the continuity of the polarizability inversion results at different times over time, and can more clearly display the slight changes in the underground medium.

[0005] In order to solve the above technical problems, the present invention provides a four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method, which mainly includes the following steps:

[0006] 1) Preprocess the observation data;

[0007] 2) Read the data, parameters and model files required for inversion;

[0008] 3) The master-slave mode of MPI parallel programming is used for computing. In the master-slave mode, all processes are divided into master processes and slave processes. The master process is responsible for allocating tasks, transmitting messages and collecting data. The assigned tasks are completed by the master and slave processes together.

[0009] 4) The forward response of the time-shifted polarizability model is calculated by MPI parallel computing. The forward response function is F(M η ), M η is a time-shifted polarizability model composed of underground polarizability models at multiple moments. η is the initial model, and the subsequent ones are the results of the previous iteration in the iterative process. After the forward response is obtained, the difference between the forward response and the observed apparent polarizability data is calculated;

[0010] 5) Construct the objective function of inversion. The objective function consists of data terms, model terms and time-shift terms. Inversion is to obtain a time-shift polarizability model that is closer to the actual underground situation by continuously minimizing the objective function. After the objective function is constructed, the objective function and its gradient are calculated in parallel using MPI;

[0011] 6) Use the limited memory quasi-Newton method to approximate the inverse matrix of the Hessian matrix;

[0012] 7) Calculate the update direction p of the time-shifted polarizability model k ;

[0013] 8) Searching for the step size of updating the time-shifted polarizability model and finding a suitable step size α in the current iteration;

[0014] 9) Update based on the current time-shifted polarizability model, using the step size α and update direction p k Multiply them and add them to the current time-shift polarizability model to get the new time-shift polarizability model after iteration.

[0015] 10) Determine whether the iteration termination condition is met, and if so, output the result; otherwise, return to step 3) and loop until the inversion result meets the iteration termination condition.

[0016] The four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method, wherein the apparent polarizability is obtained by forward modeling the response function in step 4), and the specific process is as follows:

[0017] First, according to the current continuity equation, the current density j at any position satisfies:

[0018]

[0019] Wherein, in the above formula (1), σ is the conductivity of the underground medium, u is the potential; r0 = (x0, y0, z0), is the position of the point source; r is the distance from the point source to the measurement point, I is the current intensity, δ is the Dirac function, and j is the current density;

[0020] Then, the finite difference method is used to discretize the above equation (1) to obtain a linear equation system. Solving the linear equation system can obtain the potential value u of the observation point. Assuming that the underground conductivity model is known, the surface potential value is u, and the input time-shift polarization rate model M η , time-shifted polarizability model M η is the set of underground polarizability models at multiple moments, and m η,i is the underground polarizability model at the i-th moment, then the time-shifted polarizability model M η It can be defined as:

[0021] M η =[m η,1 m η,2 ... m η,nt ] T (2);

[0022] In the above formula (2), nt represents the number of moments;

[0023] Because the polarizability of the underground medium will induce the induced polarization effect, the potential value of the ground surface changes from u to u η , then the apparent polarizability η a for:

[0024]

[0025] The four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method, wherein: the objective function in step 5) is composed of the data item Φ d (M η ), model term Φ m (M η ) and the time-shift term Φ t (M η ) and is obtained by the following formula (4):

[0026]

[0027] In the above formula (4), and are the weight factors of the model term and the time-shift constraint term respectively. By adjusting the size of the weight factor, the strength of the constraint effect in the inversion process can be controlled;

[0028] The data item Φ d (Mη ) is obtained by the following formula (5):

[0029]

[0030] In the above formula (5), D η,obs is the apparent polarizability data observed at all times, M η is the time-shifted polarizability model, F(M η ) is the forward response function, is the data covariance matrix;

[0031] The model term Φ m (M η ) is obtained by the following formula (3):

[0032]

[0033] In the above formula (6), C M,η is the model covariance matrix, M η,0 is the initial model;

[0034] The time shift term Φ t (M η ) is obtained by the following formula (7):

[0035]

[0036] In the above formula (7), m η,k is the underground polarizability model at the kth moment;

[0037] The gradient of the objective function is obtained by the following formula (8):

[0038]

[0039] In the above formula (8), J η is the Jacobian matrix, which is the forward response function F(M η ) η The partial derivative of .

[0040] The four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method, wherein the inverse matrix of the Hessian matrix is ​​approximately obtained by using the finite memory quasi-Newton method in step 6), and the calculation method is shown in the following formula (9):

[0041]

[0042] In the above formula (9), the matrix B k+1 It can be approximately considered as the inverse matrix of the Hessian matrix, where each element in formula (9) can be calculated by formula (10):

[0043]

[0044] In the above formula (10), g k is the gradient of the objective function at the kth iteration in the iterative process, M k is the time-shifted polarizability model obtained in the kth iteration. The remaining parameters can be calculated through the gradient of the objective function of the last two times and the time-shifted polarizability model.

[0045] The four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method, wherein the matrix B k+1 It is calculated using the results of the last 10 times, using matrix B k+1 And the gradient of the objective function, we can get the update direction p of the time-shifted polarizability model k , as shown in the following formula (11):

[0046] p k =-B k+1 g k (11).

[0047] The four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method, wherein the first step size α1 in step 8) is the inverse of the second-order norm of the objective function gradient, and the initial step sizes of the remaining iterations are all 1, as shown in the following formula (12):

[0048]

[0049] The four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method, wherein: the information in the observation data in step 1) includes the coordinates, apparent polarizability and current of the transmitting electrode and the receiving electrode corresponding to all measuring points.

[0050] The four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method, wherein: the preprocessing process in the step 1) includes data conversion, elimination of distortion points and adjustment of data format; the data conversion is mainly to obtain the apparent polarizability through the observed electric field calculation; the distortion points to be eliminated are observed observation data, and the distortion points whose apparent polarizability values ​​are obviously higher or lower than those of adjacent data points are eliminated; the data format adjustment includes adjusting the data to a fixed format.

[0051] The four-dimensional inversion method of the time-shift induced polarization method based on the limited memory quasi-Newton method, wherein: the files required for inversion in the step 2) include processed observation data, parameter data and initial model; the parameter data are mainly weight factors of model terms and time-shift terms, and the number of times used in the limited memory quasi-Newton method; the initial model includes a time-shift resistivity model and a time-shift polarizability model; the time-shift resistivity model is a known underground resistivity model, and the time-shift polarizability model is set to a uniform half-space with a polarizability value of 0.01.

[0052] By adopting the above technical solution, the present invention has the following beneficial effects:

[0053] The four-dimensional inversion method of the time-shift induced polarization method based on the limited memory quasi-Newton method of the present invention is reasonably conceived. Compared with the traditional method of inverting the data at different times in the monitoring data separately, the present invention synchronously inverts the apparent polarizability data at multiple times, and can simultaneously obtain the polarizability inversion results at multiple times, which belongs to the four-dimensional inversion method. The limited memory quasi-Newton method is used to approximate the inverse matrix of the Hessian matrix through the gradient and model of the recent objective function, which avoids the calculation of the Hessian matrix, reduces the calculation amount of inversion and the required memory, and the required memory is only 10% of the original. In the inversion calculation process, MPI parallel programming is introduced, and the method of multi-threaded synchronous calculation is adopted to solve the problem of too long time required for inversion, and improve the calculation efficiency. In the embodiment, the running time of this program is 7.5% of the original. In the inversion process, by adding the time-shift term in the objective function, the inversion process is constrained, the influence of external noise interference is suppressed, and it is ensured that the inversion results at different times change continuously with time, which can better reflect the real changes in the underground. BRIEF DESCRIPTION OF THE DRAWINGS

[0054] In order to more clearly illustrate the specific implementation methods of the present invention or the technical solutions in the prior art, the drawings required for use in the specific implementation methods or the description of the prior art will be briefly introduced below. Obviously, the drawings in the following description are some implementation methods of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0055] Figure 1 It is a flow chart of the four-dimensional inversion method of the time-shifted induced polarization method based on the limited memory quasi-Newton method according to the present invention;

[0056] Figure 2 This is a schematic diagram of the model constructed in the four-dimensional inversion method of the time-shift induced polarization method based on the finite memory quasi-Newton method of the present invention, which includes the models at three times T1, T2, and T3 from left to right; (Note: the resistivity of A is 100Ω·m, and the polarizability is 0.01; the resistivity of B is 10Ω·m, and the polarizability is 0.2);

[0057] Figure 3 It is a schematic diagram of the observation system involved in the four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method of the present invention;

[0058] Figure 4 A comparison chart of the polarizability inversion results of the traditional method (left) and the four-dimensional inversion method of the time-shifted induced polarization method based on the finite memory quasi-Newton method of the present invention (right) at a depth of 8m. DETAILED DESCRIPTION

[0059] The technical solution of the present invention will be described clearly and completely below in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0060] The present invention is further explained below in conjunction with specific implementation modes.

[0061] like Figure 1 As shown, this embodiment provides a four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method, which specifically includes the following steps:

[0062] S100, preprocessing the observation data, the observation data includes information such as the coordinates of the transmitting electrode and the receiving electrode corresponding to all the measuring points, the current and the apparent polarizability. The preprocessing process includes data conversion, elimination of distortion points, adjustment of data format, etc. The data conversion is mainly to obtain the apparent polarizability through the observed electric field calculation; elimination of distortion points is to observe the observation data, and the distortion points whose apparent polarizability values ​​are obviously higher or lower than the adjacent data points are eliminated; adjustment of data format includes adjusting the data to a fixed format, which is from front to back the coordinates of the transmitting electrode, the coordinates of the receiving electrode, the current, the apparent polarizability, etc.

[0063] In this embodiment, the forward response of the model is used as the observation data. Figure 2 As shown in Figure 1, the model contains data at three times: T1, T2, and T3. The polarizability of the surrounding rock A is 0.01, and the polarizability of the abnormal body B is 0.2. Figure 3 This is a schematic diagram of the observation system. Data is collected by arranging transmitting and receiving electrodes on the surface. The electrodes are arranged in a rectangular array. The triangles are the positions of the transmitting electrodes, the dots are the positions of the receiving electrodes, and the solid box is the projection of the anomaly B on the surface. When each transmitting electrode is powered, all adjacent receiving electrodes in the X direction receive simultaneously to ensure that enough information is obtained in three-dimensional space. After forward modeling, 5%, 20%, and 5% Gaussian random errors are added to the observation data at times T1, T2, and T3, respectively, to simulate the noise interference in the actual collection process. Here, a higher error level is added to the observation data at time T2 to simulate the situation where the external noise interference is large, and the data after the error is added is used as the observation data.

[0064] S200, read the files required for inversion, including the processed observation data, parameter data, initial model, etc. The parameter data is mainly the weight factors of the model term and the time-shift term. The inverse matrix of the Hessian matrix is ​​calculated in the limited memory quasi-Newton method using the objective function gradient and the number of models (in this embodiment, the objective function gradient and model for nearly 10 times). The initial model includes a time-shift resistivity model and a time-shift polarizability model. The time-shift resistivity model is generally a known underground resistivity model, and the time-shift polarizability model is set to a uniform half-space with a polarizability value of 0.01.

[0065] S300, the master process collects data, transmits messages and assigns tasks. The master process is responsible for assigning tasks, transmitting messages and collecting data. The assigned tasks are completed by the master and slave processes. In this embodiment, the number of master and slave processes is 20.

[0066] S400, calculating the forward response of the current time-shifted polarizability model in parallel through MPI, the current model is the initial model in the first iteration, and the subsequent model is the result of the previous iteration in the iteration process; the forward response function is F(M η ), M η The time-shift polarizability model is composed of multiple underground polarizability models at different moments. η ) can calculate the apparent polarizability; the implementation process is as follows:

[0067] First, according to the current continuity equation, the current density j at any position satisfies:

[0068]

[0069] Wherein, in the above formula (1), σ is the conductivity of the underground medium, u is the potential, r0 = (x0, y0, z0), is the position of the point source, r is the distance from the point source to the measurement point, I is the current intensity, δ is the Dirac function, and j is the current density;

[0070] Then, the finite difference method is used to discretize the above equation (1) to obtain a large linear equation system. The potential value u of the observation point can be obtained by solving it. Assuming that the underground conductivity model is known, the surface potential value is u, and the input time-shift polarization rate model M η , time-shifted polarizability model M η is a collection of underground polarizability models at multiple moments; let m η,i is the underground polarizability model at a certain moment, then the time-shifted polarizability model M η It can be defined as:

[0071] M η =[m η,1 m η,2 ... m η,nt ] T (2);

[0072] In the above formula (2), nt represents the number of moments;

[0073] Due to the induced polarization effect caused by the polarizability of the underground medium, the potential of the surface changes from u to u η , then the apparent polarizability is:

[0074]

[0075] The MPI multi-threaded parallel computing is used to calculate the difference between the forward response of the current time-shifted polarizability model and the observed apparent polarizability data.

[0076] S500, construct the inversion objective function, which is composed of data term, model term and time shift term. By continuously minimizing the objective function, a time shift polarizability model closer to the actual underground situation is obtained. The objective function and its gradient are calculated in parallel through MPI. The objective function is composed of data term Φ d (M η ), model term Φ m (M η ) and the time-shift term Φ t (M η ), as shown in formula (4):

[0077]

[0078] In the above formula (4), and are the weight factors of the model term and the time-shift constraint term respectively. By adjusting the size of the weight factor, the strength of the constraint effect in the inversion process can be controlled;

[0079] The above data item Φ d (M η ) is obtained by the following formula (5):

[0080]

[0081] In the above formula (5), D η,obs is the apparent polarizability data observed at all times; M η is the time-shifted polarizability model; F(M η ) is the forward response function, is the data covariance matrix.

[0082] The above model term Φ m (M η ) is obtained by the following formula (6):

[0083]

[0084] In the above formula (6), C M,η is the model covariance matrix, Mη,0 is the initial model.

[0085] The time shift term Φ mentioned above t (M η ) is obtained by the following formula (7):

[0086]

[0087] In the above formula (7), m η,k is the underground polarizability model at time k, m η,k+1 It is the polarizability model of the underground at time k+1.

[0088] The gradient of the above objective function is obtained by the following formula (8):

[0089]

[0090] In the above formula (8), J η is the Jacobian matrix, which is the forward response function F(M η ) η The partial derivative, Jacobian matrix J η The T in the upper right corner represents the transpose of the matrix.

[0091] S600. Using the limited memory quasi-Newton method, the inverse matrix of the Hessian matrix is ​​approximately obtained. The calculation method is shown in the following formula (9):

[0092]

[0093] In the above formula (9), the matrix B k+1 It can be approximately considered as the inverse matrix of the Hessian matrix, where each element in formula (9) can be calculated by formula (10):

[0094]

[0095] In the above formula (10), g k is the gradient of the objective function at the kth iteration in the iterative process, m k is the time-shifted polarizability model obtained in the kth iteration. The remaining parameters can be obtained by the gradient of the last two objective functions and the time-shifted polarizability model M k In this embodiment, the matrix B is obtained by calculating the results of nearly 10 times. k+1 .

[0096] S700, using matrix B k+1 And the gradient of the objective function, the update direction p of the time-shifted polarizability model can be calculated k , as shown in the following formula (11):

[0097] p k =-B k+1g k (11);

[0098] S800, search for the step size of the time-shifted polarizability model update, and find the appropriate step size α in the current iteration. In the present invention, the first step size α1 is the inverse of the second-order norm of the objective function gradient, and the initial step sizes of the remaining iterations are all 1, as shown in the following formula (12); this is related to the characteristics of the limited memory quasi-Newton method itself, because once the approximate Hessian matrix can fully approximate the real Hessian matrix, the step size 1 can meet the requirements:

[0099]

[0100] S900, updating the current time-shift polarization rate model to obtain a new time-shift polarization rate model after iteration, the calculation formula is:

[0101] M k+1 =M k +αp k (13).

[0102] S1000, calculating the RMS value of the time-shifted polarizability model obtained after being updated by formula (13), and the calculation method is as formula (14); judging whether the iteration termination condition is met, if it is met, outputting the time-shifted polarizability model of the final iteration as the inversion result, otherwise returning to step S300 until the termination condition of the inversion iteration is met; the termination condition of the inversion iteration is that the RMS of the data fitting error calculated at all times is less than the expected value or the number of iterations reaches the set maximum number of iterations;

[0103]

[0104] The invention has a reasonable conception, utilizes the limited memory quasi-Newton method, approximates the inverse matrix of the Hessian matrix through the gradient and model of the recent objective function, avoids the calculation of the Hessian matrix, reduces the calculation amount of inversion and the required memory; at the same time, introduces MPI parallel programming in the inversion calculation process, adopts a multi-threaded synchronous calculation method, solves the problem of too long time required for inversion, and improves the inversion efficiency; in addition, in the inversion objective function, adds a time shift term to constrain the inversion, so that the polarizability inversion results at different times maintain continuity in time, and can better reflect the slight changes of the underground medium over time.

[0105] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or replace some or all of the technical features therein with equivalents. However, these modifications or replacements do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A four-dimensional inversion method of time-shifted induced polarization method based on finite memory quasi-Newton method, characterized in that: The main steps include: 1) Preprocess the observation data; 2) Read the data, parameters and model files required for inversion; 3) The master-slave mode of MPI parallel programming is used for computing. In the master-slave mode, all processes are divided into master processes and slave processes. The master process is responsible for allocating tasks, transmitting messages and collecting data. The assigned tasks are completed by the master and slave processes together. 4) The forward response of the time-shifted polarizability model is calculated by MPI parallel computing. The forward response function is F(M η ), M η The time-shift polarizability model is composed of underground polarizability models at multiple moments. η is the initial model, and the subsequent ones are the results of the previous iteration in the iterative process. After the forward response is obtained, the difference between the forward response and the observed apparent polarizability data is calculated; 5) Construct the objective function of inversion. The objective function consists of data terms, model terms and time-shift terms. Inversion is to obtain a time-shift polarizability model that is closer to the actual underground situation by continuously minimizing the objective function. After the objective function is constructed, the objective function and its gradient are calculated in parallel using MPI; 6) Use the limited memory quasi-Newton method to approximate the inverse matrix of the Hessian matrix; 7) Calculate the update direction p of the time-shifted polarizability model k ; 8) Searching for the step size of updating the time-shifted polarizability model and finding a suitable step size α in the current iteration; 9) Update based on the current time-shifted polarizability model, using the step size α and update direction p k Multiply them and add them to the current time-shift polarizability model to get the new time-shift polarizability model after iteration. 10) Determine whether the iteration termination condition is met, and if so, output the result; otherwise, return to step 3) and loop until the inversion result meets the iteration termination condition.

2. The four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method according to claim 1, characterized in that: In step 4), the apparent polarizability is obtained by forward modeling the response function, and the specific process is as follows: First, according to the current continuity equation, the current density j at any position satisfies: Wherein, in the above formula (1), σ is the conductivity of the underground medium, u is the potential; r0 = (x0, y0, z0), is the position of the point source; r is the distance from the point source to the measurement point, I is the current intensity, δ is the Dirac function, and j is the current density; Then, the finite difference method is used to discretize the above equation (1) to obtain a linear equation system. Solving the linear equation system can obtain the potential value u of the observation point. Assuming that the underground conductivity model is known, the surface potential value is u, and the input time-shift polarization rate model M η , time-shifted polarizability model M η is the set of underground polarizability models at multiple moments, and m η,i is the underground polarizability model at the i-th moment, then the time-shifted polarizability model M η It can be defined as: M η =[m η,1 m η,2 ...m η,nt ] T (2); In the above formula (2), nt represents the number of moments; Because the polarizability of the underground medium will induce the induced polarization effect, the potential value of the ground surface changes from u to u η , then the apparent polarizability η a for:

3. The four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method according to claim 1, characterized in that: The objective function in step 5) is composed of the data term Φ d (M η ), model term Φ m (M η ) and the time-shift term Φ t (M η ) and is obtained by the following formula (4): In the above formula (4), and are the weight factors of the model term and the time-shift constraint term respectively. By adjusting the size of the weight factor, the strength of the constraint effect in the inversion process can be controlled; The data item Φ d (M η ) is obtained by the following formula (5): In the above formula (5), D η,obs is the apparent polarizability data observed at all times, M η is the time-shifted polarizability model, F(M η ) is the forward response function, is the data covariance matrix; The model term Φ m (M η ) is obtained by the following formula (3): In the above formula (6), C M,η is the model covariance matrix, M η,0 is the initial model; The time shift term Φ t (M η ) is obtained by the following formula (7): In the above formula (7), m η,k is the underground polarizability model at the kth moment; The gradient of the objective function is obtained by the following formula (8): In the above formula (8), J η is the Jacobian matrix, which is the forward response function F(M η ) η The partial derivative of .

4. The four-dimensional inversion method of time-shifted induced polarization method based on finite memory quasi-Newton method according to claim 1, characterized in that: In step 6), the inverse matrix of the Hessian matrix is ​​approximately obtained by using the limited memory quasi-Newton method, and the calculation method is shown in the following formula (9): In the above formula (9), the matrix B k+1 It can be approximately considered as the inverse matrix of the Hessian matrix, where each element in formula (9) can be calculated by formula (10): In the above formula (10), g k is the gradient of the objective function at the kth iteration in the iterative process, M k is the time-shifted polarizability model obtained in the kth iteration. The remaining parameters can be calculated through the gradient of the objective function of the last two times and the time-shifted polarizability model.

5. The four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method according to claim 4, characterized in that: The matrix B k+1 It is calculated using the results of the last 10 times, using matrix B k+1 And the gradient of the objective function, we can get the update direction p of the time-shifted polarizability model k , as shown in the following formula (11): p k =-B k+1 g k (11)。 6. The four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method according to claim 1, characterized in that: The first step size α1 in step 8) is the inverse of the second-order norm of the objective function gradient, and the initial step sizes of the remaining iterations are all 1, as shown in the following formula (12):

7. The four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method according to claim 1, characterized in that: The information in the observation data in step 1) includes the coordinates, apparent polarizability and current of the transmitting electrode and the receiving electrode corresponding to all the measuring points.

8. The four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method according to claim 1, characterized in that: The preprocessing process in step 1) includes data conversion, elimination of distortion points and adjustment of data format; the data conversion is mainly to obtain the apparent polarizability through the observed electric field calculation; the elimination of distortion points is to observe the observed data, and the distortion points whose apparent polarizability values ​​are obviously higher or lower than the adjacent data points are eliminated; the adjustment of the data format includes adjusting the data to a fixed format.

9. The four-dimensional inversion method of time-shifted induced polarization method based on limited memory quasi-Newton method according to claim 1, characterized in that: The files required for inversion in the step 2) include processed observation data, parameter data and initial model; the parameter data are mainly weight factors of model terms and time-shift terms, and the number of times used in the finite memory quasi-Newton method; the initial model includes a time-shift resistivity model and a time-shift polarizability model; the time-shift resistivity model is a known underground resistivity model, and the time-shift polarizability model is set to a uniform half-space with a polarizability value of 0.01.

Images