Joint inversion method and device based on cross gradient term weight factor

Abstract

The invention provides a joint inversion method and device based on a cross gradient term weight factor, and the method comprises the steps: endowing each grid unit in an underground grid with an initial parameter value of a uniform half space, enabling the initial parameter value to serve as an initial model of joint inversion, enabling the underground grid to be obtained through the grid division of a research region, and enabling the initial model to serve as an initial model of joint inversion; the data range of the research area is determined according to the data range of the observation area; an objective function of joint inversion is constructed, the objective function comprises a cross gradient term and a first weight factor, the first weight factor comprises a coefficient matrix and a weight constant, and each matrix element in the coefficient matrix corresponds to each grid unit in the underground grid; the matrix element values corresponding to the grid units of the underground grid in different areas are different; and utilizing the apparent resistivity data and the apparent polarizability data to perform minimization calculation on the target function so as to update the initial model and obtain a resistivity model and a polarizability model. The precision of the inversion result can be improved.

Description

Technical Field

[0001] The present application relates to the technical field of geophysical exploration. Specifically, it relates to a joint inversion method and device based on the weight factor of the cross-gradient term. Background Technique

[0002] Geophysical inversion is an important concept in geophysics. It refers to the process of using geophysical observation data (such as seismic wave records, gravity measurement data, magnetic measurement data, etc.) to infer the physical properties (such as velocity structure, density distribution, magnetic distribution, etc.) and geological structure of the Earth's interior. During the inversion process, there may be multiple different subsurface structures or models that can fit the same set of observation data, which is the non-uniqueness of geophysical inversion. The non-uniqueness of geophysical inversion stems from multiple aspects. Mathematically, both non-linearity and ill-posedness can cause non-uniqueness of the solution; in the physical model, the geological structure is complex and the properties are equivalent, and multiple combinations can produce similar responses; the observation data also has problems such as being finite, incomplete, and having insufficient precision and resolution. The interweaving of these factors leads to the fact that the inversion results often have multiple possibilities.

[0003] The cross-gradient joint inversion method is a current mainstream structural coupling joint inversion method. By multiplying the gradients of different physical properties, the boundaries of different physical properties are mutually constrained, which can improve the similarity of the inversion results of different physical properties in the spatial structure. It mainly realizes the joint inversion process by introducing a cross-gradient term into the objective function. By minimizing the objective function, the subsurface model can be obtained.

[0004] In the process of minimizing the objective function, the selection of the weight factor is particularly important. If it is too small, its constraint effect will be limited; if it is too large, it will play the opposite role and reduce the accuracy of the inversion result. However, there is currently no good method for selecting the weight factor of the cross-gradient term. Because the cross-gradient term does not have a value everywhere, and it only has a value when both models are changing and the change directions are inconsistent. Therefore, how to select the weight factor of the cross-gradient term to ensure that the cross-gradient term plays a good role in the inversion process is of great significance. Summary of the Invention

[0005] In view of this, the purpose of the present application is to provide a joint inversion method and device based on the weight factor of the cross-gradient term. By selecting an effective weight factor of the cross-gradient term, the non-uniqueness of the joint inversion is reduced, which helps to improve the accuracy of the inversion result.

[0006] In the first aspect, an embodiment of the present application provides a joint inversion method based on the weight factor of the cross-gradient term. The joint inversion method includes:

[0007] Obtain apparent resistivity data and apparent polarization rate data according to multiple survey lines arranged in advance in the observation area of the underground medium;

[0008] Assign the initial parameter values of the homogeneous half-space to each grid cell in the subsurface grid as the initial model for the joint inversion. The initial parameter values include an initial resistivity value and an initial polarizability value. The subsurface grid is obtained by meshing the study area, and the data range of the study area is determined according to the data range of the observation area.

[0009] Construct an objective function for the joint inversion. The objective function includes a cross-gradient term and a first weight factor corresponding to the cross-gradient term. The first weight factor includes a coefficient matrix and a weight constant. Each matrix element in the coefficient matrix corresponds to each grid cell in the subsurface grid, and the numerical values of the matrix elements corresponding to the grid cells in different regions of the subsurface grid are different. The weight constant is used to adjust the proportion of the cross-gradient term in the objective function.

[0010] Use the apparent resistivity data and the apparent polarizability data to perform a minimization calculation on the objective function to update the initial model and obtain a resistivity model and a polarizability model.

[0011] In an alternative embodiment, multiple survey lines are arranged in the observation area of the subsurface medium, and a transmitting electrode and a receiving electrode are provided on each survey line. The obtaining of the apparent resistivity data and the apparent polarizability data according to the multiple survey lines pre-arranged in the observation area of the subsurface medium includes:

[0012] When the transmitting electrode on each survey line is powered, the adjacent receiving electrodes on the corresponding survey line receive to obtain the apparent resistivity data and the apparent polarizability data.

[0013] In an alternative embodiment, the data range of the study area is determined by the following steps:

[0014] Obtain the data ranges in the X direction, Y direction, and Z direction of the observation area;

[0015] Set target multiples for the data ranges in the X direction, Y direction, and Z direction to obtain the data ranges of the study area in the X direction, Y direction, and Z direction; wherein, the target multiple is between 3 times and 5 times.

[0016] In an alternative embodiment, the objective function further includes a data term, a model term, and a second weight factor corresponding to the model term. The data term is used to calculate the difference between the apparent resistivity data and the resistivity data of the inversion response, and the difference between the apparent polarization rate data and the polarization rate data of the inversion response. The model term is used to calculate the smoothness of the resistivity model and the polarization rate model. The second weight factor characterizes the proportion of the model term in the objective function.

[0017] In an alternative embodiment, the objective function is expressed by the following formula:

[0018] Φ(m ρ ) = Φ d (m ρ ) + λ1Φ m (m ρ ) + γ1KΦ cg (m ρ , m η )

[0019] Φ(m η ) = Φ d (m η ) + λ2Φ m (m η ) + γ2KΦ cg (m η , m ρ )

[0020] Where m ρ represents the resistivity model of the subsurface, m η represents the polarization rate model of the subsurface, Φ(m ρ ) represents the objective function of the resistivity method in the joint inversion of the resistivity method and the induced polarization method, Φ(m η ) represents the objective function of the induced polarization method in the joint inversion of the resistivity method and the induced polarization method, Φ d (m ρ ) represents the data term of the resistivity method, Φ d (m η ) represents the data term of the induced polarization method, Φ m (m ρ ) represents the model term of the resistivity method, Φ m (m η ) represents the model term of the induced polarization method, Φ cg (, ρ , m η ) represents the cross-gradient term of the resistivity method, Φ cg (m η , m ρIt represents the cross-gradient term of the induced polarization method. λ1 represents the second weight factor of the model term in the objective function of the resistivity method. λ2 represents the second weight factor of the model term in the objective function of the induced polarization method. γ1 represents the weight constant of the cross-gradient term in the objective function of the resistivity method. γ2 represents the weight constant of the cross-gradient term in the objective function of the induced polarization method. K represents the coefficient matrix.

[0021] In an alternative embodiment, the study area includes a first sub-study area, a second sub-study area, and a third sub-study area. In the X-Y direction, the center of the first sub-study area coincides with the center of the study area. The second sub-study area is wrapped around the periphery of the first sub-study area. The third sub-study area is wrapped around the periphery of the second sub-study area. In the Z direction, the detection depths corresponding to the first sub-study area, the second sub-study area, and the third sub-study area gradually increase. The data range of the first sub-study area and the area range of the second sub-study area constitute the data range of the observation area. The grid cells in the first sub-study area correspond to the values of the first matrix elements. The grid cells in the second sub-study area correspond to the values of the second matrix elements. The grid cells in the third sub-study area correspond to the values of the third matrix elements. The value of the first matrix element is greater than the value of the second matrix element. The value of the second matrix element is greater than the value of the third matrix element.

[0022] In an alternative embodiment, the X-direction data range of the first sub-study area is 2 / 3 of the X-direction data range of the observation area. The Y-direction data range of the first sub-study area is 1 / 2 of the Y-direction data range of the observation area. The Z-direction data range of the first sub-study area is 2 / 3 of the Z-direction data range of the observation area.

[0023] In an alternative embodiment, the value of the first matrix element is 1, the value of the second matrix element is 0.1, and the value of the third matrix element is 0.

[0024] In an alternative embodiment, the set value of the weight constant makes the proportion of the cross-gradient term in the objective function between 5% and 10%.

[0025] In a second aspect, the embodiments of the present application further provide a joint inversion device based on the cross-gradient term weight factor. The joint inversion device includes:

[0026] A data acquisition module, configured to obtain apparent resistivity data and apparent polarization rate data according to multiple survey lines pre-arranged in the observation area of the underground medium.

[0027] A model construction module, configured to assign initial parameter values of a homogeneous half-space to each grid cell in an underground grid as an initial model for joint inversion. The initial parameter values include an initial resistivity value and an initial polarizability value. The underground grid is obtained by partitioning a research area, and the data range of the research area is determined according to the data range of the observation area.

[0028] A function construction module, configured to construct an objective function for joint inversion. The objective function includes a cross-gradient term and a first weight factor corresponding to the cross-gradient term. The first weight factor includes a coefficient matrix and a weight constant. Each matrix element in the coefficient matrix corresponds to each grid cell in the underground grid, and the numerical values of the matrix elements corresponding to the grid cells in different regions of the underground grid are different. The weight constant is used to adjust the proportion of the cross-gradient term in the objective function.

[0029] A model update module, configured to perform a minimization calculation on the objective function by using the apparent resistivity data and the apparent polarizability data to update the initial model, and obtain a resistivity model and a polarizability model.

[0030] In a third aspect, an embodiment of the present application further provides an electronic device, including: a processor, a memory, and a bus. The memory stores machine-readable instructions executable by the processor. When the electronic device runs, the processor communicates with the memory through the bus. When the machine-readable instructions are executed by the processor, the steps of the joint inversion method as described above are executed.

[0031] In a fourth aspect, an embodiment of the present application further provides a computer-readable storage medium. A computer program is stored on the computer-readable storage medium. When the computer program is run by a processor, the steps of the joint inversion method as described above are executed.

[0032] The embodiment of the present application provides a joint inversion method and device based on the weight factor of the cross-gradient term. The joint inversion method includes: first, obtaining apparent resistivity data and apparent polarization rate data according to multiple survey lines pre-arranged in the observation area of the underground medium; then, assigning the initial parameter values of the homogeneous half-space to each grid cell in the underground grid as the initial model for joint inversion, where the initial parameter values include the initial resistivity value and the initial polarization rate value, and the underground grid is obtained by dividing the research area into grids, and the data range of the research area is determined according to the data range of the observation area; then, constructing the objective function for joint inversion, where the objective function includes the cross-gradient term and the first weight factor corresponding to the cross-gradient term, the first weight factor includes a coefficient matrix and a weight constant, each matrix element in the coefficient matrix corresponds to each grid cell in the underground grid, and the numerical values of the matrix elements corresponding to the grid cells in different regions of the underground grid are different, and the weight constant is used to adjust the proportion of the cross-gradient term in the objective function; finally, using the apparent resistivity data and the apparent polarization rate data to perform minimization calculation on the objective function to update the initial model and obtain the resistivity model and the polarization rate model. The present application uses the coefficient matrix and the weight constant to select the weight factor of the effective cross-gradient term, so as to play a better constraint role on the resistivity model and the polarization rate model, can reduce the non-uniqueness of joint inversion, and helps to improve the accuracy of the inversion result.

[0033] In order to make the above objects, features, and advantages of the present application more obvious and understandable, the following specifically gives preferred embodiments and, in conjunction with the accompanying drawings, makes a detailed description as follows. BRIEF DESCRIPTION OF THE DRAWINGS

[0034] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the following will briefly introduce the drawings required to be used in the embodiments. It should be understood that the following drawings only show some embodiments of the present application, and therefore should not be regarded as limiting the scope. For those of ordinary skill in the art, without creative efforts, other related drawings can also be obtained based on these drawings.

[0035] Figure 1 It is a flowchart of a joint inversion method based on the weight factor of the cross-gradient term provided by an embodiment of the present application;

[0036] Figure 2 It is a schematic diagram of the distribution of survey lines provided by an embodiment of the present application;

[0037] Figure 3 It is a schematic diagram of the division of different regions of the underground grid provided by an embodiment of the present application;

[0038] Figure 4 It is the classification of the research area in the X-Y direction in an embodiment of the present application;

[0039] Figure 5 Classification of the research area in the X-Z direction in the embodiments of this application

[0040] Figure 6 Comparison chart of resistivity and polarizability inversion results of the traditional method and the method of the present invention provided in the embodiments of this application

[0041] Figure 7 Schematic diagram of the cross-gradient values calculated from resistivity and polarizability obtained by the traditional method and the method of the present invention provided in the embodiments of this application

[0042] Figure 8 Schematic structural diagram of a joint inversion device based on a cross-gradient term weight factor provided in the embodiments of this application

[0043] Figure 9 Schematic structural diagram of an electronic device provided in the embodiments of this application Detailed implementation manners

[0044] To make the objectives, technical solutions, and advantages of the embodiments of this application clearer, the technical solutions in the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings in the embodiments of this application. Apparently, the described embodiments are only some of the embodiments of this application, rather than all of them. Usually, the components of the embodiments of this application described and illustrated herein can be arranged and designed in various different configurations. Therefore, the following detailed description of the embodiments of this application provided in the drawings is not intended to limit the scope of this application that is required to be protected, but merely represents the selected embodiments of this application. Based on the embodiments of this application, every other embodiment obtained by those skilled in the art without creative efforts belongs to the scope of protection of this application.

[0045] First, the applicable application scenarios of this application are introduced. This application can be applied to the field of geophysical exploration technology. Geophysical inversion is an important concept in geophysics. It refers to the process of using geophysical observation data (such as seismic wave records, gravity measurement data, magnetic measurement data, etc.) to infer the physical properties (such as velocity structure, density distribution, magnetic distribution, etc.) and geological structure of the Earth's interior. During the inversion process, there may be multiple different subsurface structures or models that can fit the same set of observation data, which is the non-uniqueness of geophysical inversion. The non-uniqueness of geophysical inversion stems from multiple aspects. Mathematically, both nonlinearity and ill-posedness can cause non-uniqueness of solutions; in physical models, the geological structure is complex and the properties are equivalent, and multiple combinations can produce similar responses; the observation data also has problems such as being limited, incomplete, and having insufficient precision and resolution. These factors are intertwined, resulting in multiple possible inversion results.

[0046] The cross-gradient joint inversion method is a currently mainstream structural coupling joint inversion method. By multiplying the gradients of different physical properties, the boundaries of different physical properties are mutually constrained, which can improve the similarity of the inversion results of different physical properties in terms of spatial structure. It mainly realizes the joint inversion process by introducing a cross-gradient term into the objective function. By minimizing the objective function, the underground model can be obtained.

[0047] In the process of minimizing the objective function, the selection of the weight factor is particularly important. If it is too small, its constraint effect will be limited; if it is too large, it will play the opposite role and reduce the accuracy of the inversion result. However, for the weight factor of the cross-gradient term, there is currently no good selection method. Because the cross-gradient term does not have values everywhere, and the cross-gradient has values only when both models change and the change directions are inconsistent. Therefore, how to select the weight factor of the cross-gradient term to ensure that the cross-gradient term plays a good role in the inversion process is of great significance.

[0048] Based on this, the embodiments of the present application provide a joint inversion method based on the weight factor of the cross-gradient term. By selecting an effective weight factor of the cross-gradient term, the non-uniqueness of the joint inversion is reduced, which helps to improve the accuracy of the inversion result.

[0049] Please refer to Figure 1 , Figure 1 which is a flowchart of a joint inversion method based on the weight factor of the cross-gradient term provided by the embodiments of the present application. As Figure 1 shown in, the joint inversion method provided by the embodiments of the present application includes:

[0050] S101. Obtain apparent resistivity data and apparent polarization data according to multiple survey lines arranged in advance in the observation area of the underground medium;

[0051] S102. Assign the initial parameter values of the homogeneous half-space to each grid cell in the underground grid as the initial model of the joint inversion. The initial parameter values include the initial resistivity value and the initial polarization value. The underground grid is obtained by dividing the research area into grids, and the data range of the research area is determined according to the data range of the observation area;

[0052] S103. Construct the objective function of the joint inversion. The objective function includes a cross-gradient term and a first weight factor corresponding to the cross-gradient term. The first weight factor includes a coefficient matrix and a weight constant. Each matrix element in the coefficient matrix corresponds to each grid cell in the underground grid. The numerical values of the matrix elements corresponding to the grid cells in different regions of the underground grid are different. The weight constant is used to adjust the proportion of the cross-gradient term in the objective function;

[0053] S104. Use the apparent resistivity data and apparent polarization rate data to perform minimization calculation on the objective function to update the initial model, and obtain a resistivity model and a polarization rate model.

[0054] In the above steps S101 to S104, a coefficient matrix and a weight constant are used to select the weight factor of the effective cross-gradient term, so as to play a better constraint role on the resistivity model and the polarization rate model, reduce the non-uniqueness of the joint inversion, and help improve the accuracy of the inversion result.

[0055] The following uses specific embodiments to exemplarily illustrate the above steps S101 to S104:

[0056] In step S101, according to multiple survey lines pre-arranged in the observation area of the underground medium, apparent resistivity data and apparent polarization rate data are obtained.

[0057] Here, a survey line is a line arranged along a straight line at a certain scale for observing the underground medium and consists of a series of observation points. The underground medium includes but is not limited to mineral resources, oil and gas resources, groundwater, and archaeological remains. Exemplarily, taking groundwater as an example, different geological structures will affect the migration of groundwater. Arranging survey lines to obtain apparent resistivity data and apparent polarization rate data can help analyze the underground geological structures, such as the location and scale of faults, fractures, etc. Among them, the arrangement area of the survey line is defined as the observation area.

[0058] In an optional embodiment, multiple survey lines are arranged in the observation area of the underground medium, and a transmitting electrode and a receiving electrode are provided on each survey line. Step S101 specifically includes:

[0059] When the transmitting electrode on each survey line is powered, the adjacent receiving electrode on the corresponding survey line receives to obtain apparent resistivity data and apparent polarization rate data.

[0060] For example, as Figure 2 shown, a total of 5 survey lines L1 to L5 are arranged for groundwater detection in the Chaobai River Basin in Miyun District, Beijing. The triangles on the survey lines represent the positions of the transmitting electrodes, and the circles represent the positions of the receiving electrodes. Each survey line is observed using a three-electrode device. When the transmitting electrode is powered, the adjacent receiving electrode on the corresponding survey line receives. Finally, a total of 24,540 apparent resistivity data and apparent polarization rate data are obtained. The apparent resistivity data and apparent polarization rate data are sorted out to obtain the apparent resistivity value, apparent polarization rate value, observation error, transmitting point coordinates, and receiving point coordinates of each data. Among them, the receiving point coordinates range from -96m to 96m in the X direction and from -70m to 70m in the Y direction, and the detection depth is 0 to 80m.

[0061] In step S102, the initial parameter values ​​of the uniform half space are assigned to each grid cell in the underground grid as the initial model of the joint inversion. The initial parameter values ​​include initial resistivity values ​​and initial polarizability values. The underground grid is obtained by gridding the study area, and the data range of the study area is determined according to the data range of the observation area.

[0062] Here, the uniform half-space means that the medium in the half-space has uniform physical properties. For example, in terms of electrical properties, the resistivity in the entire half-space is the same, and there is no change in resistivity; in terms of polarization properties, the polarizability in the entire half-space is also the same. In the embodiment of the present application, a uniform half-space with a resistivity of 200Ω·m and a polarizability of 0.01 is selected as the initial model.

[0063] Among them, the main purpose of gridding the study area is to facilitate the mathematical description and numerical calculation of the underground space. Since the underground space is a continuous three-dimensional area, it is very difficult to analyze and calculate it directly. Through grid division, the continuous underground space is discretized into many small units, so that each unit can be analyzed and processed separately, thereby simplifying the calculation process. For example, there are 155 grid nodes in the X direction, 19 grid nodes in the Y direction, and 90 grid nodes in the Z direction in the embodiment of the present application.

[0064] After the embodiment of the present application performs a grid division operation on the study area, the obtained discrete small unit set constitutes an underground grid. Each grid unit has its specific position and number, and can be assigned corresponding physical parameters (such as resistivity, polarizability). By calculating and analyzing each grid unit, the inversion result of the underground situation of the entire study area can be finally obtained.

[0065] In other words, the underground space is divided into many grid units, and a value is first assigned to each grid unit, which is the initial model. The inversion is to continuously change the value of each grid unit by minimizing the objective function, and finally obtain a model that meets the minimization of the objective function.

[0066] Optionally, step S102 determines the data range of the study area by the following steps:

[0067] Step 1021: Obtain the data range in the X direction, the data range in the Y direction, and the data range in the Z direction of the observation area.

[0068] Here, the X-direction data range, Y-direction data range, and Z-direction data range of the observation area can be determined according to the X-direction range, Y-direction data range, and Z-direction data range of the receiving point coordinates. Figure 2Based on the data collected by the surveyed line shown, the coordinates of the receiving points range from -96m to 96m in the X direction, from -70m to 70m in the Y direction, and the detection depth is from 0 to 80m. Furthermore, the data range of the observation area in the X direction is from -96m to 96m, in the Y direction is from -70m to 70m, and in the Z direction is from 0 to 80m.

[0069] Step 1022: Set target multiples for the data ranges in the X direction, Y direction, and Z direction to obtain the data ranges of the study area in the X direction, Y direction, and Z direction; wherein, the target multiple is between 3 times and 5 times.

[0070] Here, the data ranges of the study area in the X direction, Y direction, and Z direction are determined according to the data ranges of the observation area in the X direction, Y direction, and Z direction, and the data ranges of the study area in the X direction, Y direction, and Z direction are the data ranges of the initial model in the X direction, Y direction, and Z direction. Exemplarily, Figure 2 Based on the data collected by the surveyed line shown, the data range of the study area in the X direction can be defined as from -300m to 300m, in the Y direction as from -300m to 300m, and in the Z direction as from 0 to 250m.

[0071] In the embodiments of the present application, if the entire spatial range of the study area is set too small, boundary effects will occur during the inversion process and the inversion result will be poor. Therefore, the data range of the study area is generally set to 3 to 5 times the data range of the observation area, which can obtain a more accurate inversion result, thereby reducing the boundary effect during the inversion process and improving the accuracy of the inversion result.

[0072] In step S103, construct an objective function for joint inversion. The objective function includes a cross-gradient term and a first weight factor corresponding to the cross-gradient term. The first weight factor includes a coefficient matrix and a weight constant. Each matrix element in the coefficient matrix corresponds to each grid unit in the underground grid. The numerical values of the matrix elements corresponding to the grid units in different regions of the underground grid are different. The weight constant is used to adjust the proportion of the cross-gradient term in the objective function.

[0073] In this step, the first weight factor corresponding to the cross-gradient term introduces the coefficient matrix K and constants γ1 and γ2. The number of matrix elements in the coefficient matrix K is the same as the number of grids in the underground grid, so that we can set different numerical values for each grid unit of the underground grid. γ1 and γ2 are used to adjust the proportion of the entire cross-gradient term in the objective function.

[0074] Specifically, the matrix element values corresponding to the grid cells in different regions of the underground grid are different. Each grid cell in the underground grid is assigned a value according to the region where it is located. Each grid cell in the underground grid corresponds to each matrix element in the coefficient matrix, and thus the coefficient matrix can be obtained.

[0075] Exemplarily, as Figures 3 to 5 shown, the entire cube represents the research area, which includes the first sub-research area A, the second sub-research area B, and the third sub-research area C. As Figure 4 shown, in the X-Y direction, the center of the first sub-research area A coincides with the center of the research area. The second sub-research area B is wrapped around the periphery of the first sub-research area A, and the third sub-research area C is wrapped around the periphery of the second sub-research area B; as Figure 5 shown, in the Z direction, the detection depths corresponding to the first sub-research area A, the second sub-research area B, and the third sub-research area C gradually increase; the data range of the first sub-research area A and the area range of the second sub-research area B constitute the data range of the observation area; the grid cells in the first sub-research area A correspond to the first matrix element value, the grid cells in the second sub-research area B correspond to the second matrix element value, and the grid cells in the third sub-research area C correspond to the third matrix element value. The first matrix element value is greater than the second matrix element value, and the second matrix element value is greater than the third matrix element value.

[0076] Optionally, the X-direction data range of the first sub-research area A is 2 / 3 of the X-direction data range of the observation area, the Y-direction data range of the first sub-research area A is 1 / 2 of the Y-direction data range of the observation area, and the Z-direction data range of the first sub-research area A is 2 / 3 of the Z-direction data range of the observation area.

[0077] Exemplarily, the first matrix element value can be 1, the second matrix element value can be 0.1, and the third matrix element value can be 0.

[0078] For example, taking the data collected by the survey line shown in Figure 2 as the analysis basis, and as Figures 4 to 5As shown in the figure, the underground grid is divided into three different regions A, B, and C. Among them, A represents the first research region, whose coordinate range is -60m to 60m in the X direction, -30m to 30m in the Y direction, and 0m to 60m in the Z direction. The values corresponding to the grid cells in the first research region in the coefficient matrix K are set to 1; B represents the second research region, whose coordinate range is -96m to 96m in the X direction, -70m to -30m and 30m to 70m in the Y direction, and 60m to 80m in the Z direction. The values corresponding to the grid cells in the second research region in the coefficient matrix K are set to 0.1; C represents the third research region, whose coordinate range is -300m to -96m and 96m to 300m in the X direction, -300m to -70m and 70m to 300m in the Y direction, and 80m to 250m in the Z direction. The values corresponding to the grid cells in the third research region in the coefficient matrix K are set to 0.

[0079] The above method makes the constraint effect of the cross-gradient term mainly concentrated in the first research region during the joint inversion process, avoiding the large influence of the data observation boundary on the inversion result. The resistivity and polarization rate inversion results finally obtained are also more similar in structure within the first research region.

[0080] Optionally, the set value of the weight constant makes the proportion of the cross-gradient term in the objective function between 5% and 10%.

[0081] In the embodiment of the present application, the weight constant γ1 of the resistivity model and the weight constant γ2 of the polarization rate model are set mainly to adjust the proportion of the cross-gradient term in the objective function. Through a large number of tests in the present invention, it is found that when the proportion of the entire cross-gradient term in the objective function is 5% to 10% by adjusting γ1 and γ2, a better inversion effect can be obtained. In the embodiment of the present application, based on the data collected by the Figure 2 survey line shown, γ1 can be selected as 500 and γ2 can be selected as 2000.

[0082] In an alternative embodiment, the objective function further includes a data term, a model term, and a second weight factor corresponding to the model term. The data term is used to calculate the difference between the apparent resistivity data and the resistivity data of the inversion response, and the difference between the apparent polarization rate data and the polarization rate data of the inversion response. The model term is used to calculate the smoothness of the resistivity model and the polarization rate model. The second weight factor represents the proportion of the model term in the objective function.

[0083] Here, the data term, the model term, and the second weight factor cooperate with each other in geophysical inversion and act on the objective function together to guide the optimization and adjustment of the model, and finally obtain a more accurate and reasonable underground model.

[0084] Exemplarily, the objective function is represented by the following formula:

[0085] Φ(m ρ ) = Φ d (m ρ ) + λ1Φ n (m ρ ) + γ1KΦ cg (m ρ ,, η )

[0086] Φ(, η ) = Φ d (m η ) + λ2Φ m (m η ) + γ2KΦ cg (m η ,m ρ )

[0087] Among them, m ρ represents the resistivity model of the subsurface, m η represents the polarizability model of the subsurface, Φ(m ρ ) represents the objective function of the resistivity method in the joint inversion of the resistivity method and induced polarization method, Φ(m η ) represents the objective function of the induced polarization method in the joint inversion of the resistivity method and induced polarization method, Φ d (m ρ ) represents the data term of the resistivity method, Φ d (m η ) represents the data term of the induced polarization method, Φ m (m ρ ) represents the model term of the resistivity method, Φ m (m η ) represents the model term of the induced polarization method, Φ cg (m ρ ,m η ) represents the cross-gradient term of the resistivity method, Φ cg (m η ,m ρ ) represents the cross-gradient term of the induced polarization method, λ1 represents the second weight factor of the model term in the objective function of the resistivity method, λ2 represents the second weight factor of the model term in the objective function of the induced polarization method, γ1 represents the weight constant of the cross-gradient term in the objective function of the resistivity method, γ2 represents the weight constant of the cross-gradient term in the objective function of the induced polarization method, and K represents the coefficient matrix.

[0088] Here, the inversion is achieved by continuously minimizing the objective function, and m ρ and m η will change continuously, and finally the resistivity model and polarizability model that best match the subsurface conditions are obtained.

[0089] In step S104, using the apparent resistivity data and the apparent polarizability data, the objective function is minimized to update the initial model, obtaining a resistivity model and a polarizability model.

[0090] In this step, through multiple iterations and minimization calculations of the objective function, when certain convergence conditions are met (such as the objective function value no longer significantly decreasing, the parameter change amount being less than a certain threshold, etc.), the final resistivity model and polarizability model are obtained. The resistivity model describes the distribution of the resistivity of the underground medium, and the polarizability model describes the distribution of the polarizability of the underground medium. These models can intuitively display the electrical properties at different positions underground and provide important reference bases for geological interpretation, resource exploration (such as searching for mineral resources, underground water resources, etc.), engineering construction (such as evaluating the foundation stability, etc.).

[0091] For example, as Figure 6 shown in the comparison chart of the resistivity and polarizability inversion results of the traditional method and the method of the present invention, it can be seen that compared with the traditional method, the structural similarity of the resistivity and polarizability models obtained by the method used in the present invention is significantly improved. As Figure 7 shown in the schematic diagram of the cross-gradient values calculated by the resistivity and polarizability of the traditional method and the method of the present invention, it can be seen that the cross-gradient value of the method of the present invention in the first research area is smaller, indicating that the cross-gradient term in the method of the present invention has a more obvious effect in the first research area. In the third research area, the cross-gradient value basically remains unchanged, indicating that the method of the present invention avoids the influence of the third research area on the first research area.

[0092] The joint inversion method based on the cross-gradient term weight factor provided by the embodiments of the present application uses a coefficient matrix and a weight constant to select the weight factor of the effective cross-gradient term, thereby playing a better constraint role on the resistivity model and the polarizability model, being able to reduce the non-uniqueness of the joint inversion, and contributing to improving the accuracy of the inversion result.

[0093] Based on the same inventive concept, the embodiments of the present application also provide a joint inversion device based on the cross-gradient term weight factor corresponding to the joint inversion method based on the cross-gradient term weight factor. Since the principle of solving problems by the device in the embodiments of the present application is similar to the above joint inversion method based on the cross-gradient term weight factor in the embodiments of the present application, the implementation of the device can refer to the implementation of the method, and the repeated parts will not be described again.

[0094] Please refer to Figure 8 , Figure 8 which is a schematic structural diagram of a joint inversion device based on the cross-gradient term weight factor provided by the embodiments of the present application. As Figure 8 shown in

[0095] A data acquisition module 801, configured to obtain apparent resistivity data and apparent polarization rate data according to a plurality of survey lines pre-arranged in an observation area of a subsurface medium;

[0096] A model construction module 802, configured to assign initial parameter values of a homogeneous half-space to each grid cell in an underground grid as an initial model for joint inversion, where the initial parameter values include an initial resistivity value and an initial polarization rate value, and the underground grid is obtained by dividing a research area into grids, and the data range of the research area is determined according to the data range of the observation area;

[0097] A function construction module 803, configured to construct an objective function for joint inversion, where the objective function includes a cross-gradient term and a first weight factor corresponding to the cross-gradient term, the first weight factor includes a coefficient matrix and a weight constant, each matrix element in the coefficient matrix corresponds to each grid cell in the underground grid, and the numerical values of the matrix elements corresponding to the grid cells in different areas of the underground grid are different, and the weight constant is used to adjust the proportion of the cross-gradient term in the objective function;

[0098] A model update module 804, configured to perform minimization calculation on the objective function by using the apparent resistivity data and the apparent polarization rate data to update the initial model to obtain a resistivity model and a polarization rate model.

[0099] In an optional embodiment, a plurality of survey lines are arranged in an observation area of a subsurface medium, and a transmitting electrode and a receiving electrode are arranged on each survey line. The data acquisition module 801 is specifically configured to:

[0100] When the transmitting electrode on each survey line supplies power, adjacent receiving electrodes on the corresponding survey line receive to obtain apparent resistivity data and apparent polarization rate data.

[0101] In an optional embodiment, the model construction module 802 is configured to determine the data range of the research area through the following steps:

[0102] Obtain the X-direction data range, Y-direction data range, and Z-direction data range of the observation area;

[0103] Set target multiples for the X-direction data range, the Y-direction data range, and the Z-direction data range to obtain the data ranges of the research area in the X direction, Y direction, and Z direction; where the target multiples are between 3 times and 5 times.

[0104] In an alternative embodiment, the objective function further includes a data term, a model term, and a second weight factor corresponding to the model term. The data term is used to calculate the difference between the apparent resistivity data and the resistivity data of the inversion response, and the difference between the apparent polarization rate data and the polarization rate data of the inversion response. The model term is used to calculate the smoothness of the resistivity model and the polarization rate model. The second weight factor characterizes the proportion of the model term in the objective function.

[0105] In an alternative embodiment, the objective function is expressed by the following formula:

[0106] Φ(m ρ ) = Φ d (m ρ ) + λ1Φ m (m ρ ) + γ1KΦ cg (m ρ , m η )

[0107] Φ(m η ) = Φ d (m η ) + λ2Φ m (m η ) + γ2KΦ cg (m η , m ρ )

[0108] where m ρ represents the resistivity model of the subsurface, m η represents the polarization rate model of the subsurface, Φ(m ρ ) represents the objective function of the resistivity method in the joint inversion of the resistivity method and the induced polarization method, Φ(m η ) represents the objective function of the induced polarization method in the joint inversion of the resistivity method and the induced polarization method, Φ d (m ρ ) represents the data term of the resistivity method, Φ d (m η ) represents the data term of the induced polarization method, Φ m (m ρ ) represents the model term of the resistivity method, Φ m (m η ) represents the model term of the induced polarization method, Φ cg (m ρ , m η ) represents the cross-gradient term of the resistivity method, Φ cg (m η , m ρ) represents the cross-gradient term of the induced polarization method, λ1 represents the second weight factor of the model term in the objective function of the resistivity method, λ2 represents the second weight factor of the model term in the objective function of the induced polarization method, γ1 represents the weight constant of the cross-gradient term in the objective function of the resistivity method, γ2 represents the weight constant of the cross-gradient term in the objective function of the induced polarization method, and K represents the coefficient matrix.

[0109] In an alternative embodiment, the study area includes a first sub-study area, a second sub-study area, and a third sub-study area. In the X-Y direction, the center of the first sub-study area coincides with the center of the study area, the second sub-study area is wrapped around the periphery of the first sub-study area, and the third sub-study area is wrapped around the periphery of the second sub-study area; in the Z direction, the detection depths corresponding to the first sub-study area, the second sub-study area, and the third sub-study area gradually increase; the data range of the first sub-study area and the area range of the second sub-study area constitute the data range of the observation area; the grid cells in the first sub-study area correspond to the values of the first matrix elements, the grid cells in the second sub-study area correspond to the values of the second matrix elements, the grid cells in the third sub-study area correspond to the values of the third matrix elements, the value of the first matrix element is greater than the value of the second matrix element, and the value of the second matrix element is greater than the value of the third matrix element.

[0110] In an alternative embodiment, the X-direction data range of the first sub-study area is 2 / 3 of the X-direction data range of the observation area, the Y-direction data range of the first sub-study area is 1 / 2 of the Y-direction data range of the observation area, and the Z-direction data range of the first sub-study area is 2 / 3 of the Z-direction data range of the observation area.

[0111] In an alternative embodiment, the value of the first matrix element is 1, the value of the second matrix element is 0.1, and the value of the third matrix element is 0.

[0112] In an alternative embodiment, the set value of the weight constant is such that the proportion of the cross-gradient term in the objective function is between 5% and 10%.

[0113] The joint inversion device based on the weight factor of the cross-gradient term provided by the embodiments of the present application uses the coefficient matrix and the weight constant to select the weight factor of the effective cross-gradient term, thereby playing a better constraint role on the resistivity model and the polarizability model, being able to reduce the non-uniqueness of the joint inversion, and helping to improve the accuracy of the inversion result.

[0114] Please refer to Figure 9 , Figure 9The following is a schematic structural diagram of an electronic device provided by an embodiment of the present application. As Figure 9 shown in the figure, the electronic device 900 includes a processor 901, a memory 902, and a bus 903.

[0115] The memory 902 stores machine-readable instructions executable by the processor 901. When the electronic device 900 runs, the processor 901 communicates with the memory 902 through the bus 903. When the machine-readable instructions are executed by the processor 901, the steps of the joint inversion method based on the cross-gradient term weight factor in the method embodiment as described above can be executed. For the specific implementation manner, reference can be made to the method embodiment, which will not be elaborated herein. Figure 1 shown in the figure, the steps of the joint inversion method based on the cross-gradient term weight factor in the method embodiment as described above can be executed. For the specific implementation manner, reference can be made to the method embodiment, which will not be elaborated herein.

[0116] An embodiment of the present application further provides a computer-readable storage medium. A computer program is stored on the computer-readable storage medium. When the computer program is run by a processor, the steps of the joint inversion method based on the cross-gradient term weight factor in the method embodiment as described above can be executed. For the specific implementation manner, reference can be made to the method embodiment, which will not be elaborated herein. Figure 1 shown in the figure, the steps of the joint inversion method based on the cross-gradient term weight factor in the method embodiment as described above can be executed. For the specific implementation manner, reference can be made to the method embodiment, which will not be elaborated herein.

[0117] Those skilled in the art can clearly understand that for the convenience and brevity of description, the specific working processes of the above-described systems, devices, and units can refer to the corresponding processes in the foregoing method embodiments, which will not be elaborated herein.

[0118] In several embodiments provided by the present application, it should be understood that the disclosed systems, devices, and methods can be implemented in other ways. The device embodiments described above are merely illustrative. For example, the division of the units is only a logical function division. In actual implementation, there may be other division methods. For another example, multiple units or components can be combined or integrated into another system, or some features can be ignored or not executed. Another point is that the displayed or discussed mutual coupling or direct coupling or communication connection can be through some communication interfaces. The indirect coupling or communication connection of the devices or units can be in an electrical, mechanical, or other form.

[0119] The units described as separate components may or may not be physically separated. The components displayed as units may or may not be physical units, that is, they can be located in one place, or they can be distributed to multiple network units. Some or all of the units can be selected according to actual needs to achieve the purpose of the solution of this embodiment.

[0120] In addition, in each embodiment of the present application, the functional units can be integrated into one processing unit, or each unit can exist physically alone, or two or more units can be integrated into one unit.

[0121] When the above-mentioned functions are implemented in the form of software function units and sold or used as independent products, they can be stored in a non-volatile computer-readable storage medium executable by a processor. Based on such understanding, the technical solution of the present application, in essence, or the part that contributes to the prior art, or a part of the technical solution, can be embodied in the form of a software product. The computer software product is stored in a storage medium and includes several instructions for causing a computer device (which may be a personal computer, a server, or a network device, etc.) to execute all or part of the steps of the methods described in various embodiments of the present application. The foregoing storage medium includes: various media such as USB flash drives, mobile hard disks, read-only memories (ROM), random access memories (RAM), magnetic disks, or optical discs that can store program codes.

[0122] Finally, it should be noted that: the above-mentioned embodiments are only specific implementation manners of the present application, used to illustrate the technical solutions of the present application, rather than limiting them. The protection scope of the present application is not limited thereto. Although the present application has been described in detail with reference to the foregoing embodiments, those of ordinary skill in the art should understand that: any person skilled in the art within the technical scope disclosed by the present application can still modify the technical solutions recorded in the foregoing embodiments, or can easily think of changes, or perform equivalent replacements on some of the technical features; and these modifications, changes, or replacements do not make the essence of the corresponding technical solutions deviate from the spirit and scope of the technical solutions of the embodiments of the present application, and should all be covered within the protection scope of the present application. Therefore, the protection scope of the present application shall be subject to the protection scope of the claims.

Claims

1. A joint inversion method based on cross gradient term weight factor, characterized in that: The joint inversion method comprises: Obtaining apparent resistivity data and apparent polarizability data according to a plurality of survey lines pre-arranged in an observation area of ​​the underground medium; Assigning initial parameter values ​​of a uniform half space to each grid cell in a subsurface grid as an initial model for joint inversion, wherein the initial parameter values ​​include initial resistivity values ​​and initial polarizability values, wherein the subsurface grid is obtained by gridding a study area, and the data range of the study area is determined according to the data range of the observation area; Constructing an objective function of joint inversion, wherein the objective function includes a cross-gradient term and a first weight factor corresponding to the cross-gradient term, wherein the first weight factor includes a coefficient matrix and a weight constant, wherein each matrix element in the coefficient matrix corresponds to each grid cell in the underground grid, and the values ​​of the matrix elements corresponding to the grid cells in different areas of the underground grid are different, and the weight constant is used to adjust the proportion of the cross-gradient term in the objective function; The objective function is minimized by using the apparent resistivity data and the apparent polarizability data to update the initial model and obtain a resistivity model and a polarizability model.

2. The joint inversion method according to claim 1, characterized in that: A plurality of survey lines are arranged in the observation area of ​​the underground medium, and a transmitting electrode and a receiving electrode are arranged on each survey line. The method of obtaining the apparent resistivity data and the apparent polarizability data according to the plurality of survey lines pre-arranged in the observation area of ​​the underground medium includes: When the transmitting electrode on each measuring line is powered, the adjacent receiving electrode on the corresponding measuring line receives the data to obtain the apparent resistivity data and the apparent polarizability data.

3. The joint inversion method according to claim 1, characterized in that: The data extent of the study area is determined by the following steps: Obtaining the data range in the X direction, the data range in the Y direction, and the data range in the Z direction of the observation area; A target multiple is set for the X-direction data range, the Y-direction data range and the Z-direction data range to obtain the data range of the study area in the X-direction, the Y-direction and the Z-direction; wherein the target multiple is between 3 times and 5 times.

4. The joint inversion method according to claim 1, characterized in that: The objective function also includes a data item, a model item and a second weight factor corresponding to the model item. The data item is used to calculate the difference between the apparent resistivity data and the resistivity data of the inversion response, and the difference between the apparent polarizability data and the polarizability data of the inversion response. The model item is used to calculate the smoothness of the resistivity model and the polarizability model. The second weight factor represents the proportion of the model item in the objective function.

5. The joint inversion method according to claim 4, characterized in that: The objective function is expressed by the following formula: Φ(m ρ )=Φ d (m ρ )+λ1Φ m (m ρ )+γ1KΦ cg (m ρ ,m η ) Φ(m η )=Φ d (m η )+λ2Φ m (m η )+γ2KΦ cg (m η ,m ρ ) Among them, m ρ Represents the resistivity model of the underground, m η represents the underground polarizability model, Φ(m ρ ) represents the objective function of the resistivity method in the joint inversion of the resistivity method and the induced polarization method, Φ(m η ) represents the objective function of the induced polarization method in the joint inversion of the resistivity method and the induced polarization method, Φ d (m ρ ) represents the data item of the resistivity method, Φ d (m η ) represents the data item of induced polarization method, Φ m (m ρ ) represents the model term of the resistivity method, Φ m (m η ) represents the model term of the induced polarization method, Φ cg (m ρ ,m η ) represents the cross gradient term of the resistivity method, Φ cg (m η ,m ρ ) represents the cross-gradient term of the induced polarization method, λ1 represents the second weight factor of the model term in the objective function of the resistivity method, λ2 represents the second weight factor of the model term in the objective function of the induced polarization method, γ1 represents the weight constant of the cross-gradient term in the objective function of the resistivity method, γ2 represents the weight constant of the cross-gradient term in the objective function of the induced polarization method, and K represents the coefficient matrix.

6. The joint inversion method according to claim 1, characterized in that: The research area includes a first sub-study area, a second sub-study area and a third sub-study area. In the XY direction, the center of the first sub-study area coincides with the center of the research area, the second sub-study area is wrapped around the periphery of the first sub-study area, and the third sub-study area is wrapped around the periphery of the second sub-study area. In the Z direction, the detection depths corresponding to the first sub-study area, the second sub-study area and the third sub-study area respectively increase gradually. The data range of the first sub-study area and the area range of the second sub-study area constitute the data range of the observation area. The grid cells in the first sub-study area correspond to the first matrix element values, the grid cells in the second sub-study area correspond to the second matrix element values, and the grid cells in the third sub-study area correspond to the third matrix element values. The first matrix element value is greater than the second matrix element value, and the second matrix element value is greater than the third matrix element value.

7. The joint inversion method according to claim 6, characterized in that: The data range of the first sub-study area in the X direction is 2 / 3 of the data range of the observation area in the X direction, the data range of the first sub-study area in the Y direction is 1 / 2 of the data range of the observation area in the Y direction, and the data range of the first sub-study area in the Z direction is 2 / 3 of the data range of the observation area in the Z direction.

8. The joint inversion method according to claim 6, characterized in that: The value of the first matrix element is 1, the value of the second matrix element is 0.1, and the value of the third matrix element is 0.

9. The joint inversion method according to claim 1, characterized in that: The setting value of the weight constant makes the proportion of the cross gradient term in the objective function between 5% and 10%.

10. A joint inversion device based on cross gradient term weight factor, characterized in that: The joint inversion device comprises: A data acquisition module, used to acquire apparent resistivity data and apparent polarizability data according to a plurality of survey lines pre-arranged in an observation area of ​​the underground medium; A model building module, used for assigning initial parameter values ​​of a uniform half space to each grid cell in a subsurface grid as an initial model for joint inversion, wherein the initial parameter values ​​include initial resistivity values ​​and initial polarizability values, the subsurface grid is obtained by gridding a study area, and the data range of the study area is determined according to the data range of the observation area; A function construction module, used to construct an objective function of joint inversion, wherein the objective function includes a cross-gradient term and a first weight factor corresponding to the cross-gradient term, wherein the first weight factor includes a coefficient matrix and a weight constant, wherein each matrix element in the coefficient matrix corresponds to each grid cell in the underground grid, and the values ​​of the matrix elements corresponding to the grid cells in different areas of the underground grid are different, and the weight constant is used to adjust the proportion of the cross-gradient term in the objective function; The model updating module is used to use the apparent resistivity data and the apparent polarizability data to perform a minimization calculation on the objective function to update the initial model and obtain a resistivity model and a polarizability model.

Images