Hydrogeological parameter space-time heterogeneity geostatistics inversion method

Abstract

The invention discloses a hydrogeological parameter space-time heterogeneity geostatistics inversion method. The method comprises the following steps: constructing a to-be-solved hydrogeological parameter column vector based on a geologic model and a'space-time equivalence 'principle; a priori covariance matrix of partition parameters in a time domain is assembled in a blocking mode; carrying out dimension reduction on the covariance matrix based on a principal component analysis method; carrying out parameter sensitivity matrix approximate calculation, solving a collaborative Kriging equation, and calculating an estimated value of a to-be-solved parameter by utilizing an iteration mode; and determining a target function, and performing iterative updating according to the target function to obtain a to-be-solved parameter column vector. The invention designs a partition geological statistics dimension reduction inversion method for describing the spatial-temporal heterogeneity of the hydrogeological parameters, and provides an efficient and accurate solution for describing the spatial-temporal heterogeneity of the parameters in a complex spatial-temporal process.

Description

Technical Field

[0001] The present invention relates to a geostatistical inversion method for spatio-temporal heterogeneity of hydrogeological parameters, and belongs to the technical field of hydrogeological parameter characterization. Background Art

[0002] The accurate characterization of hydrogeological parameters has an important impact on the research of underground seepage fields, solute transport, and engineering groundwater control. Existing hydrogeological parameter characterization methods mostly only consider the distribution of parameters in the spatial domain or the temporal domain. However, for complex hydrogeological spatio-temporal processes, parameters will have both spatial heterogeneity and temporal heterogeneity. Characterizing the spatio-temporal heterogeneity of parameters often involves a large number of parameters to be solved, which increases the computational cost of parameter inversion and makes it difficult for traditional parameter inversion methods to meet the higher-precision calculation requirements. Therefore, for the characterization of spatio-temporal heterogeneity of parameters under complex hydrogeological processes, there is still a lack of efficient and accurate inversion methods. Summary of the Invention

[0003] The technical problem to be solved by the present invention is to provide a geostatistical inversion method for spatio-temporal heterogeneity of hydrogeological parameters to achieve efficient and accurate characterization of spatio-temporal heterogeneity of hydrogeological parameters under complex spatio-temporal processes.

[0004] The present invention adopts the following technical solutions to solve the above technical problems:

[0005] A geostatistical inversion method for spatio-temporal heterogeneity of hydrogeological parameters includes the following steps:

[0006] Step 1: Based on the geological model of the target geological entity, the target geological entity is partitioned in the spatial domain, and the hydrogeological parameters to be solved in each partition are the same. Based on the spatio-temporal equivalence principle, the dynamic changes of the hydrogeological parameters to be solved in each partition in the temporal domain are discretized into column vectors, thereby constructing column vectors of the hydrogeological parameters to be solved in all partitions;

[0007] Step 2: For the hydrogeological parameters to be solved in each partition, a variogram in the temporal domain is constructed, and the corresponding covariance matrix is calculated according to the variogram. The prior covariance matrix of the hydrogeological parameters to be solved in all partitions is constructed by a block assembly method;

[0008] Step 3: Based on the principal component analysis method, the prior covariance matrix is dimensionally reduced to obtain a dimensionally reduced prior covariance matrix;

[0009] Step 4: Set an initial value of the hydrogeological parameters to be solved, and based on the dimensionally reduced prior covariance matrix, use an iterative method to calculate the estimated value of the hydrogeological parameters to be solved;

[0010] Step 5: Construct the objective function for the inversion of the hydrogeological parameters to be determined. Substitute the estimated values of the hydrogeological parameters to be determined obtained in each iteration into the objective function, compare the objective function values of the previous and current iterations. When the preset number of iterations is reached or the objective function values of the previous and current iterations are within the tolerance range, the iteration terminates, and the estimated value obtained in the last iteration is taken as the final estimated value;

[0011] Step 6: After the iteration terminates, calculate the posterior covariance matrix of the hydrogeological parameters to be determined, thereby obtaining the value range of the hydrogeological parameters to be determined.

[0012] Compared with the prior art, the present invention adopts the above technical solutions and has the following technical effects:

[0013] 1. The present invention constructs the spatio-temporal heterogeneous parameter vector to be determined through the zonal geological model and the "spatio-temporal equivalence" principle, constructs the prior covariance matrix and reduces its dimension based on principal component analysis, and avoids the calculation of the parameter sensitivity matrix, realizing the high-precision characterization of the spatio-temporal heterogeneity of hydrogeological parameters. Compared with the parameter inversion methods that only consider the spatial domain or only consider the time domain, the inversion accuracy is higher.

[0014] 2. The inversion idea based on dimension reduction of the present invention improves the inversion efficiency of the method and at the same time improves the applicability to the complex spatio-temporal variation process of parameters.

[0015] 3. The method of the present invention is applicable to various parameter inversion scenarios in the hydrogeological field and can provide more appropriate inversion results for different application scenarios such as engineering and scientific research. BRIEF DESCRIPTION OF THE DRAWINGS

[0016] Figure 1 is a flowchart of the geostatistical inversion method for the spatio-temporal heterogeneity of hydrogeological parameters of the present invention;

[0017] Figure 2 is a schematic diagram of the layered heterogeneous semi-pervious layer system in the embodiment of the present invention;

[0018] Figure 3 is a spatio-temporal distribution diagram of the hydraulic conductivity and storage rate of the layered heterogeneous semi-pervious layer in the embodiment of the present invention. DETAILED DESCRIPTION OF THE INVENTION

[0019] The following details the embodiments of the present invention, and the examples of the embodiments are shown in the drawings. The embodiments described below with reference to the drawings are exemplary and are only used to explain the present invention and should not be construed as limiting the present invention.

[0020] As Figure 1 shown, the present invention proposes a geostatistical inversion method for the spatio-temporal heterogeneity of hydrogeological parameters, including the following steps:

[0021] 1. Construct the column vector of the parameters to be determined: partition the geological model in the spatial domain, rely on the "time-space equivalence" to discretize the dynamic changes of the partition parameters in the time domain into column vectors, and assemble them into a comprehensive column vector of the parameters to be determined

[0022] First, the spatial domain is divided into N partitions, assuming that each partition is homogeneous and has p independent parameters to be determined. In the context of high-dimensional parameter inversion optimization, the dynamic changes of parameters in the time domain are equivalent to the one-dimensional heterogeneity of model parameters in the spatial domain, that is, the principle of "time-space equivalence". Therefore, the method of geostatistical inversion can be used to identify the dynamic change characteristics of model parameters. Specifically, for the i-th parameter s to be determined i , we can use a discrete column vector of sufficiently large size to approximate the dynamic change process of the parameters in the time domain, as follows:

[0023]

[0024] In the formula, M is the number of parameters to be determined in the time domain. According to formula (1), a column vector containing all the parameters to be determined in the study area can be constructed, as shown in formula (1). In this case, s is a column vector of MNp×1.

[0025] s=(s1,s2,…,s Np ) T (2)

[0026] 2. Assemble the prior covariance matrix: Based on the variogram, obtain the covariance matrix of each partition parameter in the time domain, and construct the prior covariance matrix of the parameters to be determined by block assembly.

[0027] For each element s in the parameter s to be found i , assuming that its distribution satisfies Gaussian continuous distribution, its variance function can be constructed using the exponential model:

[0028]

[0029] Where: var(h) represents the covariance corresponding to the distance h, represents the variance of the parameter, and d represents the correlation distance of the parameter, which means that when the distance h between two parameters exceeds d, the two parameters no longer have obvious correlation. In the time domain, the correlation length can also be expressed as the two parameters have obvious correlation within the time range d. According to the variogram, the correlation with parameter s can be calculated. i The corresponding covariance matrix Q i There is no correlation between different parameters, so the block assembly method should be adopted when assembling the covariance matrix, such as formula (4).

[0030]

[0031] 3. Dimensionality reduction of the covariance matrix: Based on the principal component analysis method, dimensionality reduction is performed on the prior covariance matrix of the parameters to be determined, and the terms related to the prior covariance matrix in the inversion model are simplified.

[0032] Since the covariance matrix Q corresponding to each parameter to be determined i is a symmetric positive definite matrix, the overall covariance matrix Q is a symmetric positive definite matrix. Therefore, the eigenvalue decomposition method can be used to decompose Q according to formula (5).

[0033]

[0034] In the formula, λ is the eigenvalue matrix, and V is the matrix composed of the eigenvectors corresponding to the eigenvalues. After the eigenvalue decomposition of the covariance matrix, the eigenvalues are sorted from large to small. The larger the eigenvalue, the greater the weight it occupies. Therefore, several eigenvalues and their corresponding eigenvectors can be taken as needed to achieve dimensionality reduction of the prior covariance matrix according to formula (5). According to formula (5), the resolution α is defined to characterize the weight of the obtained principal components. In actual use, it is usually assigned 0.90 or 0.95 according to the 90% or 95% confidence interval.

[0035]

[0036] 4. Approximate calculation of parameter sensitivity: Use the first-order Taylor expansion to approximately calculate the sensitivity matrix of the parameters to be determined, and simplify the calculation by using the dimension-reduced covariance matrix and parallel computing.

[0037] After the dimensionality reduction of the prior covariance, the operations of related parameters (H k Q, QH k and ) are also streamlined as follows:

[0038] H k Q = (H k U)U T (7)

[0039]

[0040] H k U can be calculated by using the first-order Taylor expansion:

[0041]

[0042] In the formula, U (j) represents the j-th column of the matrix U, δ represents a very small change. For double-precision calculations, generally take δ ≈ 10 -8 ; ||U (j) || and ||sk respectively represent U (j) and s k 's modulus; F(·) represents the forward model determined by the research problem and the parameter to be solved. Similarly, H k X i and H k s k can be calculated according to the formula:

[0043]

[0044] By calculating H k U, H k Q, QH k and can be indirectly calculated, thus completing the sensitivity operation of the parameter to be solved. When programming and operating using MATLAB, the parallel operation method can be adopted to speed up the calculation time-consuming in each iteration process.

[0045] 5. Solve the co-Kriging equation: Based on the calculation results of the parameter sensitivity analysis, use the co-Kriging equation to solve the necessary parameters required for the iterative solution of the parameter

[0046] In the geostatistical inversion method, assume that the estimated value of s after k iterations is s k , then the temporary estimated value after the (k + 1)-th iteration is:

[0047] s k+1 = Xβ k+1 +(H k Q) T ξ k+1 (13)

[0048]

[0049] In the formula: the diagonal element X i of X is a column vector of M*1, the elements of X i are all 1, H k is the sensitivity matrix of F(s k ) to s k , β k+1 and ξ k+1 are respectively the parameter mean and the error term after the (k + 1)-th iteration, and can be solved according to the following co-Kriging equation:

[0050]

[0051] In the formula: θ k is a dynamic stability multiplier, is a diagonal matrix composed of the diagonal elements of ; H k Xi , H k s k and can be calculated according to formulas (11), (12) and (9).

[0052] 6. Iterative solution: Construct the objective function, solve for the parameters according to the iterative calculation formula, substitute the parameter estimated values into the objective function, compare the objective function values in the previous and subsequent iterations until the number of iterations reaches the upper limit or the difference between the objective function values is within the tolerance range

[0053] The objective function of the inversion model is given by formula (15):

[0054] L k = [F obs - F(s k )] T R -1 [F obs - F(s k )] + (s k - Xβ k ) T Q -1 (s k - Xβ k ) (16)

[0055] where: F obs is a column vector composed of observed values, F(s) is a column vector composed of simulated values obtained by substituting the parameter s to be solved into the forward model F(·); R is a covariance matrix used to characterize the observation error; β k is the parameter mean value of the k-th iteration; Q is the covariance matrix of the parameter to be solved.

[0056] When the number of iterations reaches the preset upper limit, or the difference between the values of the objective function L in the previous and subsequent iterations is less than the tolerance range, the iteration terminates, and the s k at the last iteration is used as the final estimated value of the parameter s to be solved, and the posterior covariance matrix of the parameter s is calculated according to the formula:

[0057]

[0058] where: The diagonal elements of Q post are the variance values of the elements in the parameter s to be solved. According to the Gaussian normality assumption, the 95% confidence interval of the parameter to be solved can be calculated as [-1.96σ(s) 1.96σ(s)], where σ represents the standard deviation of s.

[0059] The following shows through a specific case the relationship between the hydraulic conductivity (K) and the storage rate (S s) calculation result. Consider a layered heterogeneous semi-pervious layer in the water release deformation and recovery stage. According to the lithological characteristics, the semi-pervious layer can be divided into three sub-layers (as shown in Figure 2 ), since each sub-layer has two parameters, K and S s , which need to be inverted. Therefore, the total number of parameters to be determined in the spatial domain is Np = 6. In the time domain, the dynamic changes of each parameter are discretized into 120 nodes, that is, M = 120. Therefore, the total number of parameters to be determined is MNp = 720. The forward model is given by the water release deformation model of the semi-pervious layer and is calculated using the finite element method. The deformation data of the semi-pervious layer measured by the hierarchical gauge is used as the observation data. In the inversion model, the prior parameters of the parameters to be determined are set as shown in Table 1. According to Table 1, the prior covariance matrix of the parameters to be determined can be constructed respectively and assembled into Q. Principal component analysis is used, and the resolution α is taken as 0.95 to reduce the dimension of the prior covariance matrix. The parameters to be determined are iteratively solved, the maximum number of iterations is set to 300, and the tolerance is 0.01. In each iteration, the co-Kriging equation is used to solve β k and ξ k (1 < k < 300) in the iterative formula. During the solution process, the operations related to the sensitivity matrix H k of the parameters to be determined are approximated by the first-order Taylor expansion. The objective function is constructed, and the parameter estimation values after the current iteration, β k and ξ k are substituted into the objective function to calculate the objective function value L k of the current iteration. The objective function values of the two consecutive iterations are compared. When the difference between the two values is less than the tolerance of 0.01, or the number of iterations reaches 300 times, the iteration terminates. Calculate the posterior covariance matrix of the parameter estimation values of the parameters to be determined at this time, take the diagonal elements of the posterior covariance matrix, and perform corresponding operations on them as the 95% confidence intervals of each parameter.

[0060] Table 1 Parameter settings of the inversion model

[0061]

[0062] According to the above parameter settings and operation procedures, the spatio-temporal distribution characteristics of the hydraulic conductivity and storage rate of the three sub-layers of the semi-pervious layer are calculated as shown in Figure 3As shown in the figure. It can be seen that the variation trend and amplitude of the permeability coefficient and water storage rate of the third sub-layer are the largest on the time scale. The permeability coefficient of the third sub-layer began to increase after 1994, while the water storage rate began to decrease around 1997. The permeability coefficients and water storage rates of the first and second sub-layers remain stable fluctuations on the time scale; however, the water storage rate of the second sub-layer is slightly higher than that of the first sub-layer. The above results well reflect the lithological differences of different sub-layers of the layered heterogeneous aquitard, as well as the characteristics of different sensitivities to water level recovery. At the same time, it proves the ability of the method of the present invention to characterize the spatio-temporal heterogeneity of parameters, which has important practical significance for characterizing the spatio-temporal heterogeneity of parameters in a more complex spatio-temporal evolution process.

[0063] Based on the same inventive concept, an embodiment of the present application provides a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor. When the processor executes the computer program, the steps of the aforementioned geostatistical inversion method for spatio-temporal heterogeneity of hydrogeological parameters are implemented.

[0064] Based on the same inventive concept, an embodiment of the present application provides a computer-readable storage medium storing a computer program, and when the computer program is executed by a processor, the steps of the aforementioned geostatistical inversion method for spatio-temporal heterogeneity of hydrogeological parameters are implemented.

[0065] Those skilled in the art should understand that the embodiments of the present invention can be provided as a method, a system, or a computer program product. Therefore, the present invention can take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware aspects. Moreover, the present invention can take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0066] The present invention is described with reference to the flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. It should be understood that each flow and / or block in the flowchart and / or block diagram, as well as the combination of flows and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to the processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing devices to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing devices generate means for implementing the functions specified in Figure 1 one process or multiple processes and / or blocks Figure 1 one block or multiple blocks.

[0067] These computer program instructions can also be stored in a computer-readable memory that can direct a computer or other programmable data processing apparatus to operate in a particular manner, such that the instructions stored in the computer-readable memory produce a manufacture including an instruction device that implements the functions specified in one or more of the processes and / or blocks Figure 1 in one or more of the processes and / or blocks Figure 1 specified in the block or blocks.

[0068] These computer program instructions can also be loaded onto a computer or other programmable data processing apparatus to cause a series of operational steps to be performed on the computer or other programmable apparatus to produce a computer-implemented process, whereby the instructions executed on the computer or other programmable apparatus provide steps for implementing the functions specified in one or more of the processes and / or blocks Figure 1 in one or more of the processes and / or blocks Figure 1 specified in the block or blocks.

[0069] The above embodiments are only for illustrating the technical idea of the present invention, and the protection scope of the present invention cannot be limited thereby. Any modifications made on the basis of the technical solution according to the technical idea proposed by the present invention shall fall within the protection scope of the present invention.

Claims

1. A geostatistical inversion method for spatiotemporal heterogeneity of hydrogeological parameters, characterized in that: The steps include: Step 1: partition the target geological entity in the spatial domain based on the geological model of the target geological entity. The hydrogeological parameters to be determined in each partition are the same. Based on the principle of spatiotemporal equivalence, the dynamic changes of the hydrogeological parameters to be determined in each partition in the time domain are discretized into column vectors, thereby constructing column vectors of the hydrogeological parameters to be determined in all partitions. Step 2: for the desired hydrogeological parameters of each sub-area, construct its variogram in the time domain, calculate the corresponding covariance matrix according to the variogram, and construct the prior covariance matrix of the desired hydrogeological parameters of all sub-areas by block assembly; Step 3, reducing the dimension of the prior covariance matrix based on the principal component analysis method to obtain the prior covariance matrix after dimension reduction; Step 4, setting the initial value of the hydrogeological parameter to be determined, and calculating the estimated value of the hydrogeological parameter to be determined by an iterative method based on the prior covariance matrix after dimensionality reduction; Step 5, constructing the objective function of the inversion of the hydrogeological parameter to be determined, substituting the estimated value of the hydrogeological parameter to be determined obtained in each iteration into the objective function, comparing the objective function values ​​of the two iterations, and terminating the iteration when the preset number of iterations is reached or the objective function values ​​of the two iterations are within the tolerance range, and taking the estimated value obtained in the last iteration as the final estimated value; Step 6: After the iteration is terminated, the posterior covariance matrix of the desired hydrogeological parameter is calculated to obtain the value range of the desired hydrogeological parameter.

2. The method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters according to claim 1, characterized in that: The specific process of step 1 is as follows: Based on the geological model, the target geological entity is divided into N partitions in the spatial domain, and each partition has p hydrogeological parameters to be determined. For the i-th parameter, a discrete column vector is used to represent the dynamic change process of the parameter in the time domain. The formula is as follows: In the formula, s i Represents the discrete column vector of the dynamic change of the i-th parameter in the time domain, Both represent s i , M is the length of the column vector in the time domain, and Np is the total number of hydrogeological parameters to be determined; The column vector of the hydrogeological parameters to be determined for all sub-regions is expressed as: s=(s1,s2,…,s Np ) T In the formula, s represents the column vector of the hydrogeological parameters to be determined in all partitions, and s is a column vector of MNp×1.

3. The method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters according to claim 2, characterized in that: In step 1, the target geological entity is a layered heterogeneous aquitard or other geological entity with spatial heterogeneity, and the hydrogeological parameters to be determined for each partition include but are not limited to permeability coefficient and water storage rate.

4. The method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters according to claim 2, characterized in that: The specific process of step 2 is as follows: For s i Construct the variogram: In the formula, var(h) represents s i The corresponding covariance when the distance between any two components in is h is, Indicates i The variance of , d represents the preset correlation distance; According to the variogram calculation and parameter s i The corresponding covariance matrix Q i , the prior covariance matrix Q of all the hydrogeological parameters to be determined in all partitions is constructed by block assembly as follows: In the formula, Q1, Q2, Q Np They respectively represent the covariance matrices corresponding to the first, second and Npth hydrogeological parameters to be determined, and Q is a symmetric positive definite matrix.

5. The method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters according to claim 2, characterized in that: In step 3, the prior covariance matrix is ​​subjected to eigenvalue decomposition to obtain eigenvalues ​​and corresponding eigenvectors, and the eigenvalues ​​are sorted in descending order; The resolution α is preset to characterize the weight of the selected eigenvalues, and the first r eigenvalues ​​are selected according to the following formula to reduce the dimension of the prior covariance matrix Q: In the formula, λ j is the eigenvalue, J is the number of eigenvalues, and r is the number of eigenvalues ​​selected according to the resolution; Prior covariance matrix Q after dimensionality reduction r It is expressed as: Q r =UU T Where λ is the matrix composed of eigenvalues ​​selected according to the resolution, and V is the matrix composed of eigenvectors corresponding to the eigenvalues ​​selected according to the resolution.

6. The method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters according to claim 5, characterized in that: In step 4, assume that the estimated value of s after k iterations is s k , then the estimated value s after the k+1th iteration k+1 for: s k+1 =Xβ k+1 +(H k Q) T ξ k+1 Where X is a matrix, and the diagonal elements of X are i is an M*1 column vector, X i The elements of H are all 1, 1≤i≤Np, k Q=(H k U)U T , H k F(s k ) k The sensitivity matrix of the target geological entity is F(·), β k+1 (Np×1) and ξ k+1 are the parameter mean and error term of the k+1th iteration, respectively, and are solved according to the co-kriging equation: In the formula, θ k is the dynamic stability multiplier, R is the covariance matrix used to characterize the observation error, Is The diagonal matrix composed of the diagonal elements of F obs is a column vector composed of observation values ​​representing the characteristics of the target geological entity; H k U, H k X, H k s k The calculation is performed using the first-order Taylor expansion method: Where U (j) represents the j-th column of matrix U, ||U (j) || means U (j) The modulus, j = 1, ..., r, δ represents the amount of change, ||s k || means s k The modulus of ||X i || represents X i Model.

7. The method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters according to claim 2, characterized in that: In step 5, the objective function of the inversion of the hydrogeological parameters to be determined is as follows: L k =[F obs -F(s k )] T R -1 [F obs -F(s k )]+(s k -Xβ k ) T Q -1 (s k -Xβ k ) Where, L k represents the objective function value after the kth iteration, F obs is a column vector of observations, F(·) is the forward model of the target geological entity, s k is the estimated value after the kth iteration, R is the covariance matrix used to characterize the observation error, Q is the prior covariance matrix of the hydrogeological parameters to be determined in all partitions, β k is the parameter mean of the kth iteration.

8. The method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters according to claim 2, characterized in that: In step 6, the posterior covariance matrix of the hydrogeological parameters to be determined is as follows: In the formula, Q post represents the posterior covariance matrix of the hydrogeological parameters to be determined, Q is the prior covariance matrix of the hydrogeological parameters to be determined in all partitions, and H k F(s k ) k The sensitivity matrix of the target geological entity is s. k is the estimated value after the kth iteration, X is a column vector with all elements set to 1, and Q post The diagonal elements of are the variance values ​​of each element in s. The 95% confidence interval of the hydrogeological parameter to be determined is [-1.96σ(s) 1.96σ(s)], where σ represents the standard deviation of s.

9. A computer device comprising a memory, a processor, and a computer program stored in the memory and capable of running on the processor, characterized in that: When the processor executes the computer program, the steps of the method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters as described in any one of claims 1 to 8 are implemented.

10. A computer-readable storage medium storing a computer program, characterized in that: When the computer program is executed by a processor, the steps of the method for geostatistical inversion of spatiotemporal heterogeneity of hydrogeological parameters as described in any one of claims 1 to 8 are implemented.

Images