A hybrid element local orthogonal basis random field discretization method for irregular regions

Abstract

The application discloses a hybrid element local orthogonal basis random field discretization method considering irregular areas, belongs to the technical field of geotechnical engineering random field modeling, and comprises the following steps: determining a target area; performing hybrid mesh division on a random field discretization area to obtain a set of rectangular elements and a set of triangular elements; establishing Karhunen-Loeve expansion expression and integral eigenvalue problems corresponding to the random field discretization area; constructing local orthogonal basis functions on the rectangular elements and the triangular elements, and calculating local basis function values at integral points of each element; discretizing the integral eigenvalue problems by using a Galerkin method, respectively constructing local mass matrices and local correlation matrices by using the local orthogonal basis functions of the rectangular elements and the triangular elements, assembling the local mass matrices and the local correlation matrices to form global mass matrices and global correlation matrices, establishing and solving a generalized eigenvalue equation, calculating random field discretization values in the target area, and completing random field discretization in the irregular area.

Description

Technical Field

[0001] This invention belongs to the field of geotechnical engineering random field modeling technology, specifically involving a method for discretizing hybrid unit local orthogonal basis random fields considering irregular regions. Background Technology

[0002] During the formation, deposition, weathering, and subsequent stress evolution of soil and rock masses, their physical and mechanical parameters typically exhibit significant spatial variability due to the influence of various factors such as material composition, genetic environment, groundwater conditions, and external load history. Traditional deterministic analysis methods that treat soil parameters as spatially homogeneous constants fail to accurately reflect the random distribution characteristics of soil and rock parameters in space. Therefore, using random field theory to describe the spatial variability of soil and rock parameters has become an important research direction in reliability analysis and stochastic finite element analysis in geotechnical engineering. Since continuous random fields consist of an infinite number of correlated random variables, they cannot be directly used for numerical calculations. Therefore, they must be discretized into a finite number of random variables; this process is known as random field discretization.

[0003] Existing random field discretization methods mainly include the midpoint method, local averaging method, shape function method, optimal linear estimation method, and series expansion method. Among them, the Karhunen–Loève expansion method is widely used for random field discretization of soil parameters and subsequent reliability analysis because it has superior global mean square representation capability for a given number of random variables, and the number of expansion terms usually does not depend on the number of physical model grids (Tan et al., 2020; Han and Huang, 2023; Lin et al., 2022). Furthermore, related research also shows that the Karhunen–Loève expansion method can be effectively combined with FORM, RSM, stochastic finite element method, and probabilistic analysis methods in geotechnical engineering, thus achieving a good balance between computational accuracy and efficiency.

[0004] The key to the Karhunen–Loève expansion lies in solving the second kind of Fredholm integral eigenvalue problem with the autocorrelation function as its core. For a few special covariance functions and simple regular regions, this integral eigenvalue problem has analytical solutions; however, for general regions, complex boundaries, or general forms of autocorrelation functions, analytical solutions are often difficult to obtain, and numerical methods are usually required for approximate solutions (Pranesh and Ghosh, 2016; Basmaji et al., 2023; Han and Huang, 2023). Currently, numerical methods for solving this integral eigenvalue problem mainly include the Galerkin method, interpolation methods, hierarchical matrix methods, fast multipole methods, spectral element methods, and embedding domain or bounding box methods. Among these, the Galerkin method, due to its good approximation properties and ease of integration with the finite element analysis framework, has become one of the most commonly used methods in random field discrete research (Basmaji et al., 2022; Han and Huang, 2023; Basmaji et al., 2023).

[0005] However, the traditional Galerkin method still has significant limitations in the discretization of multidimensional random fields. On the one hand, when using global basis functions, it is usually necessary to improve the approximation accuracy by increasing the polynomial order. This can easily lead to numerical instability when the correlation length is small or the spatial dimension is high, and significantly increase the cost of matrix construction and eigenvalue solving. On the other hand, for two-dimensional and three-dimensional random fields, especially in the case of non-separable autocorrelation functions, Galerkin discretization often involves high multiple integrals and high-dimensional generalized eigenvalue problems, resulting in high computational costs. Furthermore, as the dimension of the random field increases, the number of random variables and the matrix size also increase rapidly, further exacerbating the "curse of dimensionality" problem (Pranesh and Ghosh, 2016; Lin et al., 2025).

[0006] To address the aforementioned issues, existing research has improved random field discretization methods from different perspectives. For example, Lin et al. (2022) proposed the Jacobi–Lagrange–Galerkin method, which combines Jacobi orthogonal polynomials, Lagrange interpolation polynomials, and Gauss–Legendre integrals to transform high-multiplication integrals in multidimensional random field discretization into matrix multiplication, thereby improving the efficiency of multidimensional random field discretization. Gu et al. (2023) proposed a random field discretization and slope reliability analysis method based on Chebyshev-Galerkin-KL expansion, demonstrating that constructing a Galerkin-KL discretization framework using orthogonal polynomials helps improve the efficiency of slope reliability analysis. Meanwhile, Han and Huang (2023) proposed a method for solving the eigenvalue problem of KL expansion integrals based on MLS, LS, FEM, and XFEM interpolation basis functions, indicating that interpolation methods also have certain applicability in regular and partially irregular domains.

[0007] Furthermore, to adapt to complex geometric regions, some studies employ the bounding box or embedded domain approach, embedding the actual complex region within a regular rectangular or cuboid bounding domain, constructing a Karhunen-Loève expansion on the bounding domain, and then extracting the random field within the target region. Basmaji et al. (2023) systematically studied the KL expansion method based on the analytical solution of the bounding box domain, pointing out that although this type of method is simple and robust, its error variance is often more significant near the boundary of the bounding domain in the target region, and the scale of the bounding domain affects the eigenvalue decay and variance accuracy. Pranesh and Ghosh (2016) discussed the possibility of reducing the random dimension through domain decomposition from the perspective of the relationship between region size and KL eigenvalue decay, indicating that changes in region scale significantly affect the truncation efficiency of the KL expansion. The above studies show that although the regular bounding domain method can avoid the difficulty of directly solving complex regions, it may still introduce additional errors or reduce truncation efficiency near complex boundaries.

[0008] Regarding local basis function methods, Basmaji et al. (2022) proposed a Galerkin method based on discontinuous Legendre polynomials. This method constructs Legendre basis functions within a local rectangular cell domain and utilizes their orthogonality to improve the efficiency of matrix construction and assembly, thereby mitigating to some extent the problems of large matrix size and high computational cost in the traditional global Galerkin method. However, this method mainly relies on rectangular local cells for discretization, and its computational domain is not entirely consistent with the actual irregular physical region, thus still having certain limitations in adaptability to complex boundaries. Lin et al. (2025) further proposed a dual-coupled discontinuous Galerkin method for large-scale three-dimensional non-separable random fields. By dividing the random field domain into multiple internally non-overlapping subdomains, the numerical instability caused by ultra-high-order polynomials is reduced. Furthermore, by combining Lagrange interpolation, Gauss-Legendre integration, and Krylov-Schur eigenvalue solving strategies, the computational cost of discretizing non-separable three-dimensional random fields is reduced. These studies demonstrate that the Galerkin framework, based on locally discontinuous orthogonal basis functions, has significant advantages in handling high-dimensional problems and non-separable autocorrelation functions.

[0009] On the other hand, existing patents have disclosed a "method for discretizing random fields using shape function KL expansion considering irregular regions of slope soil." This method achieves random field discretization of irregular slope regions by establishing a slope model, dividing the data into discrete and shape function grids, calculating the element matrix, autocorrelation matrix, and global matrix, and solving for eigenvalues ​​and eigenfunctions. While this existing technology can handle irregular region problems to some extent, its disclosed technical approach is mainly based on shape function interpolation and quadrilateral element grid construction. There is still room for improvement in simultaneously considering computational efficiency in regular regions and geometric adaptability to complex boundary regions.

[0010] In summary, while existing random field discretization methods have made significant progress in areas such as regular regions, global basis functions, interpolation solutions, and the Galerkin method for local discontinuities, a suitable random field discretization method for irregular engineering regions is still lacking. This method can balance efficient discretization of regular regions, adaptation to complex boundary geometry, and handling of non-separable autocorrelation functions within a unified framework. In particular, for a class of regions that includes both regular regions suitable for efficient processing with rectangular elements and irregular regions requiring adaptation to complex boundaries with triangular elements, constructing local orthogonal basis functions on different types of elements within a unified Karhunen–Loève / Galerkin framework and solving the generalized eigenvalue problem through unified matrix assembly could potentially improve adaptability and computational efficiency for complex geotechnical engineering regions while maintaining discretization accuracy. Summary of the Invention

[0011] To address the aforementioned shortcomings in existing technologies, the hybrid unit local orthogonal basis random field discretization method considering irregular regions provided by this invention solves the problem that existing random field discretization methods struggle to simultaneously achieve computational efficiency in regular regions, adaptability to complex boundary geometry, and the ability to handle non-separable autocorrelation functions in irregular regions.

[0012] To achieve the aforementioned objectives, the present invention employs the following technical solution: a local orthogonal basis random field discretization method considering irregular regions of hybrid units, comprising the following steps: S1. Establish a geometric model of the discrete region of the random field and determine the target region; S2. Perform hybrid meshing on the discrete region of the random field to obtain a set of rectangular cells and a set of triangular cells, thus generating a hybrid mesh; S3. Based on the statistical parameters and autocorrelation function model of the random field to be discretized, establish the Karhunen-Loève expansion expression, integral eigenvalue problem, and Galerkin discrete solution framework based on the local orthogonal basis of the hybrid grid corresponding to the discrete region of the random field. S4. Construct local orthogonal basis functions on rectangular and triangular elements, and calculate the local basis function values ​​at the integration points of each element; S5. The Galerkin method is used to discretize the integral eigenvalue problem. Local mass matrix and local correlation matrix are constructed using the local orthogonal basis functions of rectangular and triangular elements, respectively, and then assembled to form global mass matrix and global correlation matrix. S6. Establish a generalized eigenvalue equation based on the global quality matrix and the global correlation matrix, and solve it to obtain the effective eigenvalues ​​and corresponding eigenvectors, thus obtaining the approximate eigenfunction. S7. Based on effective eigenvalues ​​and approximate eigenfunctions, construct a finite-term discrete expression for the normal random field, and calculate the discrete value of the random field within the target region to complete the discretization of the random field within the irregular region.

[0013] Furthermore: S2 specifically refers to: The discrete region of the random field is divided into rectangular units using quadtrees to generate a set of rectangular units. The remaining regions not covered by the rectangular units and the regions near irregular boundaries are divided into triangular units to generate a set of triangular units. The set of rectangular units and the set of triangular units together cover the discrete region of the random field.

[0014] Furthermore, S3 includes the following sub-steps: S31. Determine the random field to be discrete: The random field to be discretized is determined based on the discrete region of the random field, and it is expressed in the Karhunen–Loève expansion form: In the formula, For the discrete random field, The first spatial location point within the target area. For the target area, For sample space variables, For the sample space, It is a mean function. It is a function of standard deviation. For the first i 1 eigenvalue, For the first i One characteristic function, For the first i One random variable; S32. Establish a random field autocorrelation function model: Based on the spatial correlation characteristics of the random field to be discretized, the autocorrelation function of the random field is established. ,in, Indicates the second spatial location point within the target area; S33. Problem of establishing integral eigenvalues: Based on the autocorrelation function of the random field We establish the second type of Fredholm integral eigenvalue problem, the specific expression of which is: S34. Establish a Galerkin solution approach based on a hybrid mesh local orthogonal basis: The Galerkin discretization method based on local orthogonal basis of hybrid mesh is used to numerically solve the integral eigenvalue problem. Specifically, based on the hybrid mesh generated by S2, the eigenfunctions in the integral eigenvalue problem are approximated as a linear combination of global piecewise basis functions. For rectangular elements, the Legendre local orthogonal basis function is subsequently used; for triangular elements, the KD local orthogonal basis function is subsequently used.

[0015] Furthermore, S4 includes the following sub-steps: S41. Determine the element type and reference element: The elements in the hybrid mesh are identified, and corresponding reference elements are created according to the element type. When the element is a rectangular element, a reference rectangular element is used. When the element is a triangular element, a reference triangular element is used. ,in, , For local coordinates within the reference triangle element; S42. Construct Legendre local orthogonal basis functions controlled by the total degree on the rectangular element, generating local orthogonal basis functions on the rectangular element: No. Two-dimensional local orthogonal basis functions on rectangular units The specific expression is: In the formula, and For the first One-dimensional local Legendre basis functions on rectangular units and These represent the horizontal and vertical coordinate components in the physical coordinate system of the rectangular unit, respectively. p and indicators q satisfy , For the first The order of the local basis function truncation of each rectangular cell is used to limit the two-dimensional index pair. p , q The highest total number of times; S43. Construct KD local orthogonal basis functions controlled by the total degree in the triangular element to generate local orthogonal basis functions on the triangular element: No. Two-dimensional KD local orthogonal basis functions on triangular units The specific expression is: In the formula, For Jacobi polynomials in parameters The form of time corresponds to the Legendre polynomial. For parameters Jacobi polynomials, where, and For the parameters of the Jacobi polynomial, and Points inside the physical triangle unit ( x , y The local coordinates of the reference triangle element obtained through affine inverse mapping, index p and indicators q satisfy , For the first The order of the local basis function truncation of each triangular unit is used to constrain the two-dimensional index pair. p , q The highest total number of times; S44. Calculate the local basis function values ​​at each unit integration point: Based on the local orthogonal basis functions of the constructed rectangular and triangular elements, numerical integration points are selected on the corresponding reference elements, and the local basis function values ​​at each integration point are calculated.

[0016] Furthermore, S5 includes the following sub-steps: S51. Establish the mapping relationship between the local orthogonal basis function and the global orthogonal basis function numbers: For each unit in the hybrid network, establish local orthogonal basis functions and assign them a unified number. According to the unit numbering order, establish a mapping relationship between the local basis function numbers and the global basis function numbers. The mapping relationship is used to write the local quality matrix and the local correlation matrix into the corresponding positions of the global matrix. S52. Construct the local mass matrix: No. e Local mass matrix of each unit The specific expression is: In the formula, Operators that construct diagonal matrices for each component of a given vector, given the elements of the diagonal. For the first e The diagonal element vector of the local mass matrix of each unit cell; In the formula, For the first e Local orthogonal basis function value matrix over each element For the first e The total number of local orthogonal basis functions on each unit. For the first e The total number of physical integration points in each unit. This is the physical integral weight vector; S53. Construct the local correlation matrix: No. e The unit and the first f Local correlation matrix between units The specific expression is: In the formula, For the first f Local orthogonal basis function value matrix over each element It is the transpose symbol. For the first e The unit and the first f The integral point correlation weight matrix between units; S54. Assemble the global quality matrix: Based on the mapping relationship between local basis function numbers and global basis function numbers, the local quality matrix of each element is written into the global quality matrix as a diagonal block to obtain the global quality matrix; S55. Assemble the global correlation matrix: Based on the mapping relationship between local basis function numbers and global basis function numbers, the local correlation matrix between each pair of units is written into the global correlation matrix to obtain the global correlation matrix.

[0017] Furthermore, S6 includes the following sub-steps: S61. Establish the generalized eigenvalue equation: Based on the global quality matrix and the global correlation matrix, establish the generalized eigenvalue equation: In the formula, This is the global correlation matrix. The global quality matrix, Let be the eigenvector matrix, and let its eigenvector be the eigenvector matrix. Listed as the number eigenvectors , It is an eigenvalue diagonal matrix; S62. Solving for eigenvalues ​​and eigenvectors: Numerical solutions are used to solve the generalized eigenvalue equations to obtain eigenvalues ​​and corresponding eigenvectors, where the obtained eigenvalues ​​are denoted as . The corresponding eigenvector is denoted as , N The dimension of the discrete descendant number generalized eigenvalue problem; S63. Sort and filter the feature values: The obtained eigenvalues ​​are sorted in descending order, and the order of the corresponding eigenvectors is adjusted simultaneously. Furthermore, based on the numerical precision requirements, non-positive eigenvalues, minimal eigenvalues, or invalid eigenvalues ​​caused by numerical errors are eliminated, while the valid eigenvalues ​​and corresponding eigenvectors that meet the discretization requirements of the random field are retained. S64. Normalize the feature vectors: The eigenvectors are normalized, and the normalized eigenvectors are denoted as follows: ; S65. Constructing approximate characteristic functions: Calculate the approximate characteristic function based on the normalized eigenvectors. : In the formula, For the first e Indicator functions for each unit, For the firste The area occupied by each unit For the first e The first unit on p Locally orthogonal basis functions The numbering of local basis functions in the normalized eigenvectors ( e , p The corresponding component.

[0018] Furthermore, S7 includes the following sub-steps: S71. Determine the number of truncated terms in the random field expansion: The number of truncated expansion terms used in the discrete random field analysis is determined based on the effective eigenvalue sequence. M ; S72. Generating random variables: Generate a vector of random variables based on the number of terms in the truncated expansion. , Let be the random variable corresponding to the first expansion term. For the first M Random variables corresponding to each expansion term; S73. Constructing a finite-term discrete expression for a normal random field: Substituting the effective eigenvalues ​​and approximate eigenfunctions into the Karhunen–Loève expansion of the random field to be discretized, we obtain the finite-term discrete expression of the random field based on element-partitioned local basis functions. : S74. Calculate the discrete values ​​of the random field at discrete points: In the target area Select a set of discrete points for a random field, and for all discrete points... Point-by-point calculation yields the discrete value vector of the random field. : In the formula, For discrete points The discrete expression of a random field with finite terms. For discrete points The discrete expression of a random field with finite terms. For discrete points The discrete expression of a random field with finite terms. This represents the total number of discrete points.

[0019] The beneficial effects of this invention are as follows: (1) This invention uses a hybrid mesh combining rectangular and triangular elements. Rectangular elements can be used in regular regions to improve discretization efficiency, while triangular elements can be used to fit complex boundaries, so that the computational region is consistent with the actual target region and the geometric adaptability of irregular regions is improved. Legendre local orthogonal basis and KD local orthogonal basis are constructed on rectangular and triangular elements respectively. The orthogonality within the element simplifies the construction of the local mass matrix, and the truncation order of the local basis function can be adjusted according to the element scale, so as to balance discretization accuracy and computational scale. Through unified Galerkin discretization and matrix assembly, the integral eigenvalue problem is transformed into a generalized eigenvalue problem, which is applicable to separable and non-separable autocorrelation functions. It can provide a random field discretization basis for stochastic finite element analysis and reliability analysis of complex geotechnical engineering areas.

[0020] (2) The present invention can take into account both the computational efficiency of regular regions and the geometric adaptability of complex boundary regions. It is applicable to separable and non-separable autocorrelation functions and can provide a random field discretization basis for stochastic finite element analysis and reliability analysis in geotechnical engineering. Attached Figure Description

[0021] Figure 1 This is a flowchart of the local orthogonal basis random field discretization method for hybrid units considering irregular regions according to the present invention.

[0022] Figure 2 This is a schematic diagram of an irregular target area, where, Figure 2 (a) is the slope area. Figure 2 (b) is a region of inclined strata.

[0023] Figure 3 This is a schematic diagram of the hybrid mesh generation result, where, Figure 3 (a) is a mixed grid for the slope area. Figure 3 (b) is a mixed grid for inclined strata.

[0024] Figure 4 This is a schematic diagram of the discrete point distribution of a random field, where, Figure 4 (a) shows the discrete point distribution in the slope area. Figure 4 (b) shows the distribution of discrete points in the inclined strata region.

[0025] Figure 5 This is a schematic diagram illustrating the equivalent scale of the unit and the truncation order allocation of the local basis functions, where... Figure 5 (a) is the slope area. Figure 5 (b) is a region of inclined strata.

[0026] Figure 6 This represents a two-fold realization of the random field of cohesion in the slope region, where... Figure 6 (a) is the first implementation. Figure 6(b) is the second implementation.

[0027] Figure 7 This is a quadratic realization of the random field of friction angle within the slope region, where, Figure 7 (a) is the first implementation. Figure 7 (b) is the second implementation.

[0028] Figure 8 This represents a two-fold realization of the random field of cohesion in the inclined strata region, where... Figure 8 (a) is the first implementation. Figure 8 (b) is the second implementation.

[0029] Figure 9 This represents a second-order realization of the random field of friction angle within the inclined stratum region, where... Figure 9 (a) is the first implementation. Figure 9 (b) is the second implementation. Detailed Implementation

[0030] The specific embodiments of the present invention are described below to enable those skilled in the art to understand the present invention. However, it should be understood that the present invention is not limited to the scope of the specific embodiments. For those skilled in the art, various changes are obvious as long as they are within the spirit and scope of the present invention as defined and determined by the appended claims. All inventions utilizing the concept of the present invention are protected.

[0031] This invention proposes a random field discretization method under the unified Karhunen–Loève / Galerkin framework. This method uses a hybrid grid local orthogonal basis function (LOB) on the rectangular grid and a KD LOB on the triangular grid to divide irregular regions. This results in a hybrid grid LOB discretization framework suitable for complex regions, which improves the geometric adaptability, non-separable autocorrelation function handling capability, and random field discretization efficiency of irregular geotechnical engineering regions.

[0032] Example 1: like Figure 1 As shown, in one embodiment of the present invention, the method for discretizing a hybrid unit local orthogonal basis random field considering irregular regions includes the following steps: S1. Establish a geometric model of the discrete region of the random field and determine the target region; S2. Perform hybrid meshing on the discrete region of the random field to obtain a set of rectangular cells and a set of triangular cells, thus generating a hybrid mesh; S3. Based on the statistical parameters and autocorrelation function model of the random field to be discretized, establish the Karhunen-Loève expansion expression, integral eigenvalue problem, and Galerkin discrete solution framework based on the local orthogonal basis of the hybrid grid corresponding to the discrete region of the random field. S4. Construct local orthogonal basis functions on rectangular and triangular elements, and calculate the local basis function values ​​at the integration points of each element; S5. The Galerkin method is used to discretize the integral eigenvalue problem. Local mass matrix and local correlation matrix are constructed using the local orthogonal basis functions of rectangular and triangular elements, respectively, and then assembled to form global mass matrix and global correlation matrix. S6. Establish a generalized eigenvalue equation based on the global quality matrix and the global correlation matrix, and solve it to obtain the effective eigenvalues ​​and corresponding eigenvectors, thus obtaining the approximate eigenfunction. S7. Based on effective eigenvalues ​​and approximate eigenfunctions, construct a finite-term discrete expression for the normal random field, and calculate the discrete value of the random field within the target region to complete the discretization of the random field within the irregular region.

[0033] S1 is used to build a two-dimensional geometric model of the target region to be analyzed in the MATLAB environment, and represent it as a single polygonal region that can be directly called upon for subsequent hybrid mesh generation, local orthogonal basis function construction, and random field discretization calculations. S1 includes the following steps: S11: Determine the discrete region of the random field: Based on the actual geotechnical engineering problem, the target region to be discretized is determined to be a two-dimensional single polygon region. The target region includes, but is not limited to, the following two types: (1) Geotechnical engineering areas with irregular geometric shapes, such as slope areas with irregular slope surface shapes, slope shoulders or toes with broken lines. (2) Areas in geotechnical engineering where the boundaries of strata are irregularly distributed, such as irregular areas formed in a two-dimensional plane by a certain layer of soil, rock, weak interlayer or other material partitions.

[0034] In this method, only one polygonal region is selected for random field discretization of the target region in a single calculation. When there are multiple strata or multiple partitions in the project, the method of this invention can be implemented independently for each region.

[0035] S12. Establish a two-dimensional plane coordinate system: Establishing a two-dimensional rectangular coordinate system in the MATLAB environment Oxy This is used to describe the geometric location of the target region. x The axis is horizontal. y The axis is vertical, and the origin is the coordinate system. OThe reference point can be specified by the user, the lower left corner of the target area, or other locations convenient for modeling. Preferably, the coordinate unit adopts a uniform length unit format to ensure consistency in subsequent mesh generation, cell area calculation, and related length parameter input.

[0036] S13. Input the coordinates of the boundary points of the discrete region of the random field: In the two-dimensional rectangular coordinate system established in step S12, input the coordinates of the control points of the target area boundary and establish the boundary point coordinate matrix. : .

[0037] In the formula, For the first The coordinates of the boundary points , The number of boundary points for the target region is specified. Boundary points are entered sequentially according to the outer boundary of the target region, and connected end-to-end to form a closed polygon boundary. Preferably, boundary points are numbered counter-clockwise to ensure consistent boundary representation direction and facilitate subsequent region identification and element clipping. Boundary point coordinates can be directly input by the user or imported into MATLAB from external engineering geometric data.

[0038] S14. Construct the target region geometry object: Based on the boundary point coordinate matrix in S13 In MATLAB, a two-dimensional polygonal geometric object corresponding to the target region is constructed to obtain the discrete region of the random field.

[0039] The target area is represented as When the target area has an irregular slope shape, Ω represents the region of the slope soil or slope analysis object in a two-dimensional plane; when the target area has an irregular distribution of stratigraphic boundaries, Ω represents the region of a certain layer of soil, rock, or a certain independent material zone in a two-dimensional plane. In this way, irregular engineering areas from different sources can be uniformly represented as a single polygonal region, thus facilitating subsequent hybrid mesh generation and random field discretization calculations within a unified framework. The two types of irregular target areas selected in this embodiment are as follows: Figure 2 As shown, where, Figure 2 (a) is the slope area. Figure 2 (b) is a region of inclined strata.

[0040] S2 specifically refers to: The discrete region of the random field is divided into rectangular units using quadtrees to generate a set of rectangular units. The remaining regions not covered by the rectangular units and the regions near irregular boundaries are divided into triangular units to generate a set of triangular units. The set of rectangular units and the set of triangular units together cover the discrete region of the random field.

[0041] S2 is used in the MATLAB environment to perform hybrid mesh generation on the discrete region of the random field obtained in S1, generating a hybrid mesh composed of rectangular and triangular elements. The hybrid mesh uses quadtree rectangular elements as the main mesh and triangular elements to fill in the remaining areas. The shape quality of the elements is controlled to avoid excessively large differences in side length, excessively small interior angles, or excessively long and narrow shapes in the triangular elements. In specific implementation, S2 includes the following steps: S21. Set the hybrid mesh generation parameters: In MATLAB, you can set hybrid meshing parameters to control the mesh density and cell shape quality of the target region.

[0042] The hybrid meshing parameters include the quadtree maximum level parameter, rectangular cell aspect ratio control parameter, triangular cell side length difference control parameter, triangular cell shape quality control parameter, and local mesh refinement control parameter. Among them, the local mesh refinement control parameter is used to appropriately refine the triangular mesh in narrow regions or regions with significant geometric changes, while keeping it relatively sparse in regions with gentler geometry, thus allowing for significant differences in mesh size at different locations.

[0043] S22. Generate a quadtree background rectangular grid: Based on the hybrid mesh partitioning parameters set in S21, a minimum bounding rectangle region covering the target area is constructed in MATLAB, and the minimum bounding rectangle region is used as the initial root region for quadtree partitioning.

[0044] A recursive quadtree partitioning is performed on the initial root region: when a rectangular cell intersects the boundary of the target region Ω and its current level has not reached the preset maximum level, the rectangular cell is further divided into four sub-rectangular cells; when a rectangular cell is completely inside the target region Ω, it is retained as a candidate main rectangular cell; when a rectangular cell is completely outside the target region Ω, it is discarded. This yields a quadtree rectangular mesh covering the main part of the target region, denoted as the set of rectangular cells. In the formula A rectangular cell in a quadtree. The total number of rectangular cells to be retained in the quadtree; S23. Determine the geometric relationship between the rectangular unit and the target region: Check each rectangular cell of the quadtree one by one The geometric relationship with the target region Ω includes the background rectangular unit being completely outside the target region, the background rectangular unit being completely inside the target region, and the background rectangular unit intersecting the boundary of the target region.

[0045] Among them, rectangular cells that are completely located inside the target area are given priority as rectangular cells in the hybrid mesh; the regions corresponding to rectangular cells that intersect with the boundary of the target area are used as candidate regions for subsequent triangular mesh completion.

[0046] S24. Generate rectangular elements: When a rectangular cell is completely located inside the target area and its aspect ratio meets the preset shape quality requirements, the rectangular cell is directly retained as a rectangular cell in the hybrid mesh.

[0047] All rectangular elements that meet the conditions are combined into a rectangular element set. In the formula For the first k The area occupied by a rectangular cell that meets the retention criteria. This represents the total number of rectangular units in the rectangular unit set. S25. Generate triangular units: The remaining areas not covered by the rectangular main grid are divided into triangular grids to generate several triangular cells, thereby filling in the uncovered areas near the boundary of the target area and in some local areas.

[0048] During the triangular mesh generation process, the differences in side length and shape quality of the triangular units are controlled to avoid generating triangular units with excessively large differences in side length, excessively small interior angles, or excessively narrow shapes. Furthermore, based on the local mesh refinement control parameters, narrow or locally complex regions are appropriately refined to ensure that the side lengths of the generated triangular units are as uniform as possible, while allowing for significant differences in the size of the triangular mesh in different regions.

[0049] This yields the triangular unit set. In the formula For the first k The area occupied by one effective triangular unit This represents the total number of valid triangular units in the triangular unit set. S26. Form a hybrid mesh and output element information: Set of rectangular units With triangular unit set Merge to form a hybrid grid covering the target area. .

[0050] In MATLAB, the type, vertex coordinates, number, area, and shape quality information of each element are output and used as input for subsequent local orthogonal basis function construction and random field discretization calculations. Preferably, the resulting hybrid mesh satisfies the following characteristics: the main part of the target region is dominated by quadtree rectangular elements, while the remaining regions are dominated by triangular elements; significant differences in mesh size are allowed in different regions; and a hybrid mesh generation effect that balances local appropriate densification with overall density is achieved while ensuring local shape quality. The hybrid mesh generation results for slope areas and inclined strata areas are shown below. Figure 3 As shown.

[0051] for Figure 3 (a) shows the cell grid division of the slope region. Since its boundary is composed of vertical edges, horizontal edges, and slope polygons, quadtree rectangular cells (blue cells) are preferentially retained in more regular locations within the slope, while triangular cells (red cells) are used in the remaining areas, thus obtaining a mixed grid for the slope region. Figure 3 (b) shows the grid division of the inclined strata region. Since its left and right boundaries are vertical boundaries and its upper and lower boundaries are oblique boundaries, quadtree rectangular cells (blue cells) are preferentially retained in the middle of the region where it is not cut by the upper and lower oblique boundaries, while the remaining areas near the upper and lower oblique boundaries are divided into triangular cells (red cells), thus obtaining a mixed grid of the inclined strata region.

[0052] S3 is used in the MATLAB environment to establish the statistical model, autocorrelation function model, and corresponding integral eigenvalue problem of the random field to be discretized, given the discrete region of the random field determined by S1 and the hybrid grid generated by S2. It also clarifies the subsequent numerical solution approach based on the Galerkin method using local orthogonal bases in the hybrid grid, providing a theoretical foundation for the construction of local orthogonal basis functions in S4 and the matrix construction in S5. S3 includes the following steps: S31. Determine the random field to be discrete: The random field to be discretized is determined based on the discrete region of the random field, and it is expressed in the Karhunen–Loève expansion form: In the formula, For the discrete random field, The first spatial location point within the target area. For the target area, For sample space variables, For the sample space, It is a mean function. It is a function of standard deviation. For the first i 1 eigenvalue, For the firsti One characteristic function, For the first i There are several random variables; related studies have shown that random field discretization is usually based on the statistical characteristics of spatial variation parameters, such as probability distribution type, correlation function, fluctuation range and cross-correlation between different parameters to establish random field models.

[0053] S32. Establish a random field autocorrelation function model: Based on the spatial correlation characteristics of the random field to be discretized, the autocorrelation function of the random field is established. ,in, This represents the second spatial location point within the target region; the autocorrelation function describes the spatial correlation of a random field between any two points within the target region. The autocorrelation function can be either separable or non-separable. A non-separable autocorrelation function is one that cannot be expressed as a product of multiple one-dimensional autocorrelation functions. For non-separable autocorrelation functions, the corresponding Fredholm integral eigenvalue problem usually requires numerical approximation.

[0054] This invention does not limit the specific analytical expression of the non-separable autocorrelation function, but only requires that it satisfy the basic conditions required by the Karhunen-Loève expansion, that is, the autocorrelation function is bounded, symmetric and semi-positive definite in the target region Ω, so that it can be used as the kernel function of the integral eigenvalue problem.

[0055] In MATLAB, it is preferable to uniformly represent autocorrelation functions as callable function handles, numerical function objects, or discrete computation modules so that they can be called point by point and matrix assembled at subsequent integration points.

[0056] S33. Problem of establishing integral eigenvalues: Based on the autocorrelation function of the random field Establish the second type of Fredholm integral eigenvalue problem, the specific expression of which is: When autocorrelation function When the target region Ω is bounded, symmetric, and satisfies positive semidefiniteness, it can be expressed as a series expansion of eigenvalues ​​and eigenfunctions according to Mercer's theorem: Thus, the KL expansion of random fields, the autocorrelation function model, and the integral eigenvalue problem form a unified correspondence. For general geometric regions and general forms of autocorrelation functions, analytical solutions are often difficult to obtain; therefore, numerical methods are needed to discretely solve the aforementioned integral eigenvalue problem.

[0057] The KL expansion essentially requires solving the above integral eigenvalue problem and constructing a discrete expression for the random field using the obtained eigenvalues ​​and eigenfunctions. For general geometric regions and general forms of autocorrelation functions, analytical solutions are usually not available, and numerical methods such as Galerkin's method are needed to obtain approximate solutions.

[0058] S34. Establish a Galerkin solution approach based on a hybrid mesh local orthogonal basis: For general geometric regions and autocorrelation functions of a general form, the above-mentioned integral eigenvalue problem is usually difficult to solve analytically, thus requiring numerical methods. Existing research shows that the basic idea of ​​the Galerkin method is to represent the eigenfunction as a linear combination of a set of basis functions, thereby transforming the continuous integral eigenvalue problem into a discrete algebraic eigenvalue problem.

[0059] This invention preferably employs the Galerkin discretization method based on hybrid grid local orthogonal bases to numerically solve the integral eigenvalue problem. This solution approach targets complex and irregular target regions by using a hybrid partitioning method combining rectangular and triangular elements, and constructing suitable local orthogonal basis functions on different element types. Through local basis function expansion, local matrix construction, and global matrix assembly, the integral eigenvalue problem corresponding to the autocorrelation function of the random field is transformed into an algebraic generalized eigenvalue problem, thus forming a hybrid grid local orthogonal base Galerkin discretization framework suitable for discretizing random fields in complex regions.

[0060] Specifically, based on the hybrid mesh generated by S2, the characteristic functions in the integral eigenvalue problem are approximated as a linear combination of global piecewise basis functions. For rectangular elements, Legendre local orthogonal basis functions are subsequently used; for triangular elements, KD local orthogonal basis functions are subsequently used. The local orthogonal basis functions on both types of elements together constitute a hybrid discretization framework suitable for complex geometries. This framework allows local basis functions to be less continuous between elements, making it more suitable for irregular region boundaries, hybrid element partitioning, and implementing local order control at the element scale.

[0061] Compared with the traditional discretization method that uses a single global basis function on regular regions, the hybrid mesh local orthogonal basis Galerkin solution approach adopted in this invention has the following advantages: (1) It can adapt to complex and irregular regional boundaries; (2) More suitable local orthogonal basis functions can be used on rectangular and triangular units respectively; (3) It facilitates the construction of local mass matrices and local correlation matrices by combining local orthogonality; (4) It facilitates adaptive control of the truncation order of local basis functions according to the unit scale, thereby limiting the maximum total degree of local basis functions in each unit; (5) It is beneficial to improve the flexibility and applicability of random field discretization in complex regions while ensuring computational accuracy. Therefore, the role of S34 is not to directly obtain eigenvalues ​​and eigenfunctions, but to clarify the numerical solution route adopted in this invention, namely: on the basis of the hybrid grid, the Galerkin discretization method based on local orthogonal basis is used to discretize and solve the integral eigenvalue problem corresponding to the autocorrelation function of the random field.

[0062] In this embodiment, regular discrete points are arranged at 0.2 m intervals in both the slope area and the inclined stratum area, and boundary discrete points are added along the area boundaries. These discrete points are used for subsequent random field discrete value calculations, and their distribution is as follows: Figure 4 As shown, Figure 4 (a) is a schematic diagram of the distribution of discrete points in the random field of the slope area. Figure 4 (b) is a schematic diagram of the distribution of discrete points in a random field in an inclined stratum region.

[0063] Building upon the established random field model, autocorrelation function, and integral eigenvalue problem in S3, S4 further constructs a local orthogonal basis function system suitable for hybrid grid discretization. Based on the hybrid grid generated in step S2, S4 constructs Legendre local orthogonal basis functions on rectangular elements and KD local orthogonal basis functions on triangular elements. The truncation order of the local basis functions is adaptively determined according to the element size, ensuring that larger grid elements use higher-order local orthogonal basis functions and smaller grid elements use lower-order local orthogonal basis functions. The truncation order limits the maximum total number of two-dimensional index pairs, thus balancing the accuracy of random field discretization with computational efficiency. The local basis functions serve as the foundation for subsequent local matrix construction and global matrix assembly. In this embodiment, the equivalent element size is used as the basis for allocating the truncation order of the local basis functions. The equivalent element size is defined as... ,in, For the first The area of ​​each unit. For rectangular units, when... When <1m, the local basis function truncation order is 2; when 1m≤ When ≤3m, the local basis function truncation order is 3; when When the value is greater than 3m, the local basis function truncation order is set to 5. For triangular elements, when When <0.1m, the local basis function cutoff order is 1; when 0.1m ≤ When ≤0.3m, the local basis function cutoff order is 2; when When the depth is greater than 0.3m, the local basis function cutoff order is set to 3. The unit equivalent scale and local basis function cutoff order allocation results for the slope region and the inclined stratum region are as follows: Figure 5 As shown, Figure 5 (a) is a schematic diagram of the equivalent scale of the slope area unit and the order of the local basis function cutoff. Figure 5 (b) is a schematic diagram of the equivalent scale of inclined stratigraphic units and the order of local basis function truncation.

[0064] S4 includes the following steps: S41. Determine the element type and reference element: The elements in the hybrid mesh are identified, and corresponding reference elements are created according to the element type. When the element is a rectangular element, a reference rectangular element is used. When the element is a triangular element, a reference triangular element is used. ,in, and These are coordinate variables along two local coordinate directions within the reference triangle element, used to describe the position of any point within the reference triangle element.

[0065] In MATLAB, each cell is recorded with its cell type, vertex coordinates, and corresponding reference cell type, so that a mapping relationship between physical cells and reference cells can be established later.

[0066] S42. Construct Legendre local orthogonal basis functions controlled by the total degree on the rectangular element, generating local orthogonal basis functions on the rectangular element: For each rectangular cell in the hybrid mesh, a two-dimensional locally orthogonal basis function is used, formed by the product of two one-dimensional Legendre polynomials. Let the... The rectangular units are In the formula and The first Each vertex of the rectangular unit is x Minimum and maximum values ​​of direction coordinates, and The first Each vertex of the rectangular unit is y Minimum and maximum values ​​of the direction coordinates; Define a one-dimensional Legendre polynomial on the interval [-1, 1] that satisfies the recurrence relation In the formula, , and The number of times are respectivelyp +1、 p and p -1 one-dimensional Legendre polynomial, Let [the variable be] an independent variable on the standard interval [-1, 1]. for x The first-order index of the one-dimensional local Legendre basis functions along the axial direction; and satisfying the orthogonality relation: In the formula, For the number of times m One-dimensional Legendre polynomial, Using Kronecker notation, and leveraging the orthogonality of Legendre polynomials across local element fields, local basis functions can be constructed for each element, and their orthogonality can be directly utilized in subsequent matrix calculations. For the Each rectangular unit is located in x Axial direction and y Introducing standardized local coordinate transformations along the axial direction: In the formula, and The first e One rectangular unit in x Minimum and maximum coordinates in the direction, and The first e One rectangular unit in y Minimum and maximum coordinates in the direction, and Physical coordinates x and y Standardized local coordinates obtained through linear transformation The above linear transformation maps the physical unit interval to the standard interval [-1, 1] respectively: when x = x min hour, ;when x = x max hour, ; y The same applies to the direction. This transformation includes both translation and scaling.

[0067] Therefore, the first The one-dimensional local Legendre basis functions on a rectangular unit cell can be written as: In the formula, and The first e One rectangular unit in x direction and y One-dimensional local Legendre basis functions in the direction, and The degree of each interval is given by the standard interval [-1, 1]. p and q Legendre polynomials, p for x First-order index of one-dimensional local Legendre basis functions along the axis. q for y The first-order index of the one-dimensional local Legendre basis functions along the axial direction, corresponding to these local basis functions satisfying the following in their respective physical domains: in, m for x The second-order index of the one-dimensional local Legendre basis functions along the axis. n for y The second-order index of the one-dimensional local Legendre basis functions along the axis.

[0068] Therefore, the first... Two-dimensional local orthogonal basis functions on the nth rectangular unit, the th Two-dimensional local orthogonal basis functions on rectangular units The specific expression is: In the formula, and For the first One-dimensional local Legendre basis functions on rectangular units Let be the horizontal physical coordinates of any point within this rectangular unit. The vertical physical coordinates of any point within this rectangular unit, index p and indicators q satisfy , For the first The order of the local basis function truncation of each rectangular cell is used to limit the two-dimensional index pair. p , qThe highest total degree of ); the above expression shows that the two-dimensional local orthogonal basis functions on the rectangular unit are constructed by two one-dimensional Legendre polynomial factors through tensor product.

[0069] In this method, the two-dimensional index pairs ( p , q Instead of using the entire set of complete tensor products, a rule of selecting complete tensors layer by layer according to the total number of tensors is adopted, retaining only those that satisfy the condition. The basis function terms. Therefore, the first... The total number of local basis functions on a rectangular unit for: Table 1 Selection Rules for Two-Dimensional Basis Function Indices of Rectangular Elements As shown in Table 1, in MATLAB, it is preferable to first determine the local basis function truncation order based on the rectangular element scale. Then, following the above-mentioned rule of selecting each layer according to the total number of times, two-dimensional index pairs are generated. p , q The set is used to calculate the local basis function values ​​of the two-dimensional Legendre tensor product on the rectangular unit element.

[0070] S43. Construct KD local orthogonal basis functions controlled by the total degree in the triangular element to generate local orthogonal basis functions on the triangular element: For each triangular element in the hybrid mesh, a KD (Koornwinder–Dubiner) local orthogonal basis function on the reference triangle is used. Let the... Triangular unit The coordinates of the three vertices are respectively The corresponding reference triangle unit is taken as... .

[0071] The affine mapping from the reference triangle element to the physical triangle element is written as: In the above mapping, The first term is a translation term, and the other two terms are linear transformation terms; therefore, this mapping includes both translation and linear transformation. Correspondingly, any point within the physical triangle unit ( x , y All of them can be uniquely mapped to local coordinates in the reference triangle element through inverse mapping. .

[0072] set up Let be a Jacobi polynomial defined on the interval [-1, 1]. In this method, the one-dimensional Jacobi type factor used for the KD local orthogonal basis functions is taken as . ,Right now Jacobi polynomials; and ,Right now Jacobi polynomials.

[0073] in, The simplest recurrence relation can be used to write it as follows: In the formula, and Parameters α =0、 β The zeroth and first-degree Jacobi polynomials when = 0 correspond to the zeroth and first-degree Legendre polynomials. , and The number of times are respectively p +1 p and p -1 Jacobi polynomials correspond to Legendre polynomials of the same degree; for Its initial term can be written as: In the formula, and Parameters α =2 p +1 β Jacobi polynomials of degree zero and degree one when = 0; Its recursive relation can be written as: In the formula, , and Parameters α =2 p +1 β =0 times q +1 q and q Jacobi polynomials of -1; To represent the two-dimensional orthogonal basis functions on the triangular unit as a combination of two one-dimensional Jacobi-type factors, a further variable transformation is introduced in the reference triangular unit. ,in, Therefore, the first Two-dimensional KD local orthogonal basis functions on triangular units The specific expression is: In the formula, For Jacobi polynomials in parameters The form of time, The first parameter of the Jacobi polynomial. This is the second parameter of the Jacobi polynomial, corresponding to the Legendre polynomial. For parameters Jacobi polynomial, and These are the first and second local coordinates of the reference triangle element obtained by the affine inverse mapping of the physical point, respectively. p and indicators q satisfy , For the first The order of the local basis function truncation of each triangular unit is used to constrain the two-dimensional index pair. p , q The highest total number of times; The above expression shows that the two-dimensional locally orthogonal basis functions on the triangular unit are composed of a class of Jacobi polynomial factors and corresponding weight factors: one of the factors is , corresponding to the Legendre polynomial form; another factor is , for parameters Jacobi polynomials with weight factors Therefore, the two-dimensional local orthogonal basis functions on the triangular unit are not a simple tensor product of two independent one-dimensional polynomials, but a product expression formed by combining two one-dimensional Jacobi-type polynomial factors and corresponding weight factors based on the coordinate transformation of the reference triangle.

[0074] The KD basis functions on the reference triangle satisfy an orthogonality relation and can be written as: in, and For Kronecker notation, To match the indicators ( p , q The corresponding orthogonal constant.

[0075] Furthermore, in the physical triangle unit Above: in, For the first Jacobian matrix of triangular units, and For physical triangle units interior point x and y coordinates and For physical triangle units interior point x and y coordinates and For physical triangle units interior point The x and y coordinates are given. This shows that the local basis functions on the triangular element remain orthogonal on both the reference element and the physical element.

[0076] For the Given a triangular unit cell, let the order of its local basis function cutoff be... This truncation order is used to limit the two-dimensional index pair ( p , q The highest total number of times ) is calculated. In this method, the two-dimensional index is used for ( p , q The rule of selecting completely layer by layer according to the total number of times is adopted, and only those that meet the requirements are retained. The basis function terms. Therefore, the first... The total number of local basis functions on the triangular units is .

[0077] Table 2 Selection Rules for Two-Dimensional Basis Function Indices of Triangular Units Therefore, the two-dimensional KD local orthogonal basis functions on the triangular unit also embody the idea of ​​constructing two-dimensional orthogonal basis functions by combining one-dimensional orthogonal polynomial factors, and, like the rectangular unit, adopt a local order control method that selects the complete function layer by layer according to the total degree. The difference is that the local orthogonal basis functions on the rectangular unit are constructed by tensor products of one-dimensional Legendre polynomials in two independent directions, while the local orthogonal basis functions on the triangular unit are constructed by combining a class of Jacobi polynomials and corresponding weight factors based on the coordinate transformation of the reference triangle.

[0078] S44. Calculate the local basis function values ​​at each unit integration point: Based on the local orthogonal basis functions of the constructed rectangular and triangular elements, numerical integration points are selected on the corresponding reference elements, and the local basis function values ​​at each integration point are calculated.

[0079] For the local orthogonal basis functions of the rectangular and triangular elements constructed in S42 and S43, numerical integration points are selected on the corresponding reference elements, and the local basis function values ​​at each integration point are calculated for subsequent numerical integration calculations of the local mass matrix and local correlation matrix. For the... Rectangular units In the reference rectangular unit The Gaussian numerical integration method is used to determine the integration points and their corresponding integration weights. Preferably, one-dimensional Gauss-Legendre integration points can be selected in two coordinate directions respectively, and combined to form two-dimensional integration points. Let the one-dimensional integration points and integration weights be... and .

[0080] Then the two-dimensional integral point can be written as The corresponding two-dimensional integral weights are Calculate the local basis function values ​​of the rectangular element at each integration point. Then, using the local coordinate transformation relationship of the rectangular element established in S42, the integration points on the reference rectangular element are mapped to the physical rectangular element. In the process, the coordinates of the corresponding physical integration point are obtained as follows: .

[0081] For the Triangular unit Its corresponding reference triangle unit is On the reference triangular element, the integration points and their corresponding integration weights are determined using the Gaussian numerical integration method. Let the integration points and their corresponding integration weights on the reference triangular element be respectively... .in, The number of integration points on the reference triangle unit.

[0082] Then, at each integration point, the local basis function values ​​of the triangular element are calculated. Then, using the affine mapping relationship of the triangular elements established in S43, the integration points on the reference triangular element are mapped to the physical triangular element. In the process, the coordinates of the corresponding physical integration point are obtained as follows: .

[0083] For each rectangular and triangular element, the physical integral weight information corresponding to each integration point of the element is calculated and saved by combining the corresponding Jacobian determinant and integral weight.

[0084] S5, in the MATLAB environment, is used to construct the local mass matrix and local correlation matrix corresponding to the rectangular and triangular elements based on the local orthogonal basis functions, integration point coordinates, integration weights, Jacobian information, and local basis function value matrices obtained in step S4. These matrices are then further assembled to form the global mass matrix and global correlation matrix, providing input for the establishment and solution of the algebraic generalized eigenvalue problem in step S6. S5 includes the following sub-steps: S51. Establish the mapping relationship between the local orthogonal basis function and the global orthogonal basis function numbers: For each unit in the hybrid network, establish local orthogonal basis functions and assign them a unified number. According to the unit numbering order, establish a mapping relationship between the local basis function numbers and the global basis function numbers. The mapping relationship is used to write the local quality matrix and the local correlation matrix into the corresponding positions of the global matrix. Assume the hybrid mesh has a total of n e Unit 1, the 1st e The total number of local basis functions on each unit is Then the local basis functions on this unit can be uniformly denoted as .

[0085] Among them, when the first e When each unit is a rectangular unit Represent the Legendre local orthogonal basis functions on this rectangular element; when the... e When each unit is a triangular unit This represents the KD local orthogonal basis function on the triangular unit.

[0086] Local basis functions on all elements are uniformly numbered according to element numbering order, and a local basis function numbering system is established. e , i ) and global basis function number k The mapping relationship between them is used to subsequently write the local quality matrix and the local correlation matrix into the corresponding positions of the global matrix.

[0087] For a rectangular element, the number of its local basis functions is determined by the total number of local basis functions corresponding to that rectangular element in step S4. The number of local basis functions for a triangular element is determined by the total number of local basis functions corresponding to that triangular element in step S4. Decide.

[0088] S52. Construct the local mass matrix: When the unit is the first e When there are units, its local mass matrix is ​​denoted as The matrix elements of the local quality matrix Represented as: in, For the first e Each physical unit region and For the first e Local orthogonal basis functions on each unit.

[0089] Since the local basis functions constructed within each element in S4 are orthogonal, therefore when i ≠ j Sometimes, ;when i = j Sometimes, .

[0090] Therefore, the first e Local mass matrix of each unit The specific expression is: In the formula, Operators that construct diagonal matrices for each component of a given vector, given the elements of the diagonal. For the first e The diagonal term of the mass matrix corresponding to the first local orthogonal basis function on each unit. For the first e The first unit The diagonal terms of the mass matrix corresponding to each locally orthogonal basis function For the first e The total number of local orthogonal basis functions on each unit; The local mass matrix of the first diagonal elements Approximately written as: In the formula, For the first e The first unit on Locally orthogonal basis functions For the first g One physical integration point, , For the first e The total number of physical integration points in each unit. For the first g Each physical integral weight; Furthermore, the first e The local basis function values ​​at each integration point of each element are organized into a local basis function value matrix. In this matrix, rows correspond to local basis function numbers, and columns correspond to integration point numbers. The physical integration weights are written as a vector: In the formula, For the first e Local orthogonal basis function value matrix over each element For the first e The total number of local orthogonal basis functions on each unit. Multiplying matrices element by element For the first e The total number of physical integration points in each unit. Let be the physical integral weight vector; thus, the th e Local mass matrix of each unit The specific expression is: In the formula, Operators that construct diagonal matrices for each component of a given vector, given the elements of the diagonal. This is a vector of diagonal elements of the local mass matrix; For both rectangular and triangular elements, the local mass matrix is ​​calculated using the same method described above. The difference lies in the local basis functions of the rectangular element. Derived from Legendre's local orthogonal basis functions, local basis functions of triangular elements. It originates from KD local orthogonal basis functions.

[0091] Using the above method, it is not necessary to calculate the complete local mass matrix. Instead of directly calculating its diagonal elements and constructing a local mass matrix in diagonal form, it reduces computation and storage while maintaining consistency with the properties of local orthogonal basis functions.

[0092] S53. Construct the local correlation matrix: For any two units and Define the local correlation matrix between their unit pairs as Matrix elements of the local correlation matrix Represented as: In the formula, and The first e The unit and the first f Local orthogonal basis functions on each unit. In MATLAB, the preferred method for calculating the above integral is double numerical integration. Let the first... e The unit and the first f The physical integration points and weights of each unit are as follows: and .

[0093] Then the elements of the local correlation matrix Approximately written as: Preferably, the first is defined e The unit and the first f The integral point correlation weight matrix between units Its matrix elements are Then the first e The unit and the first f Local correlation matrix between units The specific expression is: In the formula, For the first f Local orthogonal basis function value matrix over each element It is the transpose symbol. For the first e The unit and the first f The integral point correlation weight matrix between units; when e = f When, the above formula represents the local correlation matrix within the unit; when e ≠ f When the autocorrelation function describes the correlation between any two points within the target region, the local correlation matrix exists not only within the same unit but also between different units.

[0094] S54. Assemble the global quality matrix: Based on the mapping relationship between local basis function numbers and global basis function numbers, the local quality matrix of each element is written into the global quality matrix as a diagonal block to obtain the global quality matrix; Based on the mapping relationship between the local basis function numbers and the global basis function numbers established in S51, the local quality matrix of each element is written into the global quality matrix as a diagonal block to obtain the global quality matrix B.

[0095] Since the local basis functions of each element in this invention are defined only within their respective elements and are not required to be continuous at the interfaces between different elements, the local basis function regions on different elements do not overlap. Therefore, the global mass matrix B can be written as a block diagonal matrix composed of the local mass matrices of each element: in, n e B represents the total number of cells in the hybrid mesh. (e) For the first eThe local mass matrix corresponding to each unit. Furthermore, since each local mass matrix B... (e) Under the influence of local orthogonal basis functions, it is a diagonal matrix, therefore the global quality matrix B has a block diagonal and highly sparse structure.

[0096] S55. Assemble the global correlation matrix: Based on the mapping relationship between local basis function numbers and global basis function numbers, the local correlation matrix between each pair of units is written into the global correlation matrix to obtain the global correlation matrix.

[0097] Global correlation matrix It can be represented in block form by unit as follows: In the formula, C (e,f) For the first e The unit and the first f The local correlation matrix between units.

[0098] In MATLAB, the local correlation matrices between each unit pair are written into the corresponding positions of the global correlation matrix according to the mapping relationship between local basis function numbers and global basis function numbers, thus completing the assembly of the global correlation matrix. Unlike the global quality matrix, the global correlation matrix... It is generally represented as a densely partitioned matrix, rather than a block diagonal matrix. This occurs when the autocorrelation function satisfies symmetry, i.e. From time to time Therefore, the global correlation matrix Theoretically, it is a symmetric matrix.

[0099] S6 includes the following steps: S61. Establish the generalized eigenvalue equation: Based on the global correlation matrix and global quality matrix obtained in S5, the corresponding algebraic generalized eigenvalue problem is expressed as follows: .

[0100] In the formula, This is the global correlation matrix. The global quality matrix, Let be the eigenvector matrix, and let its eigenvector be the eigenvector matrix. Listed as the number eigenvectors , It is an eigenvalue diagonal matrix. N Let be the dimension of the discrete descendant generalized eigenvalue problem. Equivalently, it can also be written in matrix form. The algebraic generalized eigenvalue problem originates from the integral eigenvalue problem established in S3, which is discretized by Galerkin in the locally orthogonal basis function space of the hybrid unit constructed in S4, and then obtained by constructing the local and global matrices in S5.

[0101] S62. Solving for eigenvalues ​​and eigenvectors: In MATLAB, the generalized eigenvalue equation is numerically solved to obtain eigenvalues ​​and corresponding eigenvectors. The obtained eigenvalues ​​are denoted as... The corresponding eigenvector is denoted as , N The dimension of the discrete descendant number generalized eigenvalue problem; Preferably, the above equations are solved in MATLAB using numerical methods applicable to symmetric generalized eigenvalue problems.

[0102] S63. Sort and filter the feature values: The obtained eigenvalues ​​are sorted in descending order, and the order of the corresponding eigenvectors is adjusted simultaneously. Furthermore, based on the numerical precision requirements, non-positive eigenvalues, minimal eigenvalues, or invalid eigenvalues ​​caused by numerical errors are eliminated, while valid eigenvalues ​​and corresponding eigenvectors that meet the discretization requirements of the random field are retained. The sorting and filtering steps are beneficial for determining the number of terms in the subsequent random field expansion and for extracting the main eigenvalues.

[0103] S64. Normalize the feature vectors: Since the eigenvectors obtained in S62 may not satisfy the orthogonality and normalization requirements under the mass matrix, the eigenvectors need to be normalized. The normalized eigenvectors are denoted as follows: And it satisfies the orthogonality and normalization condition in the sense of the mass matrix. In particular, when i = j Sometimes, .

[0104] After normalization, eigenvectors of different orders satisfy the orthogonality requirement in the sense of the mass matrix, which facilitates the subsequent recovery of approximate eigenfunctions and the construction of discrete expressions for random fields.

[0105] S65. Constructing approximate characteristic functions: Calculate the approximate characteristic function based on the normalized eigenvectors. : In the formula, For the first e Indicator functions for each unit, For the first e The area occupied by each unit For the first e The first unit on p The local orthogonal basis functions, where, when the... e When each unit is a rectangular unit Represent the Legendre local orthogonal basis functions on the rectangular element; when the... e When each unit is a triangular unit Denotes the KD local orthogonal basis functions on the triangular unit. The numbering of local basis functions in the normalized eigenvectors ( e , p The components corresponding to the discrete eigenvalues ​​can be obtained. The corresponding approximate characteristic function expression.

[0106] S7, in the MATLAB environment, is used to determine the number of truncation terms in the random field expansion based on the effective eigenvalues, normalized eigenvectors, and approximate eigenfunctions in the form of element partitions obtained in step S6. It then constructs a finite-term discrete expression for the random field and calculates the discrete values ​​of the random field at evaluation points within the target region, providing input for subsequent geotechnical parameter assignment, numerical analysis, reliability calculation, or random field visualization. S7 includes the following sub-steps: S71. Determine the number of truncated terms in the random field expansion: The number of truncated expansion terms used in the discrete random field analysis is determined based on the effective eigenvalue sequence. M ; S72. Generating random variables: Generate a vector of random variables based on the number of terms in the truncated expansion. , Let be the random variable corresponding to the first expansion term. For the first M Random variables corresponding to each expansion term; For a normal random field, it is preferable to generate a set of pairwise uncorrelated standard normal random variables. and satisfy .in, θ For sample space variables, For Kronecker symbols.

[0107] When the target random field is a normal random field, the standard normal random variable can be directly used in the finite-term discrete expression of the random field; when the target random field is a log-normal random field or other non-normal random field, it is preferable to first construct a normal random field, and then apply the corresponding distribution transformation according to the target probability distribution type to obtain the target random field.

[0108] S73. Constructing a finite-term discrete expression for a normal random field: Substituting the effective eigenvalues ​​and approximate eigenfunctions into the Karhunen–Loève expansion of the random field to be discretized in S31, we obtain the finite-term discrete expression of the random field based on element-partitioned local basis functions. : This formula reflects that the discrete result of a random field is determined by eigenvalues, normalized eigenvectors, local orthogonal basis functions of the unit, and random variables.

[0109] When the target random field is a log-normal random field or other non-normal random field, it is preferable to first construct a normal random field according to the above formula, and then apply the corresponding distribution transformation according to the target probability distribution type to obtain the final target random field.

[0110] S74. Calculate the discrete values ​​of the random field at discrete points: In the target area The set of discrete points in the random field is selected. These discrete points can be user-specified points, regularly arranged points, boundary supplement points, or calculation points required for subsequent numerical analysis. Let the first... a The coordinates of the discrete points are .in, N p This represents the total number of evaluation points.

[0111] If discrete points x a Located in the e a Within a unit, because If all other unit indicator functions are 0, then... The formula can be simplified to: The above formula represents the actual calculation form of the discrete values ​​of the random field at discrete points. This formula only involves the local basis function values ​​of the element to which the discrete point belongs, thus it is consistent with the piecewise expression of local orthogonal basis functions and facilitates rapid calculation in MATLAB based on the element to which the discrete point belongs.

[0112] In the target area Select a set of discrete points for a random field, and for all discrete points... Point-by-point calculation yields the discrete value vector of the random field. : In the formula, For discrete points The discrete expression of a random field with finite terms. For discrete points The discrete expression of a random field with finite terms. For discrete points The discrete expression of a random field with finite terms. This represents the total number of discrete points.

[0113] In a specific embodiment, a two-dimensional quadratic exponential autocorrelation function is used for both the slope region and the inclined stratum region, with a horizontal correlation length of 20 m and a vertical correlation length of 2 m. The mean of the cohesion random field is 20 kPa and the standard deviation is 4 kPa, while the mean of the internal friction angle random field is 30° and the standard deviation is 3°. Different random field implementations can be obtained by generating mutually independent standard normal random variable vectors. To reduce the length of the figures, the first and second implementations of each type of random field are shown. The implementation results of the cohesion random field and the internal friction angle random field in the slope region are as follows: Figure 6 and Figure 7 As shown, the implementation results of the random field of cohesion and the random field of internal friction angle in the inclined strata region are as follows: Figure 8 and Figure 9 As shown.

[0114] In the description of this invention, the above are merely preferred embodiments and are not intended to limit the scope of protection of this invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this invention should be included within the scope of protection of this invention.

Claims

1. A method for discretizing locally orthogonal basis random fields of hybrid units in irregular regions, characterized in that: Includes the following steps: S1. Establish a geometric model of the discrete region of the random field and determine the discrete region of the random field; S2. Perform hybrid meshing on the discrete region of the random field to obtain a set of rectangular cells and a set of triangular cells, thus generating a hybrid mesh; S3. Based on the statistical parameters and autocorrelation function model of the random field to be discretized, establish the Karhunen-Loève expansion expression, integral eigenvalue problem, and Galerkin discrete solution framework based on the local orthogonal basis of the hybrid grid corresponding to the discrete region of the random field. S4. Construct local orthogonal basis functions on rectangular and triangular elements, and calculate the local basis function values ​​at the integration points of each element; S5. The Galerkin method is used to discretize the integral eigenvalue problem. Local mass matrix and local correlation matrix are constructed using the local orthogonal basis functions of rectangular and triangular elements, respectively, and then assembled to form global mass matrix and global correlation matrix. S6. Establish a generalized eigenvalue equation based on the global quality matrix and the global correlation matrix, and solve it to obtain the effective eigenvalues ​​and corresponding eigenvectors, thus obtaining the approximate eigenfunction. S7. Based on effective eigenvalues ​​and approximate eigenfunctions, construct a finite-term discrete expression for the normal random field, and calculate the discrete value of the random field within the target region to complete the discretization of the random field within the irregular region.

2. The method for discretizing hybrid unit local orthogonal basis random fields considering irregular regions according to claim 1, characterized in that, S2 specifically refers to: The discrete region of the random field is divided into rectangular units using quadtrees to generate a set of rectangular units. The remaining regions not covered by the rectangular units and the regions near irregular boundaries are divided into triangular units to generate a set of triangular units. The set of rectangular units and the set of triangular units together cover the discrete region of the random field.

3. The method for discretizing hybrid unit local orthogonal basis random fields considering irregular regions according to claim 1, characterized in that, S3 includes the following steps: S31. Determine the random field to be discrete: The random field to be discretized is determined based on the discrete region of the random field, and it is expressed in the Karhunen–Loève expansion form: In the formula, For the discrete random field, The first spatial location point within the target area. For the target area, For sample space variables, For the sample space, It is a mean function. It is a function of standard deviation. For the first i 1 eigenvalue, For the first i One characteristic function, For the first i One random variable; S32. Establish a random field autocorrelation function model: Based on the spatial correlation characteristics of the random field to be discretized, the autocorrelation function of the random field is established. ,in, Indicates the second spatial location point within the target area; S33. Problem of establishing integral eigenvalues: Based on the autocorrelation function of the random field Establish the second type of Fredholm integral eigenvalue problem, the specific expression of which is: S34. Establish a Galerkin solution approach based on a hybrid mesh local orthogonal basis: The Galerkin discretization method based on local orthogonal basis of hybrid mesh is used to numerically solve the integral eigenvalue problem. Specifically, based on the hybrid mesh generated by S2, the eigenfunctions in the integral eigenvalue problem are approximated as a linear combination of global piecewise basis functions. For rectangular elements, the Legendre local orthogonal basis function is subsequently used; for triangular elements, the KD local orthogonal basis function is subsequently used.

4. The method for discretizing hybrid unit local orthogonal basis random fields considering irregular regions according to claim 3, characterized in that, S4 includes the following steps: S41. Determine the element type and reference element: The elements in the hybrid mesh are identified, and corresponding reference elements are created according to the element type. When the element is a rectangular element, a reference rectangular element is used. ; When the element is a triangular element, then a reference triangular element is used. ,in, , For local coordinates within the reference triangle element; S42. Construct Legendre local orthogonal basis functions controlled by the total degree on the rectangular element, generating local orthogonal basis functions on the rectangular element: No. Two-dimensional local orthogonal basis functions on rectangular units The specific expression is: In the formula, and For the first One-dimensional local Legendre basis functions on rectangular units and These represent the horizontal and vertical coordinate components in the physical coordinate system of the rectangular unit, respectively. p and indicators q satisfy , For the first The order of the local basis function truncation of each rectangular cell is used to limit the two-dimensional index pair. p , q The highest total number of times; S43. Construct KD local orthogonal basis functions controlled by the total degree in the triangular element to generate local orthogonal basis functions on the triangular element: No. Two-dimensional KD local orthogonal basis functions on triangular units The specific expression is: In the formula, For Jacobi polynomials in parameters The form of time corresponds to the Legendre polynomial. For parameters Jacobi polynomials, where, and For the parameters of the Jacobi polynomial, and Points inside the physical triangle unit ( x , y The local coordinates of the reference triangle element obtained through affine inverse mapping, index p and indicators q satisfy , For the first The order of the local basis function truncation of each triangular unit is used to constrain the two-dimensional index pair. p , q The highest total number of times; S44. Calculate the local basis function values ​​at each unit integration point: Based on the local orthogonal basis functions of the constructed rectangular and triangular elements, numerical integration points are selected on the corresponding reference elements, and the local basis function values ​​at each integration point are calculated.

5. The method for discretizing hybrid unit local orthogonal basis random fields considering irregular regions according to claim 4, characterized in that, S5 includes the following steps: S51. Establish the mapping relationship between the local orthogonal basis function and the global orthogonal basis function numbers: For each unit in the hybrid network, establish local orthogonal basis functions and assign them a unified number. According to the unit numbering order, establish a mapping relationship between the local basis function numbers and the global basis function numbers. The mapping relationship is used to write the local quality matrix and the local correlation matrix into the corresponding positions of the global matrix. S52. Construct the local mass matrix: No. e Local mass matrix of each unit The specific expression is: In the formula, Operators that construct diagonal matrices for each component of a given vector, given the elements of the diagonal. For the first e The diagonal element vector of the local mass matrix of each unit cell; In the formula, For the first e Local orthogonal basis function value matrix over each element For the first e The total number of local orthogonal basis functions on each unit. For the first e The total number of physical integration points in each unit. This is the physical integral weight vector; S53. Construct the local correlation matrix: No. e The unit and the first f Local correlation matrix between units The specific expression is: In the formula, For the first f Local orthogonal basis function value matrix over each element It is the transpose symbol. For the first e The unit and the first f The integral point correlation weight matrix between units; S54. Assemble the global quality matrix: Based on the mapping relationship between local basis function numbers and global basis function numbers, the local quality matrix of each element is written into the global quality matrix as a diagonal block to obtain the global quality matrix; S55. Assemble the global correlation matrix: Based on the mapping relationship between local basis function numbers and global basis function numbers, the local correlation matrix between each pair of units is written into the global correlation matrix to obtain the global correlation matrix.

6. The method for discretizing hybrid unit local orthogonal basis random fields considering irregular regions according to claim 5, characterized in that, S6 includes the following steps: S61. Establish the generalized eigenvalue equation: Based on the global mass matrix and the global correlation matrix, establish the generalized eigenvalue equation: In the formula, This is the global correlation matrix. The global quality matrix, Let be the eigenvector matrix, and let its eigenvector be the eigenvector matrix. Listed as number eigenvectors , It is an eigenvalue diagonal matrix; S62. Solving for eigenvalues ​​and eigenvectors: Numerical solutions are used to solve the generalized eigenvalue equations to obtain eigenvalues ​​and corresponding eigenvectors, where the obtained eigenvalues ​​are denoted as . The corresponding eigenvector is denoted as , N The dimension of the discrete descendant number generalized eigenvalue problem; S63. Sort and filter the feature values: The obtained eigenvalues ​​are sorted in descending order, and the order of the corresponding eigenvectors is adjusted simultaneously. Furthermore, based on the numerical precision requirements, non-positive eigenvalues, minimal eigenvalues, or invalid eigenvalues ​​caused by numerical errors are eliminated, while the valid eigenvalues ​​and corresponding eigenvectors that meet the discretization requirements of the random field are retained. S64. Normalize the feature vectors: The eigenvectors are normalized, and the normalized eigenvectors are denoted as follows: ; S65. Constructing approximate characteristic functions: Calculate the approximate characteristic function based on the normalized eigenvectors. : In the formula, For the first e Indicator functions for each unit, For the first e The area occupied by each unit For the first e The first unit on p Locally orthogonal basis functions The numbering of local basis functions in the normalized eigenvectors ( e , p The corresponding component.

7. The method for discretizing hybrid unit local orthogonal basis random fields considering irregular regions according to claim 6, characterized in that, S7 includes the following steps: S71. Determine the number of truncated terms in the random field expansion: The number of truncated expansion terms used in the discrete random field analysis is determined based on the effective eigenvalue sequence. M ; S72. Generating random variables: Generate a vector of random variables based on the number of terms in the truncated expansion. , Let be the random variable corresponding to the first expansion term. For the first M Random variables corresponding to each expansion term; S73. Constructing a finite-term discrete expression for a normal random field: Substituting the effective eigenvalues ​​and approximate eigenfunctions into the Karhunen–Loève expansion of the random field to be discretized, we obtain the finite-term discrete expression of the random field based on element-partitioned local basis functions. : S74. Calculate the discrete values ​​of the random field at discrete points: In the target area Select a set of discrete points for a random field, and for all discrete points... Point-by-point calculation yields the discrete value vector of the random field. : In the formula, For discrete points The discrete expression of a random field with finite terms. For discrete points The discrete expression of a random field with finite terms. For discrete points The discrete expression of a random field with finite terms. This represents the total number of discrete points.

Images