A method for predicting subgrade dynamic response based on physical constraint Kolmogorov-Arnold long short-term memory network

Abstract

The application provides a roadbed dynamic response prediction method based on a physical constraint Kolmogorov-Arnold long short-term memory network. The method obtains input ground motion at the bottom of the roadbed and response time history of a key response point, carries out pretreatment, sequence organization and displacement, velocity and acceleration state label construction; uses a long short-term memory network to extract time sequence characteristics, and outputs multi-state responses through a Kolmogorov-Arnold learnable activation decoding module; and trains the model according to data loss and displacement-velocity, velocity-acceleration derivative consistency loss. When applied, input the to-be-evaluated ground motion, and the acceleration time history and peak value of the top surface of the roadbed or the internal key point and other engineering indexes are obtained through response recovery, and the method can be used for rapid evaluation of highway roadbeds, airport runway foundations and the like.

Description

Technical Field

[0001] This invention belongs to the interdisciplinary fields of earthquake engineering, road engineering and artificial intelligence, and specifically relates to a method for predicting the dynamic response of roadbed based on a physically constrained Kolmogorov-Arnold long short-term memory network. Background Technology

[0002] The dynamic response of a roadbed system under seismic loading is influenced by a combination of factors, including the input ground motion at the base, the thickness of the overlying structural layers, the materials of the base and subbase layers, the stiffness and damping of the roadbed soil, boundary conditions, and soil nonlinearity. In the rapid post-earthquake assessment of highway roadbeds, airport runway foundations, urban roadbeds, and similar transportation infrastructure, it is necessary to focus on the acceleration, displacement, and related engineering parameters of the center and edge regions of the roadbed top surface and key points within the roadbed.

[0003] Current dynamic response analysis of subgrade systems typically relies on finite element method, finite difference method, or other numerical time history analysis methods. These methods can describe soil nonlinearity, structural layer interactions, and boundary conditions, but the modeling and solving of complete numerical models are costly in situations involving multiple ground motion inputs, multiple soil parameters, subgrade structure scheme comparison, and rapid post-earthquake assessment.

[0004] Deep sequence models can establish a mapping relationship between bottom-input ground motion and the response at key points of the roadbed, but purely data-driven models often lack explicit constraints on the basic kinematic relationships between displacement, velocity, and acceleration. Under conditions of strong earthquake input, weak roadbed, or insufficient sample coverage, the predicted curve may exhibit problems such as seemingly reasonable phase but mismatched derivative relationships, peak drift, or non-physical oscillations in the attenuation segment.

[0005] Physical information neural networks can improve physical consistency by introducing equation residuals. However, in roadbed engineering, material parameters, damping models, and boundary conditions for different sites and structural layers are not necessarily available in a uniform form. If the complete dynamic equilibrium equation is directly used as the loss constraint, the model may become dependent on a large number of engineering detail parameters, affecting its cross-scenario use.

[0006] Furthermore, the seismic time histories of roadbed systems are typically quite long. Directly inputting the complete long time histories into recurrent neural networks increases the time step expansion length and training difficulty; simple truncation or downsampling may disrupt the correspondence between the dominant segment, attenuation segment, and low-frequency components. Therefore, a method for predicting the dynamic response of roadbeds that balances long time histories processing efficiency, nonlinear mapping capability, and kinematic consistency is needed. Summary of the Invention

[0007] The technical problem to be solved by this invention is: to address the problems of insufficient physical consistency, limited nonlinear mapping capability, large training burden of long-term program sequences, and strong dependence on parameters of complete dynamic equations in existing fast prediction methods for subgrade system dynamic response, a subgrade dynamic response prediction method based on physically constrained Kolmogorov-Arnold long short-term memory network is proposed.

[0008] To solve the above-mentioned technical problems, the present invention adopts the following technical solution:

[0009] A method for predicting the dynamic response of roadbed based on a physically constrained Kolmogorov-Arnold long short-term memory network includes the following steps:

[0010] S1: Obtain the bottom input ground motion time history of the subgrade system and the target response time history of the key response points of the subgrade to form a training sample set;

[0011] S2: Preprocess the input ground motion time history and target response time history in the training sample set, and obtain or construct three types of response state labels for the key response points: displacement, velocity and acceleration.

[0012] S3: Perform sequence organization processing on the preprocessed input ground motion time history to obtain the input sequence for sequence modeling, and set the response recovery processing corresponding to the sequence organization processing;

[0013] S4: Construct a neural network model, which includes a long short-term memory network sequential encoder, a Kolmogorov-Arnold learnable activation decoding module, and a multi-state output module;

[0014] S5: Input the input sequence into the long short-term memory network temporal encoder to extract temporal hidden features, and generate displacement, velocity and acceleration response prediction sequences of the key response points through the Kolmogorov-Arnold learnable activation decoding module and multi-state output module;

[0015] S6: Construct a data loss based on the displacement, velocity, and acceleration response prediction sequence and the response state label, and construct a physical constraint loss based on the kinematic derivative relationship between displacement, velocity, and acceleration;

[0016] S7: Construct a joint loss function based on the data loss and physical constraint loss, and train the neural network model to obtain a trained roadbed dynamic response prediction model;

[0017] S8: Input the ground motion time history of the bottom of the subgrade system to be predicted into the trained subgrade dynamic response prediction model, and output the dynamic response prediction results of the key response points after the response recovery processing.

[0018] Preferably, the key response points include at least one of the following: the center point of the top surface of the subgrade, the edge point of the subgrade, the bottom point of the base course, the top surface point of the subgrade, or the measuring point inside the subgrade. The dynamic response prediction result includes at least the acceleration response time history of the key response points.

[0019] Preferably, in step S1, the training sample set comes from finite element numerical simulation, engineering monitoring records, or a combination of the two; the finite element numerical simulation includes the overlying structural layer, base course or subbase course, homogeneous subgrade soil, bottom seismic input boundary, and interlayer contact or continuous deformation relationship; the response state label is directly extracted from the numerical calculation results or engineering monitoring results, or constructed from the target response time history after integration, difference, or filtering correction.

[0020] Preferably, in step S3, the sequence organization processing includes window folding processing, which involves dividing the original input ground motion time history into multiple sub-sequences according to a preset window length and step size. The response recovery processing involves splicing, overlapping weighted averaging, or overlapping summation of the prediction results of the multiple sub-sequences to obtain the response prediction result of the complete time length.

[0021] Preferably, in step S3, the sequence organization processing includes sampling point rearrangement processing, which involves rearranging adjacent sampling points in the original input ground motion time history to the feature dimension, so as to shorten the time dimension and expand the single-step input feature dimension, and the response recovery processing includes inverse rearrangement processing corresponding to the sampling point rearrangement processing.

[0022] Preferably, in step S4, the Kolmogorov-Arnold learnable activation decoding module is set after the long short-term memory network temporal encoder, and is used to map the temporal latent features into multiple learnable unary functions and their linear combinations. The learnable unary functions include one or more of spline functions, piecewise polynomial functions, or parameterized basis functions.

[0023] Preferably, in step S6, the physical constraint loss takes kinematic derivative consistency as the constraint object, and includes at least a first derivative consistency constraint and a second derivative consistency constraint; the first derivative consistency constraint is used to constrain the consistency between the time derivative of the predicted displacement and the predicted velocity, the second derivative consistency constraint is used to constrain the consistency between the time derivative of the predicted velocity and the predicted acceleration, and the construction of the physical constraint loss does not require a complete mass matrix, damping matrix or stiffness matrix as necessary input.

[0024] Preferably, the time derivative is obtained by one of the following methods: central difference method, forward difference method, backward difference method, smooth difference method, or automatic differentiation method, and is supplemented at the end points of the sequence by one-sided difference or end-point extrapolation.

[0025] Preferably, in step S7, the joint loss function weights the data loss and physical constraint loss according to the response state magnitude, the training phase loss ratio, or a preset weight, and employs one or more training strategies among learning rate decay, early stopping, gradient pruning, or normalization.

[0026] Preferably, at least one of the following conditions is also met: the neural network model also receives one or more conditional features among the following: subgrade soil code, density, Poisson's ratio, initial shear modulus, damping parameter, overlying structural layer thickness, base course thickness, subgrade depth, or overlying structural layer material parameters; the dynamic response prediction results also include one or more of the following indices: acceleration response time history, peak acceleration, response spectrum, peak displacement, or duration index at multiple key response points.

[0027] In this invention, the complete input seismic motion sequence is denoted as X = {x(1), x(2), ..., x(T)}, with sampling time t(k) = (k-1)·Δt, k = 1, 2, ..., T; the m-th window when the window is collapsed is represented as X(m) = {x(s(m)), ..., x(s(m) + W - 1)}, where s(m) = 1 + (m - 1)·P, m = 1, 2, ..., M, M = floor((T - W) / P) + 1. In the formula, Δt is the sampling time interval, T is the length of the complete sequence, W is the window length, P is the window step size, m is the window number, s(m) is the starting sampling point number of the m-th window, M is the total number of windows, and floor(·) is the floor function. The above symbols are used to describe the method flow and do not require the same values ​​for all projects.

[0028] One possible construction relationship for the joint loss function is as follows:

[0029] L data = α u ·MSE(u p , u) + α v ·MSE(v p , v) + α a ·MSE(a p , a);

[0030] L phy = β1·MSE(D(u p ), v p) + β2·MSE(D(v p ), a p );

[0031] L = L data + λ·L phy .

[0032] In the formula, L data For data loss, L phy The physical constraint loss is represented by L, the joint loss function is represented by MSE, and the mean squared error function is represented by u. p v p and a p These represent the predicted displacement, predicted velocity, and predicted acceleration, respectively; u, v, and a are the corresponding response state labels; D is the time derivative operator; α u α v α a β1, β2, and λ are weighting coefficients. This set of relationships is used to limit the combinational logic between loss terms, but does not limit the specific weight values.

[0033] When the time derivative operator D adopts the central difference operator (CDM), its form can be expressed as:

[0034] D[q(k)] = [q(k + 1) - q(k - 1)] / (2·Δt), k = 2, 3, ..., T-1.

[0035] In the formula, q(k) is the response state quantity at the kth sampling time, q(k + 1) and q(k - 1) are the response state quantities at adjacent sampling times, and Δt is the sampling time interval; the sequence endpoints are supplemented by one-sided difference or endpoint extrapolation.

[0036] Compared with the prior art, the present invention has at least the following beneficial effects:

[0037] 1. This invention transforms the kinematic derivative relationship between displacement, velocity and acceleration at key response points of the subgrade into a trainable physical constraint loss. Without explicitly giving the complete mass matrix, damping matrix, stiffness matrix and soil constitutive parameters, it can still apply physical consistency constraints to the predicted dynamic response results of the subgrade system.

[0038] 2. This invention introduces a Kolmogorov-Arnold learnable activation decoding module after the Long Short-Term Memory Network (LSTM) timing encoder. By using multiple learnable univariate functions and their linear combinations, the nonlinear mapping expression capability is enhanced, which is beneficial for characterizing peak abrupt changes, frequency band shifts, soft foundation amplification, and local nonlinear waveform characteristics in the roadbed system.

[0039] 3. This invention adopts a three-state joint prediction method in the training phase and an engineering target output method in the application phase, so that the displacement and velocity states mainly serve the physical constraints and feature learning, while the acceleration time history and its engineering indicators can be used as the main output in engineering applications, thus taking into account both the sufficiency of training constraints and the simplicity of engineering use.

[0040] 4. This invention converts long-term earthquake motions into subsequences or rearranged sequences that are easy to learn through sequence organization and response recovery processing, and restores the complete response time history at the output end, which helps to reduce the training burden of long sequences while maintaining the time history correspondence between different time periods.

[0041] 5. The sample sources of this invention can be finite element numerical simulation, engineering monitoring records, or a combination of both. The physical constraints mainly rely on the universal kinematic relationship between response states and can further receive the characteristics of subgrade soil, overlying structural layers, and boundary conditions. Therefore, it can serve as an scalable modeling framework for multi-point response prediction of highway subgrade, airport runway foundation, or similar subgrade projects. Attached Figure Description

[0042] Figure 1 This is a flowchart of the overall process of the present invention, showing the processing relationships between sample construction, preprocessing, sequence organization, long short-term memory network encoding, Kolmogorov-Arnold decoding, physical constraint training, and response output.

[0043] Figure 2 This is a schematic diagram of roadbed data construction in an embodiment of the present invention, showing the relationship between the bottom input ground motion, the overlying structural layer, the homogeneous roadbed soil, the response points on the top surface of the roadbed, and the boundaries;

[0044] Figure 3 This is a schematic diagram of the physically constrained Kolmogorov-Arnold long short-term memory network structure in an embodiment of the present invention, showing the connection relationship between the sequence organization input, the long short-term memory network temporal encoder, the Kolmogorov-Arnold learnable activation decoding module, and the three-state output.

[0045] Figure 4 This is a schematic diagram of the Kolmogorov-Arnold learnable activation decoding module in an embodiment of the present invention, showing the mapping relationship between temporal latent features, learnable unary functions, linear combinations and response state outputs;

[0046] Figure 5 This is a schematic diagram of the construction of the joint loss function in an embodiment of the present invention, showing the data fitting loss, derivative consistency loss and their weighted training relationship;

[0047] Figure 6This is a schematic diagram of window folding and response recovery in an embodiment of the present invention, showing the process of subsequence prediction after long-term input is divided into windows, and the complete response time history is recovered through overlapping weighting.

[0048] Figure 7 This is a representative bottom-input ground motion time history diagram in an embodiment of the present invention, showing the change of input acceleration over time and the peak acceleration labeling method;

[0049] Figure 8 This is a comparison of the response time histories of representative roadbed top surfaces in embodiments of the present invention; wherein, Figure 8 (a) represents the low-strength input of subgrade soil A. Figure 8 (b) represents the medium strength input of subgrade soil A. Figure 8 (c) represents the high-strength input of subgrade soil A. Figure 8 (d) represents the low-strength input of subgrade soil B. Figure 8 (e) represents the medium strength input for subgrade soil B. Figure 8 (f) represents the high-strength input of subgrade soil B;

[0050] Figure 9 This is a scatter plot showing the correspondence between peak acceleration prediction points in an embodiment of the present invention; wherein, Figure 9 (a) represents the low-strength input of subgrade soil A. Figure 9 (b) represents the medium strength input of subgrade soil A. Figure 9 (c) represents the high-strength input of subgrade soil A. Figure 9 (d) represents the low-strength input of subgrade soil B. Figure 9 (e) represents the medium strength input for subgrade soil B. Figure 9 (f) represents the high-strength input of subgrade soil B;

[0051] Figure 10 This is a statistical chart of the overall time history error and correlation coefficient in an embodiment of the present invention; wherein, Figure 10 (a) shows the statistical results under condition A for the subgrade soil. Figure 10 (b) shows the statistical results under subgrade soil condition B;

[0052] Figure 11 This is a statistical chart of peak acceleration error and confidence interval coverage in an embodiment of the present invention; wherein, Figure 11 (a) shows the statistical results under condition A for the subgrade soil. Figure 11 (b) shows the statistical results under the subgrade soil condition B. Detailed Implementation

[0053] The present invention will be further described below with reference to the accompanying drawings and specific embodiments. These embodiments are only for illustrating the present invention and are not intended to limit the scope of protection of the present invention. For those skilled in the art, equivalent substitutions, simplifications, or extensions made without departing from the concept of the present invention should all fall within the scope of protection of the present invention.

[0054] Example 1: Overall Method Flow

[0055] This embodiment illustrates a method for predicting the dynamic response of roadbeds based on physically constrained Kolmogorov-Arnold long short-term memory networks. The overall process is as follows: Figure 1 As shown, it mainly includes seven stages: sample construction, data preprocessing, sequence organization, temporal coding, KAN decoding, physical constraint training, and response prediction output.

[0056] During the sample construction phase, the time histories of the bottom input ground motion and the target response time histories of key response points of the subgrade system are obtained. The bottom input ground motion is the horizontal acceleration time history; key response points are selected according to the engineering focus, including the center point P0 of the subgrade top surface, subgrade edge points, bottom points of the base course, points on the top surface of the subgrade, or measuring points inside the subgrade. The target response time histories include displacement, velocity, and acceleration; when the sample only contains acceleration response, other state labels can be supplemented through numerical integration, numerical differentiation, or numerical simulation output.

[0057] The preprocessing stage performs unified processing on the input ground motion and target response. This processing includes one or more of the following: baseline correction, filtering, resampling, normalization, zero-padding extension, or outlier removal. For long-duration ground motions, sequence organization processing uses a fixed window length and step size for window folding, converting each complete time history into several subsequences. Sequence organization processing can also employ a sampling point rearrangement method, rearranging adjacent sampling points to the feature dimension to shorten the time dimension and expand the single-step input feature dimension.

[0058] In the model building phase, a physically constrained Kolmogorov-Arnold long short-term memory network is established. This network takes the sequence-organized input sequence as input, extracts temporal latent features using the long short-term memory network, enhances the nonlinear mapping expression capability using a Kolmogorov-Arnold learnable activation decoding module, and generates three types of predicted states—displacement, velocity, and acceleration—during the training phase through a multi-state output module.

[0059] In the loss function construction phase, a data fitting loss is constructed based on three types of response states, and a kinematic derivative consistency constraint is constructed based on the relationship that the first derivative of displacement with respect to time equals velocity, and the first derivative of velocity with respect to time equals acceleration. The kinematic derivative consistency constraint does not require explicit input of mass, damping, stiffness, or restoring force terms from the complete dynamic equilibrium equations, and can serve as a physical prior that is relatively decoupled from the specific roadbed structural parameters.

[0060] After model training is complete, the ground motion to be predicted is input into the model to obtain the prediction results for each window, subsequence, or rearranged sequence. The complete time history is then recovered through splicing, overlapping weighted averaging, overlapping summation, or inverse rearrangement. In the application stage, only the acceleration time history at the center point P0 of the roadbed top surface can be output, or displacement, velocity, peak acceleration, response spectrum, or other engineering response indicators can be further output.

[0061] Example 2: Construction of Roadbed Samples

[0062] like Figure 2 As shown, the subgrade system sample consists of an overlying structural layer, a base course or subbase course, homogeneous subgrade soil, a bottom seismic input boundary, and response extraction points. The seismic motion is input from the bottom of the system in the form of horizontal acceleration time history. The overlying structural layer and the subgrade soil can have a continuous, bonded, or equivalent contact relationship to form the subgrade dynamic response sample.

[0063] Training samples can be generated from two-dimensional plane strain finite element models, axisymmetric models, equivalent layered models, or engineering monitoring data. The thickness of the overlying structure layer, the base course thickness, the subgrade depth, the width of the computational domain, and material parameters are set according to the engineering object or standard calculation example. The target response point is preferably set at the center point P0 of the top surface of the subgrade, used to characterize the key response after the earthquake motion is transmitted through the subgrade and amplified or attenuated by the overlying structure layer; alternatively, measuring points inside the subgrade can be set simultaneously for multi-point prediction.

[0064] In this embodiment, the subgrade soil is defined as two types of homogeneous soil suitable for subgrade engineering, and is represented by codes to avoid limiting the invention to specific soil names or specific examples. Subgrade soil A represents low-stiffness, high-compressibility subgrade soil, suitable for soft or moderately soft subgrades; subgrade soil B represents medium-stiffness subgrade soil, suitable for subgrades with high compaction or good bearing capacity. Subsequent descriptions will use subgrade soil A and subgrade soil B, without using specific soil names as limitations.

[0065] Ground motion records are standardized in sampling interval and duration before being input into the model, and polynomial baseline correction and bandpass filtering are performed to reduce the impact of low-frequency drift and high-frequency noise on sequence learning and derivative consistency constraints. Preprocessed bottom input accelerations and key response point state variables form sample pairs.

[0066] As an optional controlled numerical example, the subgrade system can be idealized as a shear wave propagation system consisting of an overlying structural layer, homogeneous subgrade soil, and an equivalent input base, with ground motions propagating upwards in the form of horizontal shear waves. This example can neglect three-dimensional edge effects to highlight the sequential mapping relationship between the bottom-input ground motion and the response at the center point P0 of the subgrade top surface.

[0067] The roadbed system discretized by finite element method is subjected to horizontal acceleration input a at the base. g Under the action of (t), the dynamic equilibrium relationship can be used: M·a(t) + C·v(t) + K·u(t) = -M·r·a g The expression (t) is used to describe this. Here, u(t), v(t), and a(t) represent the nodal displacement, nodal velocity, and nodal acceleration vectors, respectively; M, C, and K represent the mass matrix, damping matrix, and stiffness matrix, respectively; and r represents the input influence vector. This dynamic equilibrium relationship is used to explain the physical origin of the sample response. The training constraints of this invention do not require the explicit input of the complete M, C, and K matrices in the loss function.

[0068] Boundary conditions can be set as follows: the top surface of the roadbed is a free boundary, a horizontal seismic motion input is applied to the base end, the lateral boundary restricts normal displacement and maintains tangential freedom, or an equivalent absorbing boundary is used to reduce wave reflection. The calculation process can be divided into a static stage and a dynamic stage. In the static stage, self-weight or initial ground stress is applied, and in the dynamic stage, the time history of horizontal seismic acceleration is input and the response at P0 is extracted.

[0069] Table 1 provides an optional parameter combination for two types of homogeneous subgrade soil. Density ρ, Poisson's ratio ν, and initial shear modulus G0 are used to define the basic dynamic characteristics of the subgrade soil; damping ratio and nonlinear adjustment coefficient are used to characterize energy dissipation and the equivalent nonlinearity under strong earthquake input. The above parameters can be adjusted based on engineering surveys, laboratory tests, field wave velocity tests, or numerical inversion results, and are not used as numerical limitations on the protection range.

[0070] Table 1. Soil parameters for two types of homogeneous subgrade in the optional controlled examples.

[0071]

[0072] In this controlled simulation, the seismic ground motion records can be selected from an automatic strong ground motion record library. Horizontal acceleration time histories with peak ground acceleration (PGA) covering low, medium, and high intensity ranges and a duration not exceeding 50 s are selected. To standardize the input length and sampling interval, the seismic records can be zero-padding and extended to 50 s, and resampled to Δt=0.02 s, resulting in each record having a uniform set of 2501 sampling points. Figure 7 The table shows a representative bottom input ground motion time history. The peak acceleration label is only used to illustrate the expression method of the input sample and is not intended to limit the numerical range of the input ground motion intensity of this invention.

[0073] The aforementioned seismic motions were input into the subgrade models corresponding to subgrade soil A and subgrade soil B, respectively, resulting in two sets of input-output samples; the input was the base seismic motion a. g (t), the output is the acceleration response at P0. During the training phase, the displacement and velocity responses at point P0 are also extracted simultaneously as physical constraint labels. The processed data can be divided into training, validation, and test sets, for example, 70%, 10%, and 20%.

[0074] Example 3: Network Structure and Kolmogorov-Arnold Decoding Module

[0075] like Figure 3 As shown, the neural network model in this embodiment includes a Long Short-Term Memory (LSTM) network temporal encoder, a Kolmogorov-Arnold learnable activation decoding module, and a multi-state output module. The LSTM network temporal encoder may include one or more LSTM layers for extracting long-term dependencies and local non-stationary features from the seismic motion input sequence.

[0076] As a specific network configuration, the encoding end can adopt a three-layer long short-term memory network structure, with the three hidden dimensions set to 128, 256, and 256 respectively, while keeping the time step length of the input subsequence constant. When each complete 50-second time history is divided into five 10-second subsequences with a sampling interval of 0.02 seconds, the time step length of a single subsequence can be 500. In actual implementation, it can also be adjusted to 500 or 501 depending on whether endpoint sampling points are included.

[0077] The output of the three-layer Long Short-Term Memory (LSTM) network is a hidden state sequence H. H contains both the local waveform characteristics of the bottom input ground motion and the temporal dependency between the main vibration segment and the decay segment. After H enters the decoder, it undergoes a nonlinear transformation through a KAN decoding module containing a learnable univariate function φi(·), and is then linearly combined to form a three-state output of displacement, velocity, and acceleration. Here, H is the hidden state sequence output by the LSM network, and φi(·) is the hidden state sequence. i (·) is the i-th learnable unary function in the Kolmogorov-Arnold decoding module.

[0078] like Figure 4 As shown, the KAN decoding module can replace the conventional linear layer followed by fixed ReLU or tanh activation. Each φ i (·) can be represented by spline basis functions or piecewise polynomials, and its control points, coefficients or node parameters are updated together during model training; linear combination weights are used to combine the outputs of multiple learnable univariate functions into the target response state.

[0079] The number of output channels of the multi-state output module can be adjusted according to... × Settings, where The number of critical response points. The number of states. When only one measurement point P0 is predicted and the training phase outputs three states: displacement, velocity, and acceleration. =1, =3; When multiple measurement points are output simultaneously, each measurement point can share the same timing encoder and form a parallel channel at the output end.

[0080] Example 4: Construction of Physical Constraint Loss Function

[0081] like Figure 5 As shown, the loss function in this embodiment consists of the data fitting term ( Figure 5 (represented by L1) and physical consistency terms ( Figure 5 It consists of L2 (represented in Chinese). The data fitting term is used to constrain the consistency between the predicted displacement, predicted velocity, and predicted acceleration and the corresponding labels; the physical consistency term is used to constrain the kinematic derivative relationship between the three types of predicted states.

[0082] Let the predicted displacement, predicted velocity, and predicted acceleration be u, respectively. p v p and a p Given labels u, v, and a, the data fitting term is composed of a weighted average of three types of mean square errors. For differences in the magnitude of different response states, the losses are balanced using amplitude normalization, standard deviation normalization, or an adaptive weighting strategy.

[0083] The data fitting term can be further decomposed into displacement error L disp Speed ​​error L vel and acceleration error L acc The details are as follows:

[0084] L disp = MSE(u p , u),

[0085] L vel = MSE(v p , v),

[0086] L acc = MSE(a p , a).

[0087] For the nth sample and the kth time point, the mean square error can be calculated as the average of the squared differences between the predicted sequence and the label sequence. When there are multiple key response points, the mean can be calculated first in the time dimension, and then the mean or weighted sum can be calculated in the measurement point dimension.

[0088] The physical consistency term includes two types of constraints: first, the predicted displacement after being processed by the time derivative operator D should be consistent with the predicted velocity; second, the predicted velocity after being processed by the time derivative operator D should be consistent with the predicted acceleration. These derivative consistency constraints reduce the risk of waveforms in the predicted sequence appearing similar but with mismatched derivative relationships.

[0089] The time derivative operator can be calculated using central differencing, forward differencing, backward differencing, smooth differencing, or automatic differentiation. For window boundaries or endpoints of the complete time sequence, the derivative consistency constraint term is calculated using one-sided differencing or neighbor extrapolation.

[0090] In the joint loss function, α u α v and α a β1 and β2 are used to control the relative contributions of the three types of data fitting errors, respectively, while λ controls the relative contributions of the two types of derivative consistency errors. λ controls the overall weight of the physical constraint loss relative to the data loss. As an optional setting, in the normalized training data, candidate values ​​of β1, β2, or λ can be set between 0.01 and 1 through a validation set grid search, selecting a combination that balances fitting error and derivative residuals.

[0091] The physical constraint loss does not require the model to completely satisfy a specific subgrade soil constitutive equation, but rather constrains the universal kinematic relationship between the three types of response states. Because this relationship is relatively decoupled from the subgrade soil designation, the thickness of the overlying structural layer, and the source of the input seismic motion, it has a good migration foundation under subgrade soil A, subgrade soil B, and subsequent subgrade conditions.

[0092] Example 5: Window Collapse, Response Recovery, and Predictive Output

[0093] like Figure 6 As shown, sequence organization processing is used to process long-time ground motion sequences. Let the original input time history length be T, the window length be W, and the window step size be P. Then, the original input time history is divided into multiple subsequences of length W. Adjacent subsequences can be non-overlapping or form overlapping windows according to a preset overlap rate. As an alternative sequence organization method to window folding, adjacent sampling points can be rearranged to the feature dimension, and the original time length can be restored through inverse rearrangement.

[0094] In a specific example, the complete input time history lasts for 50 seconds with a sampling interval of 0.02 seconds, and the complete sequence contains 2501 sampling points. It can be divided into five 10-second windows as subsequences, each containing approximately 500 sampling intervals. After each subsequence is independently input into the model for prediction, the complete P0 response time history is reconstructed at the output by concatenating or overlapping them in the original time order.

[0095] When windows overlap, predicted values ​​within the overlapping intervals can be weighted and fused using rectangular windows, triangular windows, Hanning windows, or other window functions to reduce abrupt changes in window boundaries. For non-overlapping segments of equal length, they can be directly spliced ​​in chronological order; for tail segments generated by zero-padding extensions, they should be truncated according to the original effective duration before calculating engineering indicators.

[0096] The model input for the application phase is the ground motion acceleration time history at the bottom of the subgrade system to be evaluated, and the model output is the acceleration response time history at the center point P0 of the top surface of the subgrade. The output time history can be further used to calculate peak ground acceleration, response spectrum, duration, critical moment response, inter-layer response amplification factor, or other engineering indicators required for rapid post-earthquake assessment of the subgrade.

[0097] When the prediction task involves multiple key response points, the multi-state output module simultaneously outputs the response prediction results for multiple key response points. Physical constraints can be applied separately among the three states at each measurement point, or additional constraints can be constructed by combining the response consistency between measurement points.

[0098] Example 6: Training Strategy and Engineering Scope

[0099] In this embodiment, model training can employ Adam or other gradient optimization algorithms to update network parameters, and combine learning rate decay, early stopping, gradient pruning, batch normalization, or input normalization strategies to improve training stability. For different subgrade soil codes or different input intensity ranges, independent training, conditional input training, or transfer learning can be selected based on the sample distribution.

[0100] During training, each model can employ the same data folding strategy, number of training epochs, and optimized hyperparameters to ensure a consistent basis for comparison between different network structures. Optional comparison models include GRU, LSTM, Phy-LSTM, KAN-LSTM, and the PhyKAN-LSTM of this invention; where Phy-LSTM is used to examine the contribution of physical constraints alone, KAN-LSTM is used to examine the contribution of KAN decoding alone, and PhyKAN-LSTM is used to examine the synergistic effect of both.

[0101] When there are significant differences between two types of subgrade soil samples, independent models can be established for subgrade soil A and subgrade soil B respectively to avoid distribution shift caused by mixed training across subgrade soils. When a unified model is required, subgrade soil code, density, Poisson's ratio, initial shear modulus, damping ratio, nonlinear adjustment coefficient, overlying structural layer thickness, and base layer thickness can be used as conditional features input into the network.

[0102] Model training can be monitored using the validation set to track total loss, acceleration RMSE, and derivative consistency residuals. Training should be stopped when the validation set loss stops decreasing. Gradient clipping and learning rate decay can be used when gradient explosion or local oscillations occur during training. These training strategies are designed to improve numerical stability and do not limit the model to using a specific optimizer or the number of training epochs.

[0103] To avoid unnecessarily limiting this invention to a specific numerical example, the geometric dimensions, subgrade soil designation, window length, number of hidden units, weighting coefficients, and training rounds in this invention can all be adjusted as implementation parameters. The above parameter combinations are used to illustrate optional implementation parameters of this invention and do not constitute a limitation on the scope of protection.

[0104] This invention is applicable to rapid proxy prediction of dynamic response for highway subgrades, airport runway foundations, urban road subgrades, site slab foundations, or similar subgrade engineering projects. For situations such as layered foundations, groundwater, freeze-thaw effects, local subgrade damage, vehicle load coupling, or three-dimensional edge effects, corresponding engineering parameters or numerical simulation results can be introduced during the sample construction stage, while using the three-state physical constraints and learnable activation decoding framework of this invention.

[0105] Example 7: Implementation Effect Diagrams and Results Description

[0106] The trained roadbed dynamic response prediction model can be used Figure 8 (a) to Figure 8 (f) Figure 9 (a) to Figure 9 (f) Figure 10 (a) Figure 10 (b) Figure 11 (a) and Figure 11 The implementation effect is demonstrated as shown in (b). Figure 8 (a) to Figure 8 (f) is used to display the representative acceleration response time histories of subgrade soil A and subgrade soil B under different input strength ranges. Figure 9 (a) to Figure 9 (f) is used to show the scatter plot correspondence between the predicted peak acceleration values ​​and the reference values. Figure 10 (a) and Figure 10 (b) Used to display the overall RMSE and correlation coefficients of different sequence models. Figure 11 (a) and Figure 11 (b) is used to display the PGA mean error, PGA maximum error, and CI95 coverage.

[0107] Figure 8 (a) Figure 8 (b) and Figure 8(c) The time histories of the top surface response of subgrade soil A under low, medium and high strength input ranges, respectively; Figure 8 (d) Figure 8 (e) and Figure 8 (f) Corresponds to the time histories of the top surface response of subgrade soil B under low, medium, and high strength input ranges, respectively. The time history curves in each subplot may include the reference response, the conventional LSTM predicted response, and the PhyKAN-LSTM predicted response of this invention. Different input strength ranges can be divided according to the bottom input PGA, and the range thresholds can be adjusted according to the distribution of engineering ground motion samples, and are not limited to fixed values.

[0108] Combination Figure 8 (a) to Figure 8 (f) Further analysis can be performed. Under low-intensity input, the response of the top surface of the subgrade is usually dominated by near-linear transmission and small-amplitude oscillation. Under medium-intensity input, the peak value and phase of the main vibration segment are more sensitive to the stiffness, damping and input frequency components of the subgrade soil. Under high-intensity input, subgrade soils with lower stiffness, such as subgrade soil A, are more likely to exhibit response amplification, local peak abrupt change and tail oscillation. At this time, relying solely on data fitting may result in insufficient tracking of local peak values ​​or inconsistent derivative relationships.

[0109] Figure 8 (a) to Figure 8 (f) It can also be used to check the quality of the model's recovery of the complete time history. For input sequences using window folding, if there are jumps at the window boundaries, truncation of the main segment, or discontinuous zero line regression at the tail segment, it indicates that the response recovery strategy still needs to be adjusted; if the prediction curve can maintain phase continuity, stable peak position, and smooth decay segment after multiple window splicing, it indicates that there is a good synergistic effect between sequence organization processing, response recovery processing, and physical constraint training.

[0110] Figure 9 (a) to Figure 9 The scatter plot in (f) is used to compare the consistency between the predicted peak acceleration and the reference peak acceleration. Figure 9 (a) to Figure 9 (c) The low, medium, and high strength input ranges corresponding to subgrade soil A. Figure 9 (d) to Figure 9 (f) Low, medium, and high strength input ranges corresponding to subgrade soil B. The closer the scatter points are to the consistency reference line, the closer the peak index prediction is to the reference result; the degree of dispersion of scatter points under different models, different subgrade soils, and different input strength ranges can be used to identify the applicable boundaries of the model.

[0111] Combination Figure 9 (a) to Figure 9(f) It can be further seen that the peak acceleration scatter plot is used not only to evaluate the average error but also to identify extreme samples. For the low-intensity range, the scatter plot usually clusters near the consistency reference line; for the medium-to-high-intensity range, the scatter plot may deviate or show increased dispersion, indicating that the nonlinearity of the subgrade soil, the difference in input frequency band, and the sample coverage range all affect the stability of the peak prediction. This scatter plot relationship can serve as an auxiliary basis for screening high-risk prediction samples in the engineering application stage.

[0112] Figure 10 (a) and Figure 10 The RMSE and correlation coefficient in (b) are used to evaluate both error magnitude and waveform correlation simultaneously. Figure 10 (a) Statistical results of the model corresponding to subgrade soil condition A. Figure 10 (b) Statistical results of the model under the corresponding subgrade soil condition B. RMSE reflects the overall deviation between the predicted time history and the reference time history, and the correlation coefficient R reflects the degree of synchronization between the predicted waveform and the reference waveform.

[0113] Figure 10 (a) and Figure 10 The comparison of different models in (b) can be used to decompose the contributions of technical modules. GRU and LSTM can be used as baselines for conventional recurrent sequence models; Phy-LSTM is used to observe the impact of kinematic derivative consistency constraints on waveform continuity and derivative residuals; KAN-LSTM is used to observe the impact of learnable activation decoding on nonlinear mapping; PhyKAN-LSTM is used to observe the combined effect of physical constraints and learnable activation decoding. Evaluation should combine RMSE, correlation coefficient, PGA error, and derivative consistency residuals for comprehensive judgment.

[0114] Figure 11 (a) and Figure 11 (b) The average error of PGA, the maximum error of PGA, and the CI95 coverage are used to evaluate the peak performance indicators of the project. Figure 11 (a) Statistical results corresponding to subgrade soil condition A. Figure 11 (b) Statistical results under subgrade soil condition B. The average error of PGA reflects the overall peak deviation level, the maximum error of PGA reflects the error risk of extreme samples or extrapolated samples, and the CI95 coverage rate is used to characterize the ability of the prediction uncertainty interval to cover the reference peak.

[0115] Combination Figure 11 (a) and Figure 11(b) It can be further explained that the average error of PGA is suitable for evaluating the overall deviation of the model on a batch of samples, while the maximum error of PGA is more suitable for engineering risk control; CI95 coverage is used to determine whether the error interval has sufficient envelope capability. When the coverage is too low, it can be improved by expanding the training sample range, grouping training according to the subgrade soil code, introducing input strength condition features, or adjusting the uncertainty estimation method.

[0116] Under subgrade soil condition A, due to the lower stiffness of the subgrade soil and the more pronounced damping and nonlinear effects, the response at the center point P0 of the subgrade top surface may exhibit strong amplification and a long attenuation segment. Physical constraint loss is used to suppress the derivative residuals between the three-state predictions, and KAN decoding is used to enhance the nonlinear expression capability near the peak and in the attenuation segment. Under subgrade soil condition B, the system response is relatively stable, but this invention can still maintain consistency between waveform recovery and peak index through joint loss.

[0117] From an engineering application perspective, Figure 8 (a) to Figure 8 (f) Used to observe the time history recovery and response persistence characteristics under a single ground motion. Figure 9 (a) to Figure 9 (f) Used to assess whether peak indicators are systematically overestimated or underestimated. Figure 10 (a) and Figure 10 (b) Used to select or verify the agent model structure. Figure 11 (a) and Figure 11 (b) Provides peak error and confidence interval information for rapid post-earthquake assessment or scheme comparison. The above combination of maps enables the model output to be expanded from a single time history prediction to a multi-level expression of results covering time history, peak value, statistical error, and uncertainty.

[0118] Figure 8 (a) to Figure 8 (f) Figure 9 (a) to Figure 9 (f) Figure 10 (a) Figure 10 (b) Figure 11 (a) and Figure 11 (b) The results and statistical methods shown are used to illustrate the expression of the implementation effect and the parameter evaluation method of the present invention. Since the sample source, subgrade soil parameters, seismic motion distribution, normalization strategy, training rounds and data partitioning methods will all affect the specific values, the relevant curves, scatter plots and bars are only used to show the implementation effect of the controlled calculation and are not used as numerical limitations on the scope of protection of the present invention.

Claims

1. A method for predicting the dynamic response of roadbed based on a physically constrained Kolmogorov-Arnold long short-term memory network, characterized in that, Includes the following steps: S1: Obtain the bottom input ground motion time history of the subgrade system and the target response time history of the key response points of the subgrade to form a training sample set; S2: Preprocess the input ground motion time history and target response time history in the training sample set, and obtain or construct three types of response state labels for the key response points: displacement, velocity and acceleration. S3: Perform sequence organization processing on the preprocessed input ground motion time history to obtain the input sequence for sequence modeling, and set the response recovery processing corresponding to the sequence organization processing; S4: Construct a neural network model, which includes a long short-term memory network sequential encoder, a Kolmogorov-Arnold learnable activation decoding module, and a multi-state output module; S5: Input the input sequence into the long short-term memory network temporal encoder to extract temporal hidden features, and generate displacement, velocity and acceleration response prediction sequences of the key response points through the Kolmogorov-Arnold learnable activation decoding module and multi-state output module; S6: Construct a data loss based on the displacement, velocity, and acceleration response prediction sequence and the response state label, and construct a physical constraint loss based on the kinematic derivative relationship between displacement, velocity, and acceleration; S7: Construct a joint loss function based on the data loss and physical constraint loss, and train the neural network model to obtain a trained roadbed dynamic response prediction model; S8: Input the ground motion time history of the bottom of the subgrade system to be predicted into the trained subgrade dynamic response prediction model, and output the dynamic response prediction results of the key response points after the response recovery processing.

2. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, The key response points include at least one of the following: the center point of the top surface of the subgrade, the edge point of the subgrade, the bottom point of the base course, the top surface point of the subgrade, or the measuring point inside the subgrade. The dynamic response prediction results include at least the acceleration response time history of the key response points.

3. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, In step S1, the training sample set comes from finite element numerical simulation, engineering monitoring records, or a combination of the two; the finite element numerical simulation includes the overlying structural layer, base course or subbase course, homogeneous subgrade soil, bottom seismic input boundary, and interlayer contact or continuous deformation relationship; the response state label is directly extracted from the numerical calculation results or engineering monitoring results, or constructed from the target response time history after integration, difference, or filtering correction.

4. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, In step S3, the sequence organization processing includes window folding processing, which divides the original input ground motion time history into multiple sub-sequences according to a preset window length and step size. The response recovery processing performs splicing, overlapping weighted averaging, or overlapping summation on the prediction results of multiple sub-sequences to obtain the response prediction result of the complete time length.

5. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, In step S3, the sequence organization processing includes sampling point rearrangement processing, which involves rearranging adjacent sampling points in the original input ground motion time history to the feature dimension, so as to shorten the time dimension and expand the single-step input feature dimension. The response recovery processing includes inverse rearrangement processing corresponding to the sampling point rearrangement processing.

6. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, In step S4, the Kolmogorov-Arnold learnable activation decoding module is set after the long short-term memory network temporal encoder and is used to map the temporal latent features into multiple learnable unary functions and their linear combinations. The learnable unary functions include one or more of spline functions, piecewise polynomial functions, or parameterized basis functions.

7. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, In step S6, the physical constraint loss takes kinematic derivative consistency as the constraint object, and includes at least a first derivative consistency constraint and a second derivative consistency constraint; the first derivative consistency constraint is used to constrain the consistency between the time derivative of the predicted displacement and the predicted velocity, and the second derivative consistency constraint is used to constrain the consistency between the time derivative of the predicted velocity and the predicted acceleration, and the construction of the physical constraint loss does not require a complete mass matrix, damping matrix or stiffness matrix as necessary input.

8. The method for predicting the dynamic response of a roadbed according to claim 7, characterized in that, The time derivative is obtained by one of the following methods: central difference method, forward difference method, backward difference method, smooth difference method, or automatic differentiation method, and is supplemented at the end points of the sequence by one-sided difference or end-point extrapolation.

9. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, In step S7, the joint loss function weights the data loss and physical constraint loss according to the response state magnitude, the training phase loss ratio, or a preset weight, and employs one or more training strategies, such as learning rate decay, early stopping, gradient pruning, or normalization.

10. The method for predicting the dynamic response of a roadbed according to claim 1, characterized in that, It also satisfies at least one of the following conditions: the neural network model also receives one or more conditional features among the following: subgrade soil code, density, Poisson's ratio, initial shear modulus, damping parameter, overlying structural layer thickness, base course thickness, subgrade depth, or overlying structural layer material parameters; the dynamic response prediction results also include one or more of the following indicators among multiple key response points: acceleration response time history, peak acceleration, response spectrum, peak displacement, or duration index.

Images