Time-domain integral solution method for transient seepage in three-dimensional heterogeneous semi-infinite dam

Abstract

A time-domain integral solution method for transient seepage in dams with three-dimensional heterogeneous semi-infinite domains is proposed. This method includes constructing a novel semi-analytical method based on the scaled boundary finite element method (SBFEM) using scaled surfaces to study the three-dimensional unsteady seepage problem in heterogeneous semi-infinite domain foundations. A new coordinate transformation technique is introduced to accurately describe the irregular and discontinuous material region. Based on the three-dimensional diffusion equation, the transient hydraulic head control equation is systematically derived using energy functional theory and variational principles, ensuring the physical consistency modeling of the unsteady seepage process in heterogeneous media. To handle boundary conditions at infinity, the unsteady seepage matrix is ​​constructed using the weighted residual method and solved efficiently using the continued fraction expansion method. The time integral of the hydraulic response is achieved by combining the improved refined time history integration method (MPTSIM) with the fourth-order Runge-Kutta method, exhibiting higher stability and accuracy over long time intervals.

Description

Technical Field

[0001] This invention relates to the field of time-domain integral solution technology for transient seepage flow in dams containing three-dimensional heterogeneous semi-infinite domains, and particularly to a time-domain integral solution method for transient seepage flow in dams containing three-dimensional heterogeneous semi-infinite domains. Background Technology

[0002] This invention relates to the field of seepage analysis technology in geotechnical and hydraulic engineering, specifically to a time-domain integral solution method for transient seepage in dams with three-dimensional heterogeneous semi-infinite domains. This method is applicable to seepage stability assessment of dam-reservoir-foundation systems, dam seepage prevention design, and seepage response prediction under rainfall-water level fluctuation conditions, falling within the scope of high-precision numerical analysis technology for seepage in complex foundations. Its core lies in efficiently and stably solving the evolution law of the three-dimensional transient seepage field under heterogeneous material and geometrically semi-infinite extension conditions. In dam engineering, seepage analysis is a crucial step in ensuring structural safety and long-term service performance. Factors such as temperature difference, water level fluctuation, and rainfall infiltration can all cause changes in the seepage field in the dam body and its foundation. Inaccurate calculations may lead to serious consequences such as seepage failure, piping, or dam foundation instability. Therefore, accurately describing the motion law of fluid in porous media in the time domain is the foundation for conducting dam seepage safety assessment and optimization design. The Scaled Boundary Finite Element Method (SBFEM), as a semi-analytical method, constructs parameterized coordinates through the principle of "similarity center-boundary visibility," enabling analytical solutions in the radial direction and making it naturally suitable for semi-infinite domain problems. Although this method has been preliminarily applied to seepage field analysis, several key issues remain in three-dimensional unsteady-state scenarios: limited characterization capability for heterogeneous materials, low time-domain solution efficiency, and difficulty in simultaneously ensuring interface continuity and physical consistency.

[0003] Therefore, there is an urgent need for a transient seepage analysis method that can simultaneously adapt to three-dimensional heterogeneous semi-infinite domain structures and possesses high time-domain integration accuracy and high computational efficiency. This method should accurately capture the evolution characteristics of the transient seepage field in the dam-foundation system while ensuring physical continuity and numerical stability, providing highly reliable numerical support for dam safety assessment and seepage prevention design. Summary of the Invention

[0004] The main objective of this invention is to provide a time-domain integral solution method for transient seepage in dams containing three-dimensional heterogeneous semi-infinite domains. This method aims to solve the technical problems of poor adaptability to three-dimensional heterogeneous semi-infinite domains, low accuracy and transient efficiency in infinite domains, and distortion of physical meaning in existing numerical methods.

[0005] To achieve the above objectives, this invention provides a time-domain integral solution method for transient seepage in a dam containing a three-dimensional heterogeneous semi-infinite domain. The method includes the following steps:

[0006] S1. Construct similar surfaces parallel to the boundary surface, and define the scaled boundary finite element SBFEM coordinate system: circumferential coordinates. η , ζ Parallel boundary surface, radial coordinates ξ Given the perpendicular boundary surface, we provide the coordinate expressions for any point in the problem domain. We introduce shape functions to discretize the points on the boundary surface and similar surfaces, and obtain the mapping of any point in the problem domain within the surface element to the SBFEM coordinate system. We realize the transformation between the Cartesian coordinate system and the SBFEM coordinate system through the Jacobian matrix of coordinate transformation. We further obtain the inverse matrix of the Jacobian matrix and derive the differential operator expression for the transient seepage problem.

[0007] S2. Using the total energy functional characterizing the system's energy state as the variational basis, derive the energy functional expression for transient seepage. Discretize the field variables in the SBFEM coordinate system to obtain the interpolated forms of the head and gradient at any point in the problem domain. Substitute these interpolations into the total energy functional and then apply them to the circumferential coordinate system. By integrating the equations, the total energy functional expression is transformed into a discrete form that depends only on the radial coordinates. By applying the first variational principle, the first derivative of the total energy functional at the stationary point is set to zero. Then, the differential equations are extracted from the weak form equations using the integration by parts method, and the control equations in the SBFEM coordinate system with respect to the head are obtained.

[0008] S3. Based on the principle of virtual work and the weighted residual method, the internal equivalent node flux is obtained. According to the inverse relationship between the internal and external equivalent node fluxes, the relationship between the node head and the external equivalent node flux is described. Furthermore, the governing equations of the transient seepage matrix are derived. Dimensionless frequency and dimensionless seepage matrix are introduced, and the radial coordinates of the dimensionless seepage matrix are calculated. and frequency The partial derivatives are used to transform the original equation into the transient seepage matrix control equation in the frequency domain by variable substitution.

[0009] S4. Apply the continuous fraction method to... The transient seepage matrix approaching infinity is decomposed into a power series form, and the coefficient matrix is ​​obtained by solving a general eigenvalue problem. By setting the coefficients of higher-order terms to zero, and solving the generalized eigenvalue problem and the Lyapunov equation, we obtain the normalized frequency. The expression for the function to be solved is obtained by expanding the remaining lower-order terms through a power series, recursively calculating the subsequent coefficient matrix, setting a termination condition to truncate the expansion, and finally obtaining the expression for the transient seepage matrix.

[0010] S5, set frequency The governing equations for the steady-state seepage problem are obtained, and these equations are expressed as eigenvalue problems. The steady-state seepage matrix is ​​obtained by Schur decomposition, which describes the relationship between flux and head on the boundary surface of the problem domain. The inverse Fourier transform is applied to obtain the head expression in the time domain. The improved refined time history integral method MPTSIM is used, and the head value at any time step is calculated by matrix exponentiation. By executing a recursive loop, the head value at each time step on the boundary surface is obtained. A function of head with respect to the radial coordinate is defined, and the head distribution in the semi-infinite domain is solved using the fourth-order Runge-Kutta method.

[0011] The beneficial effects of this invention are:

[0012] (1) This invention studies the transient seepage problem of dams with three-dimensional heterogeneous semi-infinite domain by constructing a novel semi-analytical method based on the scaled boundary finite element method (SBFEM) of similar surfaces. The proposed method reconstructs the problem in a scaled coordinate system based on scaled surfaces, which can accurately and efficiently model complex geological heterogeneity and unbounded geometry.

[0013] (2) This invention introduces a new coordinate transformation technique for accurately describing irregular and discontinuous material regions. Based on the three-dimensional diffusion equation, the transient hydraulic head control equation is systematically derived using energy functional theory and variational principles to ensure the physical consistency modeling of the unsteady seepage process in heterogeneous media.

[0014] (3) To handle boundary conditions at infinity, this invention constructs an unsteady seepage matrix using the weighted residual method and solves it efficiently using the continued fraction expansion method. The time integral of the hydraulic response is achieved by combining the improved fine time history integration method (MPTSIM) with the fourth-order Runge-Kutta method, which has higher stability and accuracy over long time intervals. Attached Figure Description

[0015] Figure 1 This is a flowchart illustrating the time-domain integral solution method for transient seepage in a three-dimensional heterogeneous semi-infinite domain dam according to the present invention.

[0016] Figure 2 for Figure 1 The schematic diagrams shown include (a) a cross-section of a typical heterogeneous half-space; (b) a structural system with a heterogeneous half-space; (c) the near field (finite field); and (d) the far field (infinite field).

[0017] Figure 3 This is a schematic diagram of the computational model for the first embodiment – ​​a three-dimensional triangular channel layered semi-infinite domain;

[0018] Figure 4 This is a schematic diagram of the model unit discretization of the first embodiment - a three-dimensional triangular channel layered semi-infinite domain;

[0019] Figure 5 Comparison of transient seepage head in a three-dimensional triangular channel layered semi-infinite domain calculated by this invention and other methods: (a) Water level comparison when the water level rises; (b) Water level comparison when the water level falls.

[0020] Figure 6 This is a schematic diagram of the computational model for the second embodiment – ​​a semi-infinite domain foundation for a three-dimensional irregular homogeneous embankment;

[0021] Figure 7 This is a schematic diagram of the model element discretization for the second embodiment - a three-dimensional irregular homogeneous embankment semi-infinite domain foundation;

[0022] Figure 8 for Figure 4 Element mapping between similar surfaces and boundary surfaces;

[0023] Figure 9 for Figure 7 Element mapping between similar surfaces and boundary surfaces. Detailed Implementation

[0024] It should be understood that the specific embodiments described herein are for illustrative purposes only and are not intended to limit the scope of the invention. The realization of the invention's objectives, functional characteristics, and advantages will be further explained in conjunction with the embodiments and with reference to the accompanying drawings.

[0025] As shown in the attached diagram Figure 1 This is a flowchart illustrating the time-domain integral solution method for transient seepage in a three-dimensional heterogeneous semi-infinite domain dam according to the present invention. Figure 2 for Figure 1 The schematic diagram shown below; Figure 3 This is a schematic diagram of the computational model for the first embodiment – ​​a three-dimensional triangular channel layered semi-infinite domain; Figure 4 This is a schematic diagram of the model unit discretization of the first embodiment - a three-dimensional triangular channel layered semi-infinite domain; Figure 5 A comparison diagram of transient seepage head in a three-dimensional triangular channel layered semi-infinite domain calculated by the present invention and other methods; Figure 6 This is a schematic diagram of the computational model for the second embodiment – ​​a semi-infinite domain foundation for a three-dimensional irregular homogeneous embankment; Figure 7 This is a schematic diagram of the model element discretization for the second embodiment - a three-dimensional irregular homogeneous embankment semi-infinite domain foundation; Figure 8 for Figure 4 Element mapping between similar surfaces and boundary surfaces; Figure 9 for Figure 7 Element mapping between similar surfaces and boundary surfaces. The detailed method includes the following steps:

[0026] S1. Construct similar surfaces parallel to the boundary surface, and define the scaled boundary finite element SBFEM coordinate system: circumferential coordinates. η , ζ Parallel boundary surface, radial coordinates ξ Given the perpendicular boundary surface, we provide the coordinate expressions for any point in the problem domain. We introduce shape functions to discretize the points on the boundary surface and similar surfaces, and obtain the mapping of any point in the problem domain within the surface element to the SBFEM coordinate system. We realize the transformation between the Cartesian coordinate system and the SBFEM coordinate system through the Jacobian matrix of coordinate transformation. We further obtain the inverse matrix of the Jacobian matrix and derive the differential operator expression for the transient seepage problem.

[0027] S2. Using the total energy functional characterizing the system's energy state as the variational basis, derive the energy functional expression for transient seepage. Discretize the field variables in the SBFEM coordinate system to obtain the interpolated forms of the head and gradient at any point in the problem domain. Substitute these interpolations into the total energy functional and then apply them to the circumferential coordinate system. By integrating the equations, the total energy functional expression is transformed into a discrete form that depends only on the radial coordinates. By applying the first variational principle, the first derivative of the total energy functional at the stationary point is set to zero. Then, the differential equations are extracted from the weak form equations using the integration by parts method, and the control equations in the SBFEM coordinate system with respect to the head are obtained.

[0028] S3. Based on the principle of virtual work and the weighted residual method, the internal equivalent node flux is obtained. According to the inverse relationship between the internal and external equivalent node fluxes, the relationship between the node head and the external equivalent node flux is described. Furthermore, the governing equations of the transient seepage matrix are derived. Dimensionless frequency and dimensionless seepage matrix are introduced, and the radial coordinates of the dimensionless seepage matrix are calculated. and frequency The partial derivatives are used to transform the original equation into the transient seepage matrix control equation in the frequency domain by variable substitution.

[0029] S4. Apply the continuous fraction method to... The transient seepage matrix approaching infinity is decomposed into a power series form, and the coefficient matrix is ​​obtained by solving a general eigenvalue problem. By setting the coefficients of higher-order terms to zero, and solving the generalized eigenvalue problem and the Lyapunov equation, we obtain the normalized frequency. The expression for the function to be solved is obtained by expanding the remaining lower-order terms through a power series, recursively calculating the subsequent coefficient matrix, setting a termination condition to truncate the expansion, and finally obtaining the expression for the transient seepage matrix.

[0030] S5, set frequency The governing equations for the steady-state seepage problem are obtained, and these equations are expressed as eigenvalue problems. The steady-state seepage matrix is ​​obtained by Schur decomposition, which describes the relationship between flux and head on the boundary surface of the problem domain. The inverse Fourier transform is applied to obtain the head expression in the time domain. The improved refined time history integral method MPTSIM is used, and the head value at any time step is calculated by matrix exponentiation. By executing a recursive loop, the head value at each time step on the boundary surface is obtained. A function of head with respect to the radial coordinate is defined, and the head distribution in the semi-infinite domain is solved using the fourth-order Runge-Kutta method.

[0031] Further, in step S1, a similar surface parallel to the boundary surface is constructed, and the scaled boundary finite element SBFEM coordinate system is defined: circumferential coordinates. η , ζ Parallel boundary surface, radial coordinates ξ Given a perpendicular boundary surface, provide the coordinate expression for any point in the problem domain, as follows:

[0032]

[0033] in, , and Represents the coordinates of a point in the problem domain. , and Represents the coordinates of points on similar surfaces. , and Represents the coordinates of points on the boundary surface. This represents the radial coordinate of the SBFEM coordinate system.

[0034] We introduce shape functions to discretize points on the boundary surface and similar surfaces:

[0035]

[0036]

[0037] in, The shape function representing the surface element. , , These are the node coordinate vectors of points on similar surfaces within the surface element; , , This represents the nodal coordinate vector of a point on the boundary surface within a surface element.

[0038] Find the mapping of any point in the problem domain within the surface element to the SBFEM coordinate system:

[0039]

[0040] in, , , Represents the nodal coordinate vector of a point on the boundary surface within a surface element;

[0041] Achieving a Cartesian coordinate system through the Jacobian matrix of coordinate transformation ( x , y , z ) and SBFEM coordinate system ( ξ , η , ζ The transformation of the Jacobian matrix is ​​expressed as follows:

[0042]

[0043] in, Represents the Jacobian matrix, with subscripts , , Representing radial coordinates respectively and circumferential coordinates , Find the first derivative.

[0044] Find the inverse of the Jacobian matrix:

[0045]

[0046] Where the superscript -1 denotes the inverse of the matrix, It is the determinant of the Jacobian matrix. , , , , , , , , represents all elements in the inverse of the Jacobian matrix.

[0047] The differential operator expression for the transient seepage problem is derived as follows:

[0048]

[0049] in, Represents the differential operator; , , They represent respectively to x , y , z Find the partial derivative; , , They represent respectively to ξ , η , ζ Find the partial derivative; , , This represents the sub-coefficient vector.

[0050]

[0051] Furthermore, in step S2, using the total energy functional characterizing the system's energy state as the variational basis, the energy functional expression for transient seepage is derived, specifically as follows:

[0052]

[0053] in, Represents the total energy functional. This represents the energy dissipation functional. This represents a functional for storing energy. Indicates water head. Indicates the water storage coefficient. k Indicates the permeability coefficient. Indicates the hydraulic gradient. Indicates frequency, Represents the entire computational domain. dV Indicates volume integral, i Represents an imaginary number, superscript T This represents the transpose of a matrix.

[0054] Discretize the field variables in the SBFEM coordinate system to obtain the interpolated forms of the head and gradient at any point in the problem domain:

[0055]

[0056]

[0057] in, This represents the head interpolation form in the SBFEM coordinate system; The head shape function matrix represents the surface element; , The coordinate transformation matrix is ​​defined as follows:

[0058] ,

[0059] Substitute the above interpolation values ​​into the total energy functional, and then apply the results to the surface coordinates. Integrating over the coordinates transforms the total energy functional expression into a discrete form that depends only on the radial coordinates, yielding the following relationship:

[0060]

[0061] in, The first variation of the total energy functional is represented; , Indicates the starting and ending values ​​of the radial coordinate; This represents the ray direction head vector connecting the nodes of similar surface elements and the nodes of boundary surface elements; Represents the permeability coefficient matrix; This represents the integral over the radial coordinates. , , Represents the transient seepage coefficient matrix; The transient seepage mass matrix is ​​defined as follows:

[0062]

[0063]

[0064]

[0065]

[0066] in, , Represents the circumferential coordinates , integral, This indicates summation over the entire domain.

[0067] By applying the first variational principle, the first derivative of the total energy functional at its stationary point is set to zero. The specific expression is as follows:

[0068]

[0069] in, This is a variational symbol.

[0070] Differential equations were extracted from the weak form equations using the integration by parts method, resulting in the governing equations in the SBFEM coordinate system for the water head.

[0071]

[0072] in, Represents radial coordinates Find the second derivative.

[0073] Furthermore, in step S3, the internal equivalent node flux is obtained based on the principle of virtual work and the weighted residual method, specifically expressed as follows:

[0074]

[0075] in, This represents the internal equivalent node throughput.

[0076] Based on the inverse relationship between internal and external equivalent node fluxes, describe the relationship between nodal head and external equivalent node flux:

[0077]

[0078] in, This represents the flux of the external equivalent node. This is the transient seepage matrix.

[0079] The governing equations for the transient seepage matrix are derived as follows:

[0080]

[0081] Introducing dimensionless frequency:

[0082]

[0083] in, It is a dimensionless frequency. r It is the area of ​​the circumferential boundary surface, and the superscript 2 indicates square.

[0084] Introducing a dimensionless seepage matrix, the dimensionless seepage matrix and the transient seepage matrix have the following relationship:

[0085]

[0086] in, It is a dimensionless seepage matrix.

[0087] Find the dimensionless seepage matrix with respect to... and The partial derivative of is obtained as follows:

[0088]

[0089]

[0090] in, for abbreviation, For frequency Find the partial derivative.

[0091] By substituting variables, the original equation is transformed into the transient seepage matrix governing equation in the frequency domain:

[0092]

[0093] in, for The transient seepage coefficient matrix and mass matrix approaching infinity.

[0094] Furthermore, in step S4, the seepage matrix is ​​decomposed into a power series form using the continuous fraction method:

[0095]

[0096] in, express The transient seepage matrix approaching infinity; Indicates about of o Order function; , This is the coefficient matrix in the power series expansion; Represents the normalized frequency. With frequency The following expressions exist:

[0097]

[0098] The coefficient matrix is ​​obtained by solving a general eigenvalue problem. :

[0099]

[0100] in, Represents the eigenvector matrix; This represents the eigenvalue matrix.

[0101] By setting the coefficients of higher-order terms to zero, and solving the generalized eigenvalue problem and the Lyapunov equation, we obtain information about... The expression for the function to be solved:

[0102]

[0103] in, express o Initial value of the first-order function; express o The recursive calculation value of the order function; Indicates about of o +1 order function.

[0104] By expanding the remaining lower-order terms using power series and recursively calculating the subsequent coefficient matrices, and setting a termination condition to truncate the expansion, the expression for the transient seepage matrix is ​​obtained as follows:

[0105]

[0106] in, This represents the initial value of a first-order function; Represents the recursive calculated value of a first-order function; Represents the initial value of a second-order function; This represents the recursive calculated value of a second-order function.

[0107] Further, in step S5, let The governing equations for the steady-state seepage problem are obtained as follows:

[0108]

[0109] in It is the steady-state seepage matrix.

[0110] Formulate the governing equations as an eigenvalue problem:

[0111]

[0112] in, The eigenvector matrix, For the eigenvalue matrix, Let be a Hamiltonian matrix. The above equation can be explicitly expressed as:

[0113]

[0114] in, It is the identity matrix; , , , The eigenvector submatrix; It is the eigenvalue submatrix.

[0115] The steady-state seepage matrix is ​​obtained by solving the Schur decomposition:

[0116]

[0117] The relationship between flux and head at the boundary of the problem domain is described by the following expression:

[0118]

[0119] in, for abbreviation; for abbreviation; Let be the mass matrix of the boundary surface.

[0120] Then, by applying the inverse Fourier transform to the above equation, we obtain the time-domain expression for the head:

[0121]

[0122] in, For the water head in the time domain, This refers to the general form of water head.

[0123] The improved refined time history integration method MPTSIM is used to accurately calculate the head value at any time step using matrix exponentiation, as shown in the following expression:

[0124]

[0125]

[0126] in, t For any historical point in time, , The first k and k +1 corresponds to the time point in time; , The first k and k +1 time step corresponds to the head at the time point; This is the initial head; The state matrix; For external source terms; dt exp is the integral over time; exp is the exponential function.

[0127] By executing a recursive loop:

[0128]

[0129] in, This represents the number of iterations. N This represents the total number of iterations. Expand the exponential matrix using Taylor expansion.

[0130] Obtain the head values ​​at each time step on the boundary surface:

[0131]

[0132] in, It is an exponential matrix; For time step; n Expand the Taylor series; 2 n- A first-order negligible term; and These are the weights and Gaussian points of the Gaussian integral, respectively; From 1 to l Summation, l This indicates the number of integration points.

[0133] Define the head as a function of the radial coordinate:

[0134]

[0135] in, This indicates the point. The head is a function of the radial coordinate; This indicates the point. The transient seepage matrix at the location; This indicates the point. External equivalent node flux at the location.

[0136] The fourth-order Runge-Kutta method is used to solve for the head distribution in a semi-infinite domain. The specific expression is as follows:

[0137]

[0138] in, Indicates the first j obtained by step Water head; Indicates the first j -1 step obtained Water head; This is the iteration step size; , , , Corresponding to points respectively , , , The gradient, specifically expressed as:

[0139]

[0140] in, Indicates the first j walking The head is a function of the radial coordinate; , , Indicates the first j -1 step corresponds to point , , The transient seepage matrix at the location; , , Indicates the first j -1 step corresponds to point , , External equivalent node flux at the location.

[0141] Furthermore, to verify the efficiency and accuracy of the present invention, a transient seepage characteristic analysis in a three-dimensional triangular channel layered semi-infinite domain is presented in the first embodiment.

[0142] like Figure 3As shown, the computational domain consists of a 100-meter-long river channel with a double-layer foundation. The channel is 30 meters wide at the top and 20 meters deep. The first layer has a hydraulic conductivity and a storage coefficient of 1×10⁻⁶. -4 m / s and 1×10 -5 / m, and the semi-infinite field layers are 1×10. -7 m / s and 1×10 -3 / m. Two independent transient scenarios (water level rise and fall) were simulated within a 500-day period. To characterize the transient behavior, three representative observation points were selected: point A is located inside layer 1 (50, 20); point B is located at the interface between the two layers (-15, 10); and point C is located inside the semi-infinite domain (0, -20). Figure 4 and Figure 8 This paper demonstrates a discretization scheme for unsteady seepage elements in a layered semi-infinite domain of a three-dimensional triangular channel. The computational model uses three 19-node SBFEM elements, and after 2000 iterations, the results are compared and validated with those from GeoStudio SEEP / W finite element simulations. Figure 5 As shown, the calculation results of this invention (dashed line) are in excellent agreement with the SEEP / W curve (solid line), intuitively demonstrating that the deviation is negligible throughout the 500-day simulation. Furthermore, this invention achieves these results with significantly higher computational efficiency: it uses only two 19-node SBFEM elements and 2000 iterations, while SEEP / W requires 37568 tetrahedral elements. This results in a significant reduction of approximately 60% in the number of degrees of freedom, thus demonstrating the superior computational efficiency of this invention.

[0143] Furthermore, to verify the applicability of the present invention in complex systems, a second embodiment is presented – transient seepage characteristics analysis of a three-dimensional irregular homogeneous dam semi-infinite domain foundation.

[0144] like Figure 6 As shown, the calculation parameters for this embodiment are as follows: the upstream water level area is 6m long and EF=2m wide; the dam embedment depth is 3m; the dam length flush with the ground is FG=4m; the length of the bottommost dam is AB=1m; the initial water levels at the upstream and downstream boundaries are 1.0m and 0.0m, respectively; the upstream water level gradually rises to 4.0m over time. Furthermore, the hydraulic conductivity and water storage coefficient of the semi-infinite domain foundation are 1.02×10⁻⁶. -3 m / min and 1×10 -3 / m. Figure 7 and Figure 9This paper presents a discretization scheme for unsteady seepage elements in a three-dimensional irregular homogeneous dam. The computational model uses three 16-node SBFEM elements and is solved through 2000 iterations. For comparative analysis, two points (A and B) are selected for comparison, and the calculated water heads at these points are compared with those obtained by the finite element method (FEM) and by Yan et al. (2025) using the polyhedral scaled boundary finite element method (PSBFEM).

[0145] As shown in Table 1, in terms of accuracy, the results of this invention show excellent agreement with finite element numerical simulations, with an average error of less than 2%, thus verifying its correctness in solving complex seepage problems. Furthermore, the results are very close to those of PSBFEM, with a difference within 1%, further confirming its reliability. Notably, this invention has a significant advantage in computational efficiency. The proposed method requires a total CPU time of 16.15 s, compared to 36.90 s for PSBFEM and 45.98 s for FEM, representing reductions in computation time of 56% and 69%, respectively. These results highlight the superior performance of this invention in achieving high accuracy and computational efficiency in solving transient seepage problems in semi-infinite domains with complex boundaries.

[0146] Table 1: Comparison of head and CPU time calculated by this invention and other methods

[0147] .

Claims

1. A time-domain integral solution method for transient seepage in a dam containing a three-dimensional heterogeneous semi-infinite domain, characterized in that, The method includes the following steps: S1. Construct similar surfaces parallel to the boundary surface, and define the scaled boundary finite element SBFEM coordinate system: circumferential coordinates η and ζ are parallel to the boundary surface, and radial coordinates ξ are perpendicular to the boundary surface. Give the coordinate expression of any point in the problem domain. Introduce shape functions to discretize the points on the boundary surface and similar surfaces, and obtain the mapping of any point in the problem domain in the surface element to the SBFEM coordinate system. Realize the transformation between the Cartesian coordinate system and the SBFEM coordinate system through the Jacobian matrix of coordinate transformation. Further obtain the inverse matrix of the Jacobian matrix and derive the differential operator expression of the transient seepage problem. S2. Using the total energy functional characterizing the system's energy state as the variational basis, derive the energy functional expression for transient seepage. Discretize the field variables in the SBFEM coordinate system to obtain the interpolated forms of the head and gradient at any point in the problem domain. Substitute these interpolations into the total energy functional and then apply them to the circumferential coordinate system. By integrating the equations, the total energy functional expression is transformed into a discrete form that depends only on the radial coordinates. By applying the first variational principle, the first derivative of the stationary point of the total energy functional is made zero. The differential equations are extracted from the weak form equations using the integration by parts method, and the control equations in the SBFEM coordinate system with respect to the head are obtained. S3. Based on the principle of virtual work and the weighted residual method, the internal equivalent node flux is obtained. According to the inverse relationship between the internal and external equivalent node fluxes, the relationship between the node head and the external equivalent node flux is described. Furthermore, the governing equations of the transient seepage matrix are derived. Dimensionless frequency and dimensionless seepage matrix are introduced, and the radial coordinates of the dimensionless seepage matrix are calculated. and frequency The partial derivatives of the equation are used to transform the original equation into the transient seepage matrix control equation in the frequency domain by variable substitution. S4. Apply the continuous fraction method to... The transient seepage matrix approaching infinity is decomposed into a power series form, and the coefficient matrix is ​​obtained by solving a general eigenvalue problem. By setting the coefficients of higher-order terms to zero, and solving the generalized eigenvalue problem and the Lyapunov equation, we obtain the normalized frequency. The expression of the function to be solved is obtained, and then the remaining lower-order terms are expanded by power series. The subsequent coefficient matrix is ​​calculated recursively, and the expansion is truncated by setting a termination condition to obtain the expression of the transient seepage matrix. S5, set frequency The governing equations for the steady-state seepage problem are obtained, and these equations are expressed as eigenvalue problems. The steady-state seepage matrix is ​​obtained by Schur decomposition, which describes the relationship between flux and head on the boundary surface of the problem domain. The inverse Fourier transform is applied to obtain the head expression in the time domain. The improved refined time history integral method MPTSIM is used, and the head value at any time step is calculated by matrix exponentiation. By executing a recursive loop, the head value at each time step on the boundary surface is obtained. A function of head with respect to the radial coordinate is defined, and the head distribution in the semi-infinite domain is solved using the fourth-order Runge-Kutta method.

2. The time-domain integral solution method for transient seepage in a dam containing a three-dimensional heterogeneous semi-infinite domain, as described in claim 1, is characterized in that... In step S1, a similar surface parallel to the boundary surface is constructed, and a scaled boundary finite element (SBFEM) coordinate system is defined: circumferential coordinates η and ζ are parallel to the boundary surface, and radial coordinates ξ are perpendicular to the boundary surface. The coordinate expression for any point in the problem domain is given as follows: ； in, , and Represents the coordinates of a point in the problem domain. , and Represents the coordinates of points on similar surfaces. , and Represents the coordinates of points on the boundary surface. Represents the radial coordinates of the SBFEM coordinate system; We introduce shape functions to discretize points on the boundary surface and similar surfaces: ； ； in, The shape function representing the surface element. , , These are the node coordinate vectors of points on similar surfaces within the surface element; , , Represents the nodal coordinate vector of a point on the boundary surface within a surface element; Find the mapping of any point in the problem domain within the surface element to the SBFEM coordinate system: ； in, , , Represents the nodal coordinate vector of a point on the boundary surface within a surface element; The transformation between the Cartesian coordinate system (x, y, z) and the SBFEM coordinate system (ξ, η, ζ) is achieved using the Jacobian matrix for coordinate transformation. The expression for the Jacobian matrix is ​​as follows: ； in, Represents the Jacobian matrix, with subscripts , , Representing radial coordinates respectively and circumferential coordinates , Find the first derivative; Find the inverse of the Jacobian matrix: ； Where the superscript -1 denotes the inverse of the matrix, It is the determinant of the Jacobian matrix. , , , , , , , , For all elements in the inverse of the Jacobian matrix; The differential operator expression for the transient seepage problem is derived as follows: ； in, Represents the differential operator; , , Let x, y, and z represent the partial derivatives with respect to x, y, and z, respectively. , , Let denote the partial derivatives with respect to ξ, η, and ζ, respectively; , , Represents the sub-coefficient vector; 。 3. The time-domain integral solution method for transient seepage in a dam containing a three-dimensional heterogeneous semi-infinite domain, as described in claim 2, is characterized in that... In step S2, the energy functional describing the system's energy state is used as the variational basis to derive the energy functional expression for transient seepage. The specific expression is as follows: ； in, Represents the total energy functional. This represents the energy dissipation functional. This represents the energy storage functional. Indicates water head. The storage coefficient is represented by k, and the permeability coefficient is represented by k. Indicates the hydraulic gradient. Indicates frequency, dV represents the entire computational domain, i represents the volume integral, and the superscript T represents the transpose of the matrix. Discretize the field variables in the SBFEM coordinate system to obtain the interpolated forms of the head and gradient at any point in the problem domain: ； ； in, This represents the head interpolation form in the SBFEM coordinate system; The head shape function matrix represents the surface element; , The coordinate transformation matrix is ​​defined as follows: ， ； Substitute the above interpolation values ​​into the total energy functional, and then apply the results to the surface coordinates. Integrating over the coordinates transforms the total energy functional expression into a discrete form that depends only on the radial coordinates, yielding the following relationship: ； in, The first variation of the total energy functional is represented; , Indicates the starting and ending values ​​of the radial coordinate; This represents the ray direction head vector connecting the nodes of similar surface elements and the nodes of boundary surface elements; Represents the permeability coefficient matrix; Indicates the integral over the radial coordinates; , , Represents the transient seepage coefficient matrix; The transient seepage mass matrix is ​​defined as follows: ； ； ； ； in, , Represents the circumferential coordinates , integral, This indicates summation over the entire field; By applying the first variational principle, the first derivative of the total energy functional at its stationary point is set to zero. The specific expression is as follows: ； in, For variational notation; Differential equations were extracted from the weak form equations using the integration by parts method, and the governing equations in the SBFEM coordinate system with respect to the water head were obtained. ； in, Represents radial coordinates Find the second derivative.

4. The time-domain integral solution method for transient seepage in a dam containing a three-dimensional heterogeneous semi-infinite domain, as described in claim 3, is characterized in that... In step S3, the internal equivalent node flux is obtained based on the principle of virtual work and the weighted residual method, and the specific expression is as follows: ； in, This represents the internal equivalent node throughput. Based on the inverse relationship between internal and external equivalent node fluxes, describe the relationship between nodal head and external equivalent node flux: ； in, This represents the flux of the external equivalent node. This is the transient seepage matrix; The governing equations for the transient seepage matrix are derived as follows: ； Introducing dimensionless frequency: ； in, It is a dimensionless frequency, r is the area of ​​the circumferential boundary surface, and the superscript 2 indicates the square. Introducing a dimensionless seepage matrix, the dimensionless seepage matrix and the transient seepage matrix have the following relationship: ； in, The matrix is ​​a dimensionless seepage matrix; Find the dimensionless seepage matrix with respect to... and The partial derivative of is given by the following expression: ； ； in, for abbreviation, For frequency Find the partial derivative; By substituting variables, the original equation is transformed into the transient seepage matrix governing equation in the frequency domain: ； in, for The transient seepage coefficient matrix approaching infinity; for The mass matrix that approaches infinity.

5. The time-domain integral solution method for transient seepage in a dam containing a three-dimensional heterogeneous semi-infinite domain, as described in claim 4, is characterized in that... In step S4, the seepage matrix is ​​decomposed into a power series form using the continuous fraction method: ； in, express The transient seepage matrix approaching infinity; Indicates about an o-order function; , This is the coefficient matrix in the power series expansion; Represents the normalized frequency. With frequency The following expressions exist: ； The coefficient matrix is ​​obtained by solving a general eigenvalue problem. : ； in, Represents the eigenvector matrix; Represents the eigenvalue matrix; By setting the coefficients of higher-order terms to zero, and solving the generalized eigenvalue problem and the Lyapunov equation, we obtain information about... The expression for the function to be solved: ； in, This represents the initial value of a function of order 0; This represents the recursive calculated value of an o-order function; Indicates about A function of order o+1; By expanding the remaining lower-order terms using power series and recursively calculating the subsequent coefficient matrices, and setting a termination condition to truncate the expansion, the expression for the transient seepage matrix is ​​obtained as follows: ； in, This represents the initial value of a first-order function; Represents the recursive calculated value of a first-order function; This represents the initial value of a second-order function; This represents the recursive calculated value of a second-order function.

6. The time-domain integral solution method for transient seepage in a dam containing a three-dimensional heterogeneous semi-infinite domain, as described in claim 5, is characterized in that... In step S5, let The governing equations for the steady-state seepage problem are obtained as follows: ； in It is the steady-state seepage matrix; Formulate the governing equations as an eigenvalue problem: ； in, The eigenvector matrix, For the eigenvalue matrix, Let be a Hamiltonian matrix; the above equation is explicitly expressed as: ； in, It is the identity matrix; , , , The eigenvector submatrix; It is an eigenvalue submatrix; The steady-state seepage matrix is ​​obtained by solving the Schur decomposition: ； The relationship between flux and head at the boundary of the problem domain is described by the following expression: ； in, for abbreviation; for abbreviation; The mass matrix of the boundary surface; Then, by applying the inverse Fourier transform to the above equation, we obtain the time-domain expression for the head: ； in, For the water head in the time domain, This refers to the general form of water head; The improved refined time history integration method MPTSIM is used to accurately calculate the head value at any time step using matrix exponentiation, as shown in the following expression: ； ； Where t is any historical time point, , These are the time points corresponding to the kth and k+1th time steps, respectively; , These are the water head at the time points corresponding to the kth and k+1th time steps, respectively; This is the initial head; The state matrix; dt is the external source term; exp is the integral over time; By executing a recursive loop: ； in, The iteration number is denoted as N; the total number of iterations is denoted as N. Expand the exponent matrix using Taylor series; Obtain the head values ​​at each time step on the boundary surface: ； in, It is an exponential matrix; n is the time step; n is the Taylor series. It is a high-order negligible term of degree 2n-1; and These are the weights and Gaussian points of the Gaussian integral, respectively; The summation is performed from 1 to l, where l represents the number of integration points; Define the head as a function of the radial coordinate: ； in, This indicates the point. The head is a function of the radial coordinate; This indicates the point. The transient seepage matrix at the location; This indicates the point. External equivalent node flux at the location; The fourth-order Runge-Kutta method is used to solve for the head distribution in a semi-infinite domain. The specific expression is as follows: ； in, This indicates the result obtained in step j. Water head; This indicates the result obtained in step j-1. Water head; This is the iteration step size; , , , Corresponding to points respectively , , , The gradient, specifically expressed as: ； in, Indicates the time at step j The head is a function of the radial coordinate; , , This indicates the point corresponding to the (j-1)th step. , , The transient seepage matrix at the location; , , This indicates the point corresponding to the (j-1)th step. , , External equivalent node flux at the location.

Images