CN115758820A - Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment - Google Patents

Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment Download PDF

Info

Publication number
CN115758820A
CN115758820A CN202211424242.4A CN202211424242A CN115758820A CN 115758820 A CN115758820 A CN 115758820A CN 202211424242 A CN202211424242 A CN 202211424242A CN 115758820 A CN115758820 A CN 115758820A
Authority
CN
China
Prior art keywords
temperature
transient
matrix
thermal conductivity
prediction model
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202211424242.4A
Other languages
Chinese (zh)
Inventor
杨健
胡祥云
黄国疏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China University of Geosciences
Original Assignee
China University of Geosciences
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China University of Geosciences filed Critical China University of Geosciences
Priority to CN202211424242.4A priority Critical patent/CN115758820A/en
Publication of CN115758820A publication Critical patent/CN115758820A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02EREDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E10/00Energy generation through renewable energy sources
    • Y02E10/10Geothermal energy

Landscapes

  • Geophysics And Detection Of Objects (AREA)

Abstract

The invention provides a three-dimensional transient prediction method, a device and electronic equipment for a heat conductivity coefficient of a deep stratum, wherein the method uses an underground heat flow value, thermal physical property parameters of each stratum, a fluid migration rule and the like as prior information, carries out transient finite element temperature numerical simulation and constructs a regularization objective function, constructs a Hessian matrix and a gradient vector based on a Gauss-Newton algorithm, solves the product of the Jacobian matrix and an arbitrary vector by utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule, avoids the solution and storage of a large dense Jacobian matrix and obtains a model correction amount, searches for an optimal model step length to update a model parameter by using a Wolfe criterion, circularly predicts so that a data fitting difference is smaller than a preset fitting difference, and outputs a current prediction model parameter as an optimal heat conductivity coefficient prediction result. The method can intuitively characterize the heat conductivity coefficient distribution characteristic of the deep stratum medium in the time domain, and has high prediction precision and strong practicability.

Description

Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment
Technical Field
The invention relates to the field of thermal conductivity prediction, in particular to a deep stratum thermal conductivity three-dimensional transient prediction method and device and electronic equipment.
Background
Geothermal resources are attracting more and more attention as a new type of energy source because of their unique advantages, such as low cost, sustainable utilization, environmental protection, etc., which are incomparable with other energy sources. The development and utilization of geothermal resources are greatly promoted, the energy structure is improved, and the method has important significance for solving the increasingly serious global environment problem and realizing the double-carbon target. The thermal conductivity coefficient, the thermal diffusion coefficient, the specific heat capacity and the like are very important rock thermophysical parameters in a deep stratum, and are the basis for quantitatively evaluating temperature change, and accurate calculation of temperature is very important for clarifying an underground heat transfer mode, a geothermal causation mechanism and delineating a target area of geothermal resources.
The thermal conductivity is one of the most important parameters in the rock thermophysical parameters, and rock properties are greatly different due to different geological structures and evolution processes in the earth, so that the thermal conductivity of the rock shows obvious heterogeneity in spatial distribution, and the regional thermal characteristics are changed. Currently, predictions of rock thermal conductivity are typically characterized by laboratory measurements on drill cuttings or core samples or by geophysical logs. However, the methods only represent the steady-state rock heat conductivity coefficients of point and line scales, neglect the influence of the flow of underground hot fluid on the rock heat conductivity coefficients, and limit the accurate prediction of the rock heat conductivity coefficients of regional time domains, which has strong limitations on the research in the fields of high-temperature geothermal target area site selection, thermal cause analysis and the like of a research area.
Disclosure of Invention
The invention provides a method, a device and electronic equipment for predicting the three-dimensional transient state of the heat conductivity of a deep stratum, and aims to solve the technical problems that the existing rock heat conductivity coefficient prediction method can only represent the steady-state rock heat conductivity of a point and a line scale and neglects the influence of the flow of underground hot fluid on the rock heat conductivity. The method comprises the steps of using an underground heat flow value, each stratum thermophysical property parameter, a fluid migration rule and the like as prior information, carrying out transient finite element temperature numerical simulation and constructing a regularization target function, constructing a Hessian matrix and a gradient vector based on a Gauss-Newton algorithm, solving the product of the Jacobian matrix and an arbitrary vector by using a Jacobian-free Krylov subspace technology based on a time chain derivation rule, avoiding solving and storing the large dense Jacobian matrix and obtaining model correction quantity, searching an optimal model step length by using a Wolfe criterion to update model parameters, circularly predicting to enable the data fitting difference to be smaller than a preset fitting difference, outputting a current prediction model parameter as an optimal heat conductivity coefficient prediction result, and realizing three-dimensional transient prediction of the heat conductivity coefficient of the deep medium.
According to a first aspect of the invention, a deep formation thermal conductivity three-dimensional transient prediction method is provided, which comprises the following steps:
acquiring an underground water migration rule in a uniform half-space research area and a response mechanism thereof to a heat conductivity coefficient, constructing a heat conductivity coefficient abnormal body in the research area, giving boundary conditions, carrying out transient finite element temperature numerical simulation by combining the distribution of thermophysical parameters in the research area and the boundary conditions, and acquiring an underground space temperature field T (x, y, z) = d at a given moment obs X, y and z respectively represent the directions of x, y and z axes;
according to the temperature field d of the underground space obs Establishing an objective function of transient temperature forward modeling by utilizing Tikhonov regularization, acquiring a linear equation set consisting of a Hessian matrix and a gradient vector by adopting a Gauss-Newton algorithm, solving the Jacobian matrix of a first-order derivative of temperature numerical simulation response to a prediction model parameter in the Hessian matrix and the product of the transposition of the Jacobian matrix and an arbitrary vector by utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule, and solving to obtain a correction quantity delta m of the prediction model;
acquiring the optimal prediction model search step length meeting the Wolfe criterion according to the correction quantity delta m of the prediction model, and updating the parameters of the prediction model;
carrying out transient finite element temperature numerical simulation according to the updated prediction model parameters to obtain an underground space temperature field under the model parameters, calculating a data fitting difference between the model parameters and actual measurement data, and if the condition is smaller than a preset fitting difference condition, determining the current prediction model parameters as an optimal heat conductivity coefficient prediction result; otherwise, returning to the prediction method and solving a new prediction model correction quantity delta m.
Preferably, the groundwater migration law and the response mechanism thereof to the thermal conductivity include:
setting the migration rates of the movement of underground water in the directions x, y and z of a research area as v x 、v y 、v z And the fluid is transported and immersed into the surrounding medium, so that the heat conductivity coefficient of the medium is changed, and the following conditions are met:
k=(1-φ)·k s +φ·k f
wherein, k s 、k f The real thermal conductivity, the medium thermal conductivity and the fluid thermal conductivity are respectively, and the unit is W/m.K; phi is the porosity of the medium and is dimensionless.
Preferably, the boundary condition satisfies:
Figure BDA0003941296360000031
wherein, alpha and beta are constants, alpha =1, beta = T up The time is the upper boundary condition; α =0, β = q down The lower boundary condition is adopted; α =0 and β =0 is a four-week boundary condition. T is temperature, in degrees Celsius up Is the surface temperature in units of DEG C, n is the normal direction of the boundary, K is the thermal conductivity in units of W/m.K, q down Is a heat flow value in W/m 2
Preferably, the transient finite element temperature numerical simulation adopts a hexahedral structure discrete mode, and the temperature T, the heat conductivity k and the fluid heat capacity mu in the hexahedral unit f Medium heat capacity [ mu ] s The fluid transport velocity v and the heat source Q all adopt linear interpolation:
Figure BDA0003941296360000032
Figure BDA0003941296360000033
wherein i is a node label, T i 、k i 、μ fi 、μ si 、v i 、Q i The temperature, the heat conductivity coefficient, the fluid heat capacity, the medium heat capacity, the fluid flow velocity and the heat source of each node in the unit are respectively, and the unit is W/m. K, J/(m) 3 ·K)、J/(m 3 ·K)、m/s、W/m 3 ;N i Is a shape function, satisfies:
Figure BDA0003941296360000034
wherein ξ i 、η i 、ζ i The coordinate of the node i on the parent unit is dimensionless; xi, eta and zeta in the parent unit and coordinates x, y and z in the child unit are in a relational formula:
Figure BDA0003941296360000035
wherein x is 0 、y 0 、z 0 Respectively being the middle points of the subunits in the directions of x, y and z axes; a. b and c are the side lengths of the subunits in the directions of x, y and z respectively, and the unit is m.
Preferably, the step of performing a transient finite element temperature numerical simulation in combination with the distribution of the thermophysical parameters in the investigation region and the boundary conditions to obtain the temperature field of the subsurface space at a given time comprises:
the differential equation for a given transient heat transfer subsurface space temperature field is:
Figure BDA0003941296360000041
wherein K is the thermal conductivity coefficient, the unit W/m.K, T is the temperature, the unit degree centigrade, mu f The unit J/(m) is the heat capacity of the fluid 3 ·K),μ s The unit J/(m) is the heat capacity of the medium 3 K), v is the fluid flow rate, m/s, Q is the heat source, W/m 3 T is the time, in units of s,
Figure BDA0003941296360000042
in order to be the divergence degree,
Figure BDA0003941296360000043
is a gradient.
Multiplying and integrating the differential equation of the temperature field of the underground space by an impulse function delta T:
Figure BDA0003941296360000044
wherein, omega is the research area, gamma is the boundary, and delta T is the impulse function.
According to the Hamiltonian operation rule in the field theory and On Gao Gongshi, and substituting the boundary condition, the first term in the above formula is rewritten as:
Figure BDA0003941296360000045
therefore, the integral equation for the transient heat transfer subsurface space temperature field is:
Figure BDA0003941296360000046
carrying out discrete subdivision and linear interpolation on the integral equation of the underground space temperature field, and solving to obtain the unit integral of each item in the integral equation:
Figure BDA0003941296360000047
Figure BDA0003941296360000048
Figure BDA0003941296360000049
Figure BDA00039412963600000410
Figure BDA00039412963600000411
Figure BDA00039412963600000412
wherein T is time, unit s, delta T e For arrays of impulse functions within a cell, T e For temperature arrays within the cell, K 1e ,K 2e ,K 3e ,G e Is a stiffness matrix within the cell, P 1e ,P 2e For the source vector in the unit, solving to obtain a rigidity matrix and a source vector internal coefficient:
Figure BDA0003941296360000051
Figure BDA0003941296360000052
Figure BDA0003941296360000053
Figure BDA0003941296360000054
Figure BDA0003941296360000055
Figure BDA0003941296360000056
Figure BDA0003941296360000057
wherein k is 1ij ,k 2ij ,k 3ij ,g ij ,p 1ij ,p 2ij Respectively a stiffness matrix K 1e ,K 2e ,K 3e ,G e Source vector P 1e ,P 2e The elements in the table, i, j, l are node labels, a, b, c are the side length of the subunit in x, y, z direction, unit is m, xi i 、η i 、ζ i Is the coordinate of node i on the parent cell.
And expanding the number of each node of each unit to the element position corresponding to the research area to obtain a linear equation set of the synthesized total stiffness matrix and the source vector:
Figure BDA0003941296360000058
where δ T is the impulse function array of the study region, T is the temperature array within the study region, K 1 ,K 2 ,K 3 G is respectively a stiffness matrix in the investigation region, P 1 ,P 2 Omitting δ T for source vectors within the region of interest T After simplification, the method comprises the following steps:
Figure BDA0003941296360000061
wherein K = K 1 +K 2 +K 3 ,P=P 1 +P 2
Differentiating the time T in the transient underground heat transfer underground space temperature field, and after delta T, the temperature field is changed from T 1 Change to T 2 Then, then
Figure BDA0003941296360000062
Additionally, the interpolated temperature within the time of Δ t
Figure BDA0003941296360000063
Instead of T in the above equation, the linear equation set for the transient temperature field is:
Figure BDA0003941296360000064
by solving the above formula, after finite element numerical simulation is realized, the temperature field T (x, y, z) = d of the underground space at a given moment is obtained obs And x, y and z represent x, y and z axis directions, respectively.
Preferably, the objective function of forward modeling of transient temperature constructed by the Tikhonov regularization is as follows:
Φ=||D[d(m)-d obs ]|| 2 +λ||W(m-m apr )|| 2
wherein D is a data weighting matrix; d obs Is observed data, namely a temperature field T (x, y, z) of the underground space; m is a parameter vector of the current prediction model; d (m) is simulation data under the current prediction model parameters; w is a smooth constraint matrix; λ is a damping factor; m is apr Is predictive model prior information;
based on a linear equation set of a total rigidity matrix and a source vector in a finite element, the heat conductivity coefficient k in the linear equation set is derived to obtain:
Figure BDA0003941296360000065
Figure BDA0003941296360000066
Figure BDA0003941296360000067
where K, G is the stiffness matrix in the region of interest, P 1 ,P 2 For the source vector in the region of interest, Δ T is the time interval, and the units s, K are the thermal conductivity array, and the units W/m.K, T i 、T i+1 Temperature arrays at time i and i +1, respectively. At the initial time (t = 0), the above equation is rewritten as:
Figure BDA0003941296360000068
wherein, T 0 Is an array of temperatures at the initial time. In order to ensure that the temperature value and the heat conductivity coefficient are not negative, the data parameters and the model parameters are set as the logarithm of the temperature and the logarithm of the heat conductivity coefficient in the prediction process: d = lnT, m = lnk; according to the chain derivation rule, the elements in the Jacobian matrix are:
Figure BDA0003941296360000071
wherein, i =1,2, …, N, N is the number of data parameters, j =1,2, …, M, M is the number of model parameters, d i For the data parameter of the ith node in the investigation region, m j Model parameter, T, representing the jth node in the study area i Denotes the temperature, k, of the ith node in the investigation region j And representing the heat conductivity coefficient of the jth node in the research area, and writing the Jacobian matrix as follows:
J=RGR'
wherein R and R' are diagonal matrices,
Figure BDA0003941296360000072
R'=diag(k 1 ,k 2 ,…,k M ) G is a matrix composed of elements G (i, j);
for an arbitrary vector consisting of dimensionless elements v = (v) 1 ,v 2 ,…,v M ) T And the conversion formula based on the finite element rigidity matrix is as follows:
Figure BDA0003941296360000073
k (v) is a rigidity matrix when the heat conductivity coefficient array is v; and when t =0, there are
Figure BDA0003941296360000081
Wherein e is 1 ,e 2 ,...,e M Since all vectors are unit vectors, the product of Jacobian matrix J and arbitrary vector v can be obtained.
Similarly, the transposition of the Jacobian matrix and an arbitrary vector y = (y) composed of dimensionless elements are calculated 1 ,y 2 ,…,y N ) T The product of (a);
according to a Gauss-Newton algorithm, neglecting a second-order information item, obtaining the following linear equation system to solve the model correction quantity delta m:
HΔm=-g
wherein H is Hessian matrix, H = (J) T D T DJ+λW T W), g is gradient vector, g = [ J = T D T D(d-d obs )+λW T W(m-m apr )]。
Preferably, the optimal prediction model search step length satisfying the Wolfe criterion is obtained according to the prediction model correction amount Δ m, and the following inequality is satisfied:
Figure BDA0003941296360000082
wherein alpha is a search step length and is dimensionless, and c is a dimensionless number; in the process of obtaining the search step length, the initial value of alpha is 1, and if the inequality is not satisfied, the initial value of alpha is made to be 1
α=ρα
Wherein rho is a decimal larger than 0 and smaller than 1, and the process is circulated until the Wolfe criterion is met; thus, the update formula for the model parameters is:
m k+1 =m k +αΔm k
wherein m is k 、m k+1 Model parameters, Δ m, at kth and (k + 1) th iterations, respectively k The model correction quantity obtained in the k iteration is obtained.
Preferably, the calculation of the data fitting difference is as follows:
Figure BDA0003941296360000091
wherein, the RMS is root mean square and is used for measuring the deviation between the measured value and the predicted value, and the closer to 0, the smaller the deviation between the measured value and the predicted value is, the more reliable the prediction result is; d (m) i 、d obs i The ith predicted value and the ith measured value are respectively, n is the number of observed data, and the dimension is not existed;
if the data fitting difference is smaller than the preset fitting difference, stopping iteration and outputting the current prediction model parameters as the optimal prediction result; otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
According to a second aspect of the invention, the invention provides a deep formation thermal conductivity three-dimensional transient prediction device, which comprises the following modules:
the transient finite element numerical simulation module is used for acquiring the underground water migration rule in a uniform half-space research area and the response mechanism of the underground water migration rule to the heat conductivity coefficient, constructing a heat conductivity coefficient abnormal body in the research area, setting boundary conditions, carrying out transient finite element temperature numerical simulation by combining the distribution of the thermophysical parameters in the research area and the boundary conditions, and acquiring the underground space temperature field T (x, y, z) = d at a given moment obs X, y and z respectively represent the directions of x, y and z axes;
a model correction amount calculation module for calculating the correction amount according to the temperature field d of the underground space obs The method comprises the steps of utilizing Tikhonov regularization to construct an objective function of transient temperature forward modeling, adopting a Gauss-Newton algorithm to obtain a linear equation set consisting of a Hessian matrix and gradient vectors, and utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule to solve the Hessian matrixThe temperature numerical simulation responses to the product of the Jacobian matrix of the first derivative of the prediction model parameters and the transpose of the Jacobian matrix and any vector, and the correction quantity delta m of the prediction model is obtained through solving;
the model parameter updating module is used for acquiring the optimal prediction model searching step length meeting the Wolfe criterion according to the correction quantity delta m of the prediction model and updating the parameters of the prediction model;
the thermal conductivity coefficient prediction module is used for carrying out transient finite element temperature numerical simulation according to the updated prediction model parameters to obtain an underground space temperature field under the model parameters, calculating the data fitting difference between the data and the actually measured data, and outputting the current prediction model parameters as the optimal thermal conductivity coefficient prediction result if the current prediction model parameters are smaller than the preset fitting difference condition; otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
According to a third aspect of the invention, there is provided an electronic device comprising a memory, a processor and a computer program stored on the memory and executable on the processor, the processor implementing the steps of the deep formation thermal conductivity three-dimensional transient prediction method when executing the program.
The technical scheme provided by the invention has the following beneficial effects:
the invention provides a three-dimensional transient prediction method, a three-dimensional transient prediction device and electronic equipment for the thermal conductivity of a deep stratum, which are combined with thermophysical parameters and a ground water migration rule, realize accurate and rapid prediction of the thermal conductivity of a deep medium in a time domain by using a Jacobian-free Krylov subspace technology based on a time chain derivation rule, and have the advantages of strong practicability, wide and deep prediction range and result accuracy of 88.97-94.19%.
Drawings
The invention will be further described with reference to the accompanying drawings and examples, in which:
FIG. 1 is a flow chart of a deep formation thermal conductivity three-dimensional transient prediction method of the present invention;
FIG. 2 is a plot of groundwater migration rates for a study area according to the present invention;
fig. 3 is a diagram of a distribution of a real thermal conductivity model according to the present invention, in which fig. 3 (a) is a diagram of a thermal conductivity distribution of an XZ direction cross section (y =2.5 km), fig. 3 (b) is a diagram of a thermal conductivity distribution of an XY direction cross section (z =1.1 km), and fig. 3 (c) is a diagram of a thermal conductivity distribution of an XY direction cross section (z =3.0 km);
FIG. 4 is a graph of a heat source parameter profile of the present invention;
FIG. 5 is a temperature field profile of a subsurface space (y =2.5 km) of a study area of the present invention;
fig. 6 is a temperature field profile under a priori information (y =2.5 km) of the present invention;
fig. 7 is a temperature field profile of the present invention at optimal prediction results (y =2.5 km);
fig. 8 is a graph showing the best prediction result of the underground thermal conductivity of the study region according to the present invention, in which fig. 8 (a) is a thermal conductivity distribution graph of an XZ direction section (y =2.5 km), fig. 8 (b) is a thermal conductivity distribution graph of an XY direction section (z =1.1 km), and fig. 8 (c) is a thermal conductivity distribution graph of an XY direction section (z =3.0 km);
FIG. 9 is a graph of data fitting difference as a function of iteration number for the present invention;
fig. 10 is a block diagram of a deep formation thermal conductivity three-dimensional transient prediction apparatus according to the present invention.
Detailed Description
For a more clear understanding of the technical features, objects and effects of the present invention, embodiments of the present invention will now be described in detail with reference to the accompanying drawings.
Referring to fig. 1, fig. 1 is a flow chart of a deep formation thermal conductivity three-dimensional transient prediction method according to the present invention; the method comprises the following steps:
s1: acquiring an underground water migration rule in a uniform half-space research area and a response mechanism thereof to a heat conductivity coefficient, constructing a heat conductivity coefficient abnormal body in the research area, giving boundary conditions, carrying out transient finite element temperature numerical simulation by combining the distribution of thermophysical parameters in the research area and the boundary conditions, and acquiring an underground space temperature field T (x, y, z) = d at a given moment obs X, y and z respectively represent the directions of x, y and z axes;
in the present embodimentThe step S1 specifically includes: setting the migration rate of groundwater movement in the x, y and z directions of the research area (refer to figure 2), and the heat conductivity coefficient k of the fluid f =0.622, the fluid migration dips into the surrounding medium, causing the medium thermal conductivity to change, satisfying:
k=(1-φ)·k s +φ·k f
wherein k, k s 、k f The real thermal conductivity, the medium thermal conductivity and the fluid thermal conductivity are respectively, and the unit is W/m.K; phi is the porosity of the medium and is dimensionless.
Constructing a thermal conductivity anomaly in the uniform half-space study area, wherein the real thermal conductivity model distribution is shown in fig. 3, wherein fig. 3 (a) is a thermal conductivity distribution diagram of an XZ direction section (y =2.5 km), fig. 3 (b) is a thermal conductivity distribution diagram of an XY direction section (z =1.1 km), and fig. 3 (c) is a thermal conductivity distribution diagram of an XY direction section (z =3.0 km);
referring to FIG. 4, given a heat source parameter profile, the heat capacity μ of the fluid f =4.2×10 6 J/(m 3 K) heat capacity of the medium μ s =2×10 6 J/(m 3 K) and given the boundary conditions:
Figure BDA0003941296360000111
wherein α =1, β =15 ℃ is the upper boundary condition; α =0, β =0.05W/m 2 The lower boundary condition is adopted; α =0 and β =0 is a four-week boundary condition. T is temperature, unit ℃, n is boundary normal direction, and is dimensionless, K is heat conductivity coefficient, unit W/m.K.
Performing transient finite temperature numerical simulation under the distribution of the thermophysical parameters and the boundary conditions to obtain a temperature of 5 × 10 -5 Underground space temperature field T (x, y, z) = d after Ma obs
As an optional implementation mode, the transient finite element temperature numerical simulation adopts a discrete mode of a hexahedral structure, and the temperature T, the heat conductivity k and the fluid heat capacity mu in a hexahedral unit f Medium heat capacity μ s Fluid transport velocity v, heatThe source Q adopts linear interpolation:
Figure BDA0003941296360000121
Figure BDA0003941296360000122
wherein i is a node label, T i 、k i 、μ fi 、μ si 、v i 、Q i The temperature, the heat conductivity coefficient, the fluid heat capacity, the medium heat capacity, the fluid flow velocity and the heat source of each node in the unit are respectively, and the unit is W/m. K, J/(m) 3 ·K)、J/(m 3 ·K)、m/s、W/m 3 ;N i Is a shape function, satisfies:
Figure BDA0003941296360000123
wherein ξ i 、η i 、ζ i The coordinate of the node i on the parent unit is dimensionless; xi, eta and zeta in the parent unit and coordinates x, y and z in the child unit are in a relational formula:
Figure BDA0003941296360000124
wherein x is 0 、y 0 、z 0 Respectively being the middle points of the subunits in the directions of x, y and z axes; a. b and c are the side lengths of the subunits in the directions of x, y and z respectively, and the unit is m.
The differential equation for a given transient heat transfer subsurface space temperature field is:
Figure BDA0003941296360000125
wherein K is the thermal conductivity coefficient, the unit W/m.K, T is the temperature, the unit degree centigrade, mu f Is the heat capacity of the fluidUnit J/(m) 3 ·K),μ s The unit J/(m) is the heat capacity of the medium 3 K), v is the fluid flow rate, m/s, Q is the heat source, W/m 3 T is the sum of time, unit s,
Figure BDA0003941296360000126
in order to be the divergence degree,
Figure BDA0003941296360000127
is a gradient.
Multiplying and integrating the differential equation of the temperature field of the underground space by an impulse function delta T:
Figure BDA0003941296360000128
wherein, omega is the research area, gamma is the boundary, and delta T is the impulse function.
According to the Hamiltonian operation rule in the field theory and On Gao Gongshi, and substituting the boundary condition, the first term in the above formula is rewritten as:
Figure BDA0003941296360000129
therefore, the integral equation for the transient heat transfer subsurface space temperature field is:
Figure BDA0003941296360000131
carrying out discrete subdivision and linear interpolation on the integral equation of the underground space temperature field, and solving to obtain the unit integral of each item in the integral equation:
Figure BDA0003941296360000132
Figure BDA0003941296360000133
Figure BDA0003941296360000134
Figure BDA0003941296360000135
Figure BDA0003941296360000136
Figure BDA0003941296360000137
wherein T is time, unit s, delta T e As an array of impulse functions within a cell, T e For temperature arrays within the cell, K 1e ,K 2e ,K 3e ,G e Respectively a stiffness matrix within the cell, P 1e ,P 2e For the source vector in the unit, solving to obtain a rigidity matrix and a source vector internal coefficient:
Figure BDA0003941296360000138
Figure BDA0003941296360000139
Figure BDA00039412963600001310
Figure BDA00039412963600001311
Figure BDA0003941296360000141
Figure BDA0003941296360000142
Figure BDA0003941296360000143
wherein k is 1ij ,k 2ij ,k 3ij ,g ij ,p 1ij ,p 2ij Respectively a stiffness matrix K 1e ,K 2e ,K 3e ,G e Source vector P 1e ,P 2e The elements in the table, i, j, l are node labels, a, b, c are the side lengths of the subunits in x, y, z directions, respectively, and the unit is m, xi i 、η i 、ζ i Is the coordinate of node i on the parent cell.
And expanding the number of each node of each unit to the element position corresponding to the research area to obtain a linear equation set of the synthesized total stiffness matrix and the source vector:
Figure BDA0003941296360000144
where δ T is the impulse function array of the study region, T is the temperature array within the study region, K 1 ,K 2 ,K 3 G is the stiffness matrix in the region of interest, P 1 ,P 2 Omitting δ T for source vectors within the region of interest T After simplification, the method comprises the following steps:
Figure BDA0003941296360000145
wherein K = K 1 +K 2 +K 3 ,P=P 1 +P 2
Differentiating the time T in the transient heat transfer underground space temperature field, and after delta T, the temperature field is changed from T 1 Change to T 2 Then, then
Figure BDA0003941296360000146
Additionally, the interpolated temperature within the time of Δ t
Figure BDA0003941296360000147
Instead of T in the above equation, the linear equation set for the transient temperature field is:
Figure BDA0003941296360000148
by solving the above formula, after finite element numerical simulation is realized, the temperature field T (x, y, z) = d of the underground space at a given moment is obtained obs (refer to fig. 5), x, y, and z represent x, y, and z-axis directions, respectively.
S2: according to the temperature field d of the underground space obs Constructing an objective function of transient temperature forward modeling by utilizing Tikhonov regularization, acquiring a linear equation set consisting of a Hessian matrix and a gradient vector by adopting a Gauss-Newton optimization algorithm, solving the Jacobian matrix of a first derivative of a temperature numerical simulation response to a prediction model parameter in the Hessian matrix and the product of the transposition of the Jacobian matrix and an arbitrary vector by utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule, and solving to obtain a correction quantity delta m of the prediction model;
in this embodiment, the step S2 specifically includes: setting a priori information of a study region to be uniform half-space (m) ref =3.0W/m · K), substituting into the thermophysical parameter distribution and the boundary condition to perform transient finite element temperature numerical simulation, and obtaining an underground space temperature field of the research region under prior information (refer to fig. 6);
the objective function of the forward modeling of transient temperature is constructed by Tikhonov regularization and is as follows:
Φ=||D[d(m)-d obs ]|| 2 +λ||W(m-m apr )|| 2
wherein D is a data weighting matrix; d obs Is observed data, namely a temperature field T (x, y, z) of the underground space; m is a parameter vector of the current prediction model; d (m) is at the current prediction modelSimulation data under parameters; w is a smooth constraint matrix; λ is a damping factor; m is apr Is predictive model prior information;
based on a linear equation set of a total rigidity matrix and a source vector in a finite element, the heat conductivity coefficient k in the linear equation set is derived to obtain:
Figure BDA0003941296360000151
Figure BDA0003941296360000152
Figure BDA0003941296360000153
where K, G is the stiffness matrix in the region of interest, P 1 ,P 2 For the source vector in the region of interest, Δ T is the time interval, and the units s, K are the thermal conductivity array, and the units W/m.K, T i 、T i+1 Temperature arrays at time i and i +1, respectively. At the initial time (t = 0), the above equation is rewritten as:
Figure BDA0003941296360000154
wherein, T 0 Is an array of temperatures at the initial time. In order to ensure that the temperature value and the heat conductivity coefficient are not negative, the data parameters and the model parameters are set as the logarithm of the temperature and the logarithm of the heat conductivity coefficient in the prediction process: d = ln T, m = ln k; according to the chain derivation rule, the elements in the Jacobian matrix are:
Figure BDA0003941296360000155
wherein, i =1,2, …, N, N is the number of data parameters, j =1,2, …, M, M is the number of model parameters, d i For the number of i-th nodes in the study areaAccording to the parameter, m j Model parameter, T, representing the jth node in the study area i Representing the temperature, k, of the ith node in the region of interest j Representing the thermal conductivity of the jth node in the study area, the Jacobian matrix is written as:
J=RGR'
wherein R and R' are diagonal matrices,
Figure BDA0003941296360000161
R'=diag(k 1 ,k 2 ,…,k M ) G is a matrix composed of elements G (i, j);
for an arbitrary vector consisting of dimensionless elements v = (v) 1 ,v 2 ,…,v M ) T And the conversion formula based on the finite element rigidity matrix is as follows:
Figure BDA0003941296360000162
k (v) is a rigidity matrix when the heat conductivity coefficient array is v; and when t =0, there are
Figure BDA0003941296360000163
Wherein e is 1 ,e 2 ,...,e M Since all vectors are unit vectors, the product of Jacobian matrix J and arbitrary vector v can be obtained.
Similarly, the transposition of the Jacobian matrix and an arbitrary vector y = (y) composed of dimensionless elements are calculated 1 ,y 2 ,…,y N ) T The product of (a);
according to a Gauss-Newton algorithm, neglecting a second-order information item, obtaining the following linear equation system to solve the model correction quantity delta m:
HΔm=-g
wherein H is Hessian matrix, H = (J) T D T DJ+λW T W), g is gradient vector, g = [ J = T D T D(d-d obs )+λW T W(m-m apr )]。
S3: acquiring the optimal prediction model search step length meeting the Wolfe criterion according to the correction quantity delta m of the prediction model, and updating the parameters of the prediction model;
in this embodiment, the step S3 specifically includes: the Wolfe criterion needs to satisfy the following inequality:
Figure BDA0003941296360000171
wherein alpha is a search step length and is dimensionless, and c is a dimensionless number; in the process of obtaining the search step length, the initial value of alpha is set to be 1, and if the inequality is not satisfied, the initial value of alpha is set to be 1
α=ρα
Wherein rho is a decimal larger than 0 and smaller than 1, and the process is circulated until the Wolfe criterion is met; thus, the update formula for the model parameters is:
m k+1 =m k +αΔm k
wherein m is k 、m k+1 Model parameters, Δ m, at the kth and k +1 th iterations, respectively k The model correction quantity obtained in the k iteration is obtained.
S4: carrying out transient finite element temperature numerical simulation according to the updated prediction model parameters to obtain an underground space temperature field under the model parameters, calculating a data fitting difference between the model parameters and actual measurement data, if the data fitting difference is smaller than a preset fitting difference condition, determining the current prediction model parameters to be optimal thermal conductivity coefficient prediction results, and calculating temperature field distribution under the current model parameters to be used as a temperature field under the optimal prediction results (refer to fig. 7); otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
In this embodiment, the step S4 specifically includes: the calculation formula of the data fitting difference is as follows:
Figure BDA0003941296360000172
wherein RMS is root mean square and is used to measureThe closer to 0, the smaller the deviation between the measured value and the predicted value is, the more reliable the prediction result is; d (m) i 、d obs i The ith predicted value and the ith measured value are respectively, and n is the number of observed data and is dimensionless.
The preset fitting difference is 0.07, if the data fitting difference is smaller than the preset fitting difference, the iteration is stopped, and the current prediction model parameters are output as the optimal prediction result (refer to fig. 8); otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m. Referring to FIG. 9, FIG. 9 is a graph of data fitting difference versus number of predicted iterations for the present invention; in this embodiment, the RMS value between the measured value and the predicted value in each prediction process is calculated. The experimental result shows that the RMS between the measured value and the predicted value is converged from 1.1930 to 0.0693, and the fitting degree is good.
Referring to fig. 3 and 8, wherein fig. 3 is a graph of the thermal conductivity of the real model of the present invention, and fig. 8 is a graph of the thermal conductivity of the optimal prediction model of the present invention; in this embodiment, an Abnormal Volume Overlap Ratio (AVOR) is defined to evaluate the accuracy of the prediction results:
Figure BDA0003941296360000181
wherein, V o To predict the volume of the portion of the anomaly that intersects the true anomaly, m 3 ;V p To predict the volume of the anomaly, m 3 ;V t To predict an anomaly, m 3 . The experimental result shows that AVOR of the upper thin layer anomalous body reaches 88.97%; AVOR of the lower left side abnormal body reaches 94.19%; the AVOR of the lower right outliers reached 93.68%. Therefore, the optimal prediction model is high in fitting degree with the real model.
Referring to fig. 10, fig. 10 is a block diagram of a deep formation thermal conductivity three-dimensional transient prediction apparatus according to the present invention, which includes the following modules:
the transient finite element numerical simulation module 1 is used for acquiring the underground water migration rule in the uniform half-space research area andconstructing a heat conductivity coefficient abnormal body in the research area and setting boundary conditions according to the response mechanism of the heat conductivity coefficient abnormal body, carrying out transient finite element temperature numerical simulation by combining the distribution of the thermophysical parameters in the research area and the boundary conditions, and acquiring an underground space temperature field T (x, y, z) = d at a given moment obs X, y and z respectively represent x, y and z axis directions;
a model correction quantity calculation module 2 for calculating the correction quantity according to the temperature field d of the underground space obs Establishing an objective function of transient temperature forward modeling by utilizing Tikhonov regularization, acquiring a linear equation set consisting of a Hessian matrix and a gradient vector by adopting a Gauss-Newton algorithm, solving the Jacobian matrix of a first-order derivative of temperature numerical simulation response to a prediction model parameter in the Hessian matrix and the product of the transposition of the Jacobian matrix and an arbitrary vector by utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule, and solving to obtain a correction quantity delta m of the prediction model;
the model parameter updating module 3 is used for acquiring the optimal prediction model search step length meeting the Wolfe criterion according to the correction quantity delta m of the prediction model and updating the parameters of the prediction model;
the thermal conductivity coefficient prediction module 4 is used for carrying out transient finite element temperature numerical simulation according to the updated prediction model parameters to obtain an underground space temperature field under the model parameters, calculating a data fitting difference between the model parameters and actually measured data, and outputting the current prediction model parameters as an optimal thermal conductivity coefficient prediction result if the current prediction model parameters are smaller than a preset fitting difference condition; otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
Furthermore, as a preferred embodiment, the present invention provides an electronic device, comprising a memory, a processor and a computer program stored on the memory and executable on the processor, wherein the processor when executing the program implements the steps of the deep formation thermal conductivity three-dimensional transient prediction method, for example: acquiring the underground water migration rule in a uniform half-space research area and a response mechanism thereof to the heat conductivity coefficient, constructing a heat conductivity coefficient abnormal body in the research area, setting boundary conditions, and combining the research area and the heat conductivity coefficient abnormal bodyCarrying out transient finite element temperature numerical simulation on the distribution of internal thermal physical parameters and the boundary conditions to obtain an underground space temperature field T (x, y, z) = d at a given moment obs X, y and z respectively represent the directions of x, y and z axes; according to the temperature field d of the underground space obs Establishing an objective function of transient temperature forward modeling by utilizing Tikhonov regularization, acquiring a linear equation set consisting of a Hessian matrix and a gradient vector by adopting a Gauss-Newton algorithm, solving the Jacobian matrix of a first-order derivative of temperature numerical simulation response to a prediction model parameter in the Hessian matrix and the product of the transposition of the Jacobian matrix and an arbitrary vector by utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule, and solving to obtain a correction quantity delta m of the prediction model; acquiring the optimal prediction model search step length meeting the Wolfe criterion according to the correction quantity delta m of the prediction model, and updating the parameters of the prediction model; carrying out transient finite element temperature numerical simulation according to the updated prediction model parameters to obtain an underground space temperature field under the model parameters, calculating a data fitting difference between the model parameters and actual measurement data, and if the data fitting difference is smaller than a preset fitting difference condition, outputting the current prediction model parameters as an optimal heat conductivity coefficient prediction result; otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
It should be noted that, in this document, the terms "comprises," "comprising," or any other variation thereof, are intended to cover a non-exclusive inclusion, such that a process, method, article, or system that comprises a list of elements does not include only those elements but may include other elements not expressly listed or inherent to such process, method, article, or system. Without further limitation, an element defined by the phrase "comprising a … …" does not exclude the presence of another identical element in a process, method, article, or system that comprises the element.
The above-mentioned serial numbers of the embodiments of the present invention are only for description, and do not represent the advantages and disadvantages of the embodiments. In the unit claims enumerating several means, several of these means may be embodied by one and the same item of hardware. The use of the words first, second, third and the like do not denote any order, but rather the words first, second and the like may be interpreted as indicating any order.
The above description is only a preferred embodiment of the present invention, and not intended to limit the scope of the present invention, and all modifications of equivalent structures and equivalent processes, which are made by using the contents of the present specification and the accompanying drawings, or directly or indirectly applied to other related technical fields, are included in the scope of the present invention.

Claims (10)

1. A three-dimensional transient prediction method for the heat conductivity coefficient of a deep stratum is characterized by comprising the following steps:
acquiring an underground water migration rule in a uniform half-space research area and a response mechanism thereof to a heat conductivity coefficient, constructing a heat conductivity coefficient abnormal body in the research area, setting boundary conditions, carrying out transient finite element temperature numerical simulation by combining the distribution of thermophysical parameters in the research area and the boundary conditions, and acquiring an underground space temperature field T (x, y, z) = d at a given moment obs X, y and z respectively represent x, y and z axis directions;
according to the temperature field d of the underground space obs Establishing an objective function of transient temperature forward modeling by utilizing Tikhonov regularization, acquiring a linear equation set consisting of a Hessian matrix and a gradient vector by adopting a Gauss-Newton algorithm, solving the Jacobian matrix of a first-order derivative of temperature numerical simulation response to a prediction model parameter in the Hessian matrix and the product of the transposition of the Jacobian matrix and an arbitrary vector by utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule, and solving to obtain a correction quantity delta m of the prediction model;
acquiring the optimal prediction model search step length meeting the Wolfe criterion according to the correction quantity delta m of the prediction model, and updating the parameters of the prediction model;
carrying out transient finite element temperature numerical simulation according to the updated prediction model parameters to obtain an underground space temperature field under the model parameters, calculating a data fitting difference between the model parameters and actual measurement data, and if the data fitting difference is smaller than a preset fitting difference condition, outputting the current prediction model parameters as an optimal heat conductivity coefficient prediction result; otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
2. The deep formation thermal conductivity three-dimensional transient prediction method of claim 1, wherein the law of groundwater migration and its response mechanism to thermal conductivity comprises:
setting the migration rates of the movement of underground water in the directions x, y and z of a research area as v x 、v y 、v z And the fluid is transported and immersed into the surrounding medium, so that the heat conductivity coefficient of the medium is changed, and the following conditions are met:
k=(1-φ)·k s +φ·k f
wherein, k s 、k f The real thermal conductivity, the medium thermal conductivity and the fluid thermal conductivity are respectively, and the unit is W/m.K; phi is the porosity of the medium and is dimensionless.
3. The deep formation thermal conductivity three-dimensional transient prediction method of claim 1, in which the boundary conditions comprise:
Figure FDA0003941296350000011
wherein, both alpha and beta are constants, alpha =1, beta = T up The time is the upper boundary condition; α =0, β = q down The lower boundary condition is adopted; a boundary condition of four weeks is defined when α =0 and β = 0; t is temperature, in degrees Celsius up Is the surface temperature in units of DEG C, n is the normal direction of the boundary, K is the thermal conductivity in units of W/m.K, q down Is a heat flow value in W/m 2
4. The deep formation thermal conductivity three-dimensional transient prediction method of claim 1, wherein the transient finite element temperature numerical simulation is performed in a discrete manner by using a hexahedral structure, and the internal temperature T, the thermal conductivity k and the fluid heat capacity μ of the hexahedral unit f Medium heat capacity [ mu ] s Fluid transport rate v, heat source QLinear interpolation is adopted:
Figure FDA0003941296350000021
Figure FDA0003941296350000022
wherein i is a node label, T i 、k i 、μ fi 、μ si 、v i 、Q i The temperature, the heat conductivity coefficient, the fluid heat capacity, the medium heat capacity, the fluid flow rate and the heat source of each node in the unit are respectively, and the unit is W/m. K, J/(m) 3 ·K)、J/(m 3 ·K)、m/s、W/m 3 ;N i Is a shape function, satisfies:
Figure FDA0003941296350000023
wherein ξ i 、η i 、ζ i The coordinate of the node i on the parent unit is dimensionless; xi, eta and zeta in the parent unit and coordinates x, y and z in the child unit are in a relational formula:
Figure FDA0003941296350000024
wherein x is 0 、y 0 、z 0 Respectively being the middle points of the subunits in the directions of x, y and z axes; a. b and c are the side lengths of the subunits in the x, y and z directions respectively, and the unit is m.
5. The deep formation thermal conductivity three-dimensional transient prediction method of claim 1, wherein the step of performing a transient finite element temperature numerical simulation in combination with the distribution of thermophysical parameters in the study zone and the boundary conditions to obtain the subsurface spatial temperature field at a given time comprises:
the differential equation for a given transient heat transfer subsurface space temperature field is:
Figure FDA0003941296350000025
wherein K is the thermal conductivity coefficient, the unit W/m.K, T is the temperature, the unit degree centigrade, mu f The unit J/(m) is the heat capacity of the fluid 3 ·K),μ s The unit J/(m) is the heat capacity of the medium 3 K), v is the fluid flow rate, m/s, Q is the heat source, W/m 3 T is the time, in units of s,
Figure FDA0003941296350000031
in order to be the divergence degree,
Figure FDA0003941296350000032
is a gradient;
multiplying and integrating the differential equation of the temperature field of the underground space by an impulse function delta T:
Figure FDA0003941296350000033
wherein, omega is a research region, gamma is a boundary, and delta T is an impulse function;
according to the Hamiltonian operation rule in the field theory and On Gao Gongshi, and substituting the boundary condition, the first term in the above formula is rewritten as:
Figure FDA0003941296350000034
the integral equation of the transient heat transfer underground space temperature field is obtained as follows:
Figure FDA0003941296350000035
carrying out discrete subdivision and linear interpolation on the integral equation of the temperature field of the transient heat transfer underground space, and solving to obtain the unit integrals of all terms in the integral equation:
Figure FDA0003941296350000036
Figure FDA0003941296350000037
Figure FDA0003941296350000038
Figure FDA0003941296350000039
Figure FDA00039412963500000310
Figure FDA00039412963500000311
wherein T is time, unit s, delta T e For arrays of impulse functions within a cell, T e For temperature arrays within the cell, K 1e ,K 2e ,K 3e ,G e Is a stiffness matrix within the cell, P 1e ,P 2e For the source vector in the unit, solving to obtain a rigidity matrix and a source vector internal coefficient:
Figure FDA00039412963500000312
Figure FDA0003941296350000041
Figure FDA0003941296350000042
Figure FDA0003941296350000043
Figure FDA0003941296350000044
Figure FDA0003941296350000045
Figure FDA0003941296350000046
wherein k is 1ij ,k 2ij ,k 3ij ,g ij ,p 1ij ,p 2ij Respectively a stiffness matrix K 1e ,K 2e ,K 3e ,G e And a source vector P 1e ,P 2e The elements in the table, i, j, l are node labels, a, b, c are the side length of the subunit in x, y, z direction, unit is m, xi i 、η i 、ζ i Coordinates of the node i on the parent unit;
and expanding the number of each node of each unit to the element position corresponding to the research area to obtain a linear equation set of the synthesized total stiffness matrix and the source vector:
Figure FDA0003941296350000047
where δ T is the impulse function array of the study region, T is the temperature array within the study region, K 1 ,K 2 ,K 3 G is respectively a stiffness matrix in the investigation region, P 1 ,P 2 Omitting δ T for source vectors within the region of interest T After simplification, the following steps are included:
Figure FDA0003941296350000048
wherein K = K 1 +K 2 +K 3 ,P=P 1 +P 2
Differentiating the time T in the transient heat transfer underground space temperature field, and after delta T, the temperature field is changed from T 1 Change to T 2 Then, then
Figure FDA0003941296350000051
Additionally, the interpolated temperature within the time of Δ t
Figure FDA0003941296350000052
Instead of T in the above equation, the linear equation set for the transient temperature field is:
Figure FDA0003941296350000053
by solving the above formula, after finite element numerical simulation is realized, the temperature field T (x, y, z) = d of the underground space at a given moment is obtained obs And x, y and z represent x, y and z axis directions, respectively.
6. The deep formation thermal conductivity three-dimensional transient prediction method of claim 5, wherein the objective function of the Tikhonov regularization-constructed transient temperature forward modeling is as follows:
Φ=||D[d(m)-d obs ]|| 2 +λ||W(m-m apr )|| 2
wherein D is a data weighting matrix; d obs For observing data, i.e.A subsurface space temperature field T (x, y, z); m is a parameter vector of the current prediction model; d (m) is simulation data under the current prediction model parameters; w is a smooth constraint matrix; λ is a damping factor; m is apr Is predictive model prior information;
based on a linear equation set of a total rigidity matrix and a source vector in a finite element, the heat conductivity coefficient k in the linear equation set is derived to obtain:
Figure FDA0003941296350000054
Figure FDA0003941296350000055
Figure FDA0003941296350000056
where K, G is the stiffness matrix in the region of interest, P 1 ,P 2 For the source vector in the region of interest, Δ T is the time interval, and the units s, K are the thermal conductivity array, and the units W/m.K, T i 、T i+1 Temperature arrays at the time of i and i +1 respectively; at the initial time (t = 0), the above equation is rewritten as:
Figure FDA0003941296350000057
wherein, T 0 Is a temperature array at the initial time; setting data parameters and model parameters as logarithm of temperature and logarithm of thermal conductivity in the prediction process: d = lnT, m = lnk; according to the chain derivation rule, the elements in the Jacobian matrix are:
Figure FDA0003941296350000061
wherein i=1,2, …, N is the number of data parameters, j =1,2, …, M is the number of model parameters, d i For the data parameter of the ith node in the investigation region, m j Model parameter, T, representing the jth node in the study area i Denotes the temperature, k, of the ith node in the investigation region j Representing the thermal conductivity of the jth node in the study area, the Jacobian matrix is written as:
J=RGR'
wherein R and R' are diagonal matrices,
Figure FDA0003941296350000062
R'=diag(k 1 ,k 2 ,…,k M ) G is a matrix composed of elements G (i, j);
for an arbitrary vector consisting of dimensionless elements v = (v) 1 ,v 2 ,…,v M ) T And the conversion formula based on the finite element rigidity matrix is as follows:
Figure FDA0003941296350000063
k (v) is a rigidity matrix when the heat conductivity coefficient array is v; and when t =0, there are
Figure FDA0003941296350000064
Wherein e is 1 ,e 2 ,...,e M All the vectors are unit vectors, and the product of the Jacobian matrix J and any vector v is obtained;
similarly, the transposition of the Jacobian matrix and an arbitrary vector y = (y) composed of dimensionless elements are calculated 1 ,y 2 ,…,y N ) T The product of (a);
according to a Gauss-Newton algorithm, neglecting a second-order information item, obtaining the following linear equation system to solve the model correction quantity delta m:
HΔm=-g
wherein H is Hessian matrix, H = (J) T D T DJ+λW T W), g is gradient vector, g = [ J = T D T D(d-d obs )+λW T W(m-m apr )]。
7. The method of claim 6, wherein the step size of the optimal prediction model search satisfying Wolfe's criteria is obtained based on the prediction model correction Δ m by satisfying the following inequality:
Figure FDA0003941296350000071
wherein alpha is a search step length and is dimensionless, and c is a dimensionless number; in the process of obtaining the search step length, the initial value of alpha is 1, and if the inequality is not satisfied, the initial value of alpha is made to be 1
α=ρα
Wherein rho is a decimal larger than 0 and smaller than 1, and the process is circulated until the Wolfe criterion is met; the updating formula of the model parameters is as follows:
m k+1 =m k +αΔm k
wherein m is k 、m k+1 Model parameters, Δ m, at the kth and k +1 th iterations, respectively k The model correction quantity obtained in the k iteration is obtained.
8. The deep formation thermal conductivity three-dimensional transient prediction method of claim 1, wherein the data fitting difference is calculated as:
Figure FDA0003941296350000072
wherein, the RMS is root mean square and is used for measuring the deviation between the measured value and the predicted value; d (m) i 、d obs i The ith predicted value and the ith measured value are respectively, n is the number of observed data, and the dimension is not existed;
if the data fitting difference is smaller than the preset fitting difference, stopping iteration and outputting the current prediction model parameters as the optimal prediction result; otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
9. The three-dimensional transient prediction device for the deep stratum heat conductivity coefficient is characterized by comprising the following modules:
the transient finite element numerical simulation module is used for acquiring the underground water migration rule in a uniform half-space research area and the response mechanism of the underground water migration rule to the heat conductivity coefficient, constructing a heat conductivity coefficient abnormal body in the research area, setting boundary conditions, carrying out transient finite element temperature numerical simulation by combining the distribution of the thermophysical parameters in the research area and the boundary conditions, and acquiring the underground space temperature field T (x, y, z) = d at a given moment obs X, y and z respectively represent the directions of x, y and z axes;
a model correction amount calculation module for calculating the correction amount according to the temperature field d of the underground space obs Establishing an objective function of transient temperature forward modeling by utilizing Tikhonov regularization, acquiring a linear equation set consisting of a Hessian matrix and a gradient vector by adopting a Gauss-Newton algorithm, solving the Jacobian matrix of a first-order derivative of temperature numerical simulation response to a prediction model parameter in the Hessian matrix and the product of the transposition of the Jacobian matrix and an arbitrary vector by utilizing a Jacobian-free Krylov subspace technology based on a time chain derivation rule, and solving to obtain a correction quantity delta m of the prediction model;
the model parameter updating module is used for acquiring the optimal prediction model searching step length meeting the Wolfe criterion according to the correction quantity delta m of the prediction model and updating the parameters of the prediction model;
the thermal conductivity coefficient prediction module is used for carrying out transient finite element temperature numerical simulation according to the updated prediction model parameters to obtain an underground space temperature field under the model parameters, calculating the data fitting difference between the data and the actually measured data, and outputting the current prediction model parameters as the optimal thermal conductivity coefficient prediction result if the current prediction model parameters are smaller than the preset fitting difference condition; otherwise, updating the underground space temperature field d (m) according to the model parameters and solving a new prediction model correction quantity delta m.
10. An electronic device comprising a memory, a processor, and a computer program stored on the memory and executable on the processor, wherein the processor when executing the program performs the steps of the deep formation thermal conductivity three-dimensional transient prediction method of any of claims 1-8.
CN202211424242.4A 2022-11-14 2022-11-14 Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment Pending CN115758820A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211424242.4A CN115758820A (en) 2022-11-14 2022-11-14 Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211424242.4A CN115758820A (en) 2022-11-14 2022-11-14 Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment

Publications (1)

Publication Number Publication Date
CN115758820A true CN115758820A (en) 2023-03-07

Family

ID=85370777

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211424242.4A Pending CN115758820A (en) 2022-11-14 2022-11-14 Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment

Country Status (1)

Country Link
CN (1) CN115758820A (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116379793A (en) * 2023-06-02 2023-07-04 青岛智控菲特软件科技有限公司 Submerged arc furnace short-net regulation and control data processing method

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116379793A (en) * 2023-06-02 2023-07-04 青岛智控菲特软件科技有限公司 Submerged arc furnace short-net regulation and control data processing method
CN116379793B (en) * 2023-06-02 2023-08-15 青岛智控菲特软件科技有限公司 Submerged arc furnace short-net regulation and control data processing method

Similar Documents

Publication Publication Date Title
Chua et al. An inverse ocean modeling system
Oliver et al. Markov chain Monte Carlo methods for conditioning a permeability field to pressure data
CN114707385A (en) Deep stratum heat conductivity coefficient three-dimensional prediction method and device based on Krylov subspace
Ardito et al. Diagnostic analysis of concrete dams based on seasonal hydrostatic loading
Solonen et al. Simulation-based optimal design using a response variance criterion
Wang et al. A model of surface heat fluxes based on the theory of maximum entropy production
Yang et al. Multiscale computation for transient heat conduction problem with radiation boundary condition in porous materials
Wen et al. An on-line extended Kalman filtering technique for reconstructing the transient heat flux and temperature field in two-dimensional participating media
CN115758820A (en) Deep stratum heat conductivity coefficient three-dimensional transient prediction method and device and electronic equipment
Li et al. Efficient geostatistical inverse methods for structured and unstructured grids
Xu et al. Estimation of turbulent heat fluxes by assimilation of land surface temperature observations from GOES satellites into an ensemble Kalman smoother framework
Ishitsuka et al. Bayesian and neural network approaches to estimate deep temperature distribution for assessing a supercritical geothermal system: evaluation using a numerical model
Alpar et al. Estimation of soils thermophysical characteristics in a nonlinear inverse heat transfer problem
Liu et al. Reservoir automatic history matching: Methods, challenges, and future directions.
Wang et al. Two-dimensional steady-state boundary shape inversion of CGM-SPSO algorithm on temperature information
US20230176030A1 (en) Systems and Methods for Determining Ground Water-Surface Water Interactions
Liu et al. Observing system simulation experiments in geomagnetic data assimilation
Penazzi et al. Transfer function estimation with SMC method for combined heat transfer: insensitivity to detail refinement of complex geometries
Chang et al. Optimal observation network design for parameter structure identification in groundwater modeling
Azam et al. Parametric PGD model used with orthogonal polynomials to assess efficiently the building's envelope thermal performance
D’Antona et al. Analysis of the sensor placement for optimal temperature distribution reconstruction
Degen et al. Uncertainty quantification for basin-scale conductive models
Nho et al. Scaling of thermal conductivity of helium confined in pores
Butler et al. Simulation-based optimal experimental design: A measure-theoretic perspective
Azimi et al. Structured multiblock body-fitted grids solution of transient inverse heat conduction problems in an arbitrary geometry

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination