METHOD AND SYSTEM FOR IDENTIFYING TIME-VARYING CHARACTERISTICS OF HEAVY-LOAD VEHICLE SUSPENSION

Abstract

A method and system are provided for identifying time-varying suspension characteristics of heavy-load vehicles. The method includes collecting sequential control state data of a mining truck using sensors, predicting parameter-related factors through a deep learning network, estimating suspension stiffness and damping coefficients via a linear dynamic model considering longitudinal-vertical coupling, and predicting future system states through a nonlinear dynamic model based on the estimated parameters and learned factors. According to the method, a deep learning network is integrated into a physical model of the mining truck, an accurate longitudinal-vertical dynamical model of the mining truck is established, accurate suspension parameters are identified, the stiffness damping time-varying characteristics of the suspension of the mining truck are given through a physical model-data driving method, and the model has certain interpretability and generalization; the rigidity and damping of the four suspensions can be obtained only through sprung information.

Claims

1. A method for identifying time-varying characteristics of heavy-load vehicle suspension, comprising a non-transitory computer readable medium operable on a computer with memory for the method for identifying time-varying characteristics of heavy-load vehicle suspension, and comprising program instructions for executing the following steps of: using sensors to collect the sequence state of mining truck operations; using a deep learning network to predict parameter-related factors based on the sequence state, and meanwhile; using a linear dynamics model considering longitudinal-vertical coupling effects roughly estimates the stiffness and damping coefficients as key suspension parameters; using a nonlinear dynamics model considering longitudinal-vertical coupling effects then predicts the state at the next time step based on the parameter-related factors and key suspension parameters; the state loss Lp between the predicted state and the target state detected by the sensors serves as the loss error term for the deep learning network to update the network parameters; the linear dynamics model considering longitudinal-vertical coupling effects includes a body with mass m.sub.c and moment of inertia Ic, and front and rear axles with unsprung masses m.sub.tf and m.sub.tr respectively, and the suspension forces transmitted to the body from the front and rear axles are F.sub.f and F.sub.r, and the road excitations for the front and rear tires are z.sub.qf and z.sub.qr, and the degrees of freedom (DoF) of the model include the vertical displacement of the center of gravity (CoG)z.sub.c, the pitch angle .sub.c, and the vertical displacements of the front and rear unsprung masses z.sub.tf and z.sub.tr, and F.sub.c represents the inertial force acting on the CoG due to acceleration a and velocity v along the x-axis, and M.sub.c represents the moment of inertia, and the suspension stiffness and damping coefficients are k and c, with subscripts f and r representing the front and rear suspensions, and subscript t representing the tire, and the distances from the CoG to the front and rear axles are a.sub.c and b.sub.c respectively, and the vertical distance between the CoG and the pitch center PC is h; and improving vehicle stability under heavy-duty and off-road conditions while prolonging suspension system longevity based on the method for identifying time-varying characteristics of heavy-load vehicle suspension through a deep learning-physics hybrid-driven network.

2. The method for identifying time-varying characteristics of heavy-load vehicle suspension according to claim 1, characterized in that the vertical acceleration and velocity at the CoG position are calculated by common calculation apparatus based on the nonlinear mapping relationship of the IMU $\{\begin{array}{c}{\stackrel{.}{z}}_c=A\left({}^i{u}^i\right)\\ {\stackrel{.Math.}{z}}_c=A\left({}^i{u}^i+{}^j{}^i{u}^i\right)\end{array},$ comprises the following steps: i) according to the Euler theorem of rigid body kinematics, the motion of a point in the absolute coordinate system r is decomposed into the translational displacement R in the relative coordinate system plus the vector sum of rotation around the base point: $r=R+u$ $u=A{u}^i,u^i={\left[x^i,y^i,z^i\right]}^T,$ where: A is the coordinate rotation matrix, u.sub.i is the projection on the coordinate axis after rotation, and OXYZ is the inertial reference frame (absolute reference frame) fixed on the ground; O.sub.iX.sub.mY.sub.mZ.sub.m is the coordinate system (relative coordinate system) fixed on the IMU, with the three-axis direction the same as the inertial reference frame; O.sub.iX.sub.iY.sub.iZ.sub.i represents the coordinate system definition of the IMU itself, and P is the point to be measured, u is the vector from the coordinate axis reference point to the measured point; ii) according to the actual situation, reasonably determine the rotation order of the coordinate axis and calculate the rotation matrix: according to vehicle kinematics knowledge, the roll angle and pitch angle are relatively small compared to the yaw angle, and therefore, the rotation order when selecting the coordinate axis is: first rotate around the Z axis (corresponding to the yaw angle .sub.c); then rotate around the Y axis (corresponding to the pitch angle .sub.c); finally rotate around the X axis (corresponding to the roll angle .sub.c), and Obtain the rotation matrix expression A (.sub.c, .sub.c, .sub.c)=A.sub.A.sub.A.sub., where: the positive direction of rotation is the right-hand spiral direction, and the expressions of the three rotation matrices are $A_{}=\left(\begin{array}{ccc}\cos {}_c& -\sin {}_c& 0\\ \sin {}_c& \cos {}_c& 0\\ 0& 0& 1\end{array}\right);$ $A_{}=\left(\begin{array}{ccc}\cos {}_c& 0& \sin {}_c\\ 0& 1& 0\\ -\sin {}_c& 0& cos{}_c\end{array}\right);$ $A_{}=\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos {}_c& -\sin {}_c\\ 0& \sin {}_c& \cos {}_c\end{array}\right)$ iii) calculate the speed of point P, and take the derivative of both sides of the expressions of the three rotation matrices with respect to time, and the rotation matrix is regular, and its derivative with respect to time is calculated through linear algebra to obtain {dot over (A)}={tilde over ()}A, where: the angular velocity matrix of each axis is an antisymmetric matrix $\stackrel{~}{}=\left(\begin{array}{ccc}0& -{}_z& {}_y\\ {}_z& 0& -{}_x\\ -{}_y& {}_x& 0\end{array}\right),$ obtain the speed expression of the measured point $\begin{array}{c}\stackrel{.}{r}=\stackrel{}{R}+\stackrel{}{A}{u}^i=\stackrel{}{R}+\stackrel{~}{}A{u}^i\\ =\stackrel{.}{R}+u=\stackrel{.}{R}+A\left({}^i{u}^i\right)\end{array}$ and the acceleration expression of the measured point (the angular acceleration of each axis) $\begin{array}{c}\stackrel{.Math.}{r}=\stackrel{.Math.}{R}+au+u\\ =\stackrel{.Math.}{R}+A\left({}^i{u}^i+{}^i{}^i{u}^i\right)\end{array},$ that is, calculate the acceleration (including angular acceleration) of any point rigidly connected to the IMU installation position, and thus, the motion of the vehicle body is obtained.

3. The method for identifying time-varying characteristics of heavy-load vehicle suspension according to claim 1, characterized in that the key suspension parameters, namely the suspension stiffness and damping coefficients, are calculated by common calculation apparatus in the following steps: take the front suspension as an example, under the input of force F.sub.f, obtain the vertical displacement of the sprung mass z.sub.f according to the transfer relationship from the suspension force to the suspension displacement in the dynamics system of the mining truck, and use the System identification method to calculate the suspension stiffness, and calculate the displacement of the unsprung position, specifically including: i) according to the Laplace transform relationship between force and displacement $\{\begin{array}{c}{F}_f=\left(k+cs\right)\left(Z_f-Z_{\mathrm{tf}}\right)\\ F_f=\left(k_{\mathrm{tf}}+m_{\mathrm{tf}}{s}^2\right)Z_{\mathrm{tf}}\end{array},$ where: the above are all defined, and the script is the Laplace transform of the corresponding variable, and ii) by eliminating the vertical displacement of the unsprung mass, the transfer relationship from the suspension force to the suspension displacement in the dynamics system of the mining truck $\frac{Z_f}{F_f}=\frac{m_{\mathrm{tf}}{s}^2+c_{f}s+\left(k_{\mathrm{tf}}+k_f\right)}{m_{\mathrm{tf}}{c}_f{s}^3+m_{\mathrm{tf}}{k}_f{s}^2+k_{\mathrm{tf}}{c}_{f}s+k_{\mathrm{tf}}{k}_f},$ where: s is the time domain differential operator, and According to the measured real-time timing signal of the mining truck, use the continuous-time system simplified refined instrumental variable method (SRIVC) for system identification.sup.[22], and rewrite the system transfer equation into the standard form with limited continuous-time white noise e(t) $Z_f=G\left(s\right)F_f+e\left(t\right)$ $G\left(s\right)=\frac{B\left(s\right)}{A\left(s\right)}=\frac{b_2{s}^2+b_{1}s+b_0}{a_3{s}^3+a_2{s}^2+a_{1}s+1},$ where: the parameter vector to be identified is =[a.sub.1, a.sub.2, a.sub.3, b.sub.0, b.sub.1, b.sub.2].sup.T; iii) using a continuous-time state variable filter to perform low-pass filtering on the entire equation to obtain the pre-filtered time derivatives of the input and output required for identification, and Since the filtered regression vector .sub.f is related to the noise vector, the identification result is asymptotically biased, and to eliminate the influence of the noise vector, define the instrumental vector at step $$\begin{gathered} {}_f\left(t_k\right)={\frac{1}{A_j\left(s\right)}\left[-s{Z}_f,.Math.,-s^3{Z}_f,F_f,.Math.,s^2{F}_f\right]}_{f_k}^T \\ {}_f\left(t_k\right){\frac{F_f\left(t_k\right)}{A_j\left(s\right)}\left[\frac{-s{B}_j\left(s\right)}{A_j\left(s\right)},.Math.,\frac{-s^3{B}_j\left(s\right)}{A_j\left(s\right)},1,.Math.,s^2\right]}^T, \end{gathered}$$ $Z_f^f\left(t_k\right)=\frac{1}{A_j\left(p\right)}{Z}_f\left(t_k\right)$ calculate the parameter estimate at step j+1 ${}_{j+1}={\left[\underset{k=1}{\overset{N}{.Math.}}{}_f\left(t_k\right){}_f^T\left(t_k\right)\right]}^{-1}\left[\underset{k=1}{\overset{N}{.Math.}}{}_f\left(t_k\right)Z_f^f\left(t_k\right)\right];$ continue iterating until the error between the identification result and the previous identification result is small enough, that is $\frac{{.Math.{}_{j+1}-{}_j.Math.}_2}{{.Math.{}_j.Math.}_2}<;$ iv) according to the identified system transfer function and the original input-output relationship, use the normalized root mean square error (NRMSE) and the final prediction error (FPE) as evaluation indicators to evaluate the accuracy of the model.

4. The method for identifying time-varying characteristics of heavy-load vehicle suspension according to claim 1, characterized in that the state at the next time step is identified by combining data-driven and neural network methods to identify suspension parameters, which including the following steps: i) take the front suspension as an example, perform backward differentiation on $0=k_{\mathrm{tf}}^k\left(z_{\mathrm{qf}}^k-z_{\mathrm{tf}}^k\right)-F_f^k-m_{\mathrm{tf}}\frac{{\stackrel{.}{z}}_{\mathrm{tf}}^k-{\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}}{t^k}$ to obtain the discrete system expression of the suspension at time k, where: the superscripts k and k1 represent the current state and the last state respectively, and the sampling time difference of each state; since the displacement of the unsprung mass is obtained through parameter identification, its differentiation is not accurate, and it needs to be combined with other equations for identification, and the speed of the unsprung mass at time k is $F_f^k=k_{f,\mathrm{new}}\left(z_{\mathrm{tf}}^k-z_c^k+a_c{}_c^k\right)+c_{f,\mathrm{new}}\left({\stackrel{.}{z}}_{\mathrm{tf}}^k-{\stackrel{.}{z}}_c^k+a_c{\stackrel{.}{}}_c^k\right),$ $F_f^{k-1}=k_{f,\mathrm{new}}\left(z_{\mathrm{tf}}^{k-1}-z_c^{k-1}+a_c{}_c^{k-1}\right)+c_{f,\mathrm{new}}\left({\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}-{\stackrel{.}{z}}_c^{k-1}+a_c{\stackrel{.}{}}_c^{k-1}\right),$ where: the subscript new represents the parameters modified by the neural network, and When the suspension force at time k1 is introduced, the entire algorithm will iterate, causing the accumulation of errors each time, and Introduce the forward differentiation of the body DoF to represent the speed of the unsprung mass at time $k-1F_f^{k-1}=\frac{m_c{\stackrel{.Math.}{z}}_c^{k-1}{b}_c-I_c{\stackrel{.Math.}{}}_c^{k-1}-h.Math.m_c{a}^{k-1}}{L_c},$ ${\stackrel{.Math.}{}}_c^{k-1}=\frac{{\stackrel{.}{}}_c^k-{\stackrel{.}{}}_c^{k-1}}{t^{k-1}},{\stackrel{.Math.}{z}}_c^{k-1}=\frac{{\stackrel{.}{z}}_c^k-{\stackrel{.}{z}}_c^{k-1}}{t^{k-1}},$ obtain the speed of the unsprung mass and then obtain $\{\begin{array}{c}{\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}=\frac{F_f^{k-1}-k_{f,\mathrm{new}}\left(z_{\mathrm{tf}}^{k-1}-z_c^{k-1}+a_c{}_c^{k-1}\right)}{c_{f,\mathrm{new}}}+{\stackrel{.}{z}}_c^{k-1}-a_c{\stackrel{.}{}}_c^{k-1}\\ {\stackrel{.}{z}}_{\mathrm{tf}}^k=\left[k_{\mathrm{tf}}^k\left(z_{\mathrm{qf}}^k-z_{\mathrm{tf}}^k\right)-F_f^k\right]\frac{t^k}{m_{\mathrm{tf}}}+{\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}\end{array},$ and solve the function to obtain the force Ff transmitted by the front suspension at time k:Ff $f\left(F_f^k\right)=0;$ ii) as to the rear suspension, the equation $\{\begin{array}{c}{z}_{\mathrm{tr}}^{k-1}=\frac{F_r^{k-1}-k_{r,\mathrm{new}}\left(z_{\mathrm{tr}}^{k-1}-z_c^{k-1}-b_c{}_c^{k-1}\right)}{c_{r,\mathrm{new}}}+{\stackrel{.}{z}}_c^{k-1}+b_c{\stackrel{.}{}}_c^{k-1}\\ z_{\mathrm{tr}}^k=\left[k_{\mathrm{tr}}^k\left(z_{qr}^k-z_{\mathrm{tr}}^k\right)-F_r^k\right]\frac{t^k}{m_{\mathrm{tr}}}+z_{\mathrm{tr}}^{k-1}\\ F_r^k=k_{r,\mathrm{new}}\left(z_{\mathrm{tr}}^k-z_c^k-b_c{}_c^k\right)+c_{r,\mathrm{new}}\left({\stackrel{.}{z}}_{\mathrm{tr}}^k-{\stackrel{.}{z}}_c^k-b_c{\stackrel{.}{}}_c^k\right)\end{array}$ can be obtained, and the F.sub.r at time k is calculated, and therefore, the motion of the vehicle body is $\{\begin{array}{c}{\stackrel{.Math.}{}}_c^k=\frac{F_r^k{b}_c-F_f^k{a}_c-h.Math.m_c{a}^k}{I_c}\\ {\stackrel{.Math.}{z}}_c^k=\frac{F_f^k+F_r^k}{m_c}\end{array},$ and according to the forward differentiation method, the predicted state of the vehicle body is summarized as $\{\begin{array}{c}{\stackrel{.}{}}_c^{k+1}={\stackrel{.Math.}{}}_c^k{t}^k+{\stackrel{.}{}}_c^k\\ {\stackrel{.}{z}}_c^{k+1}={\stackrel{.Math.}{z}}_c^k{t}^k+{\stackrel{.}{z}}_c^k\end{array},$ obtaining the iterative method of the mining truck body parameters in discrete time, and iii) as the front and rear suspension parameters of the mining truck show nonlinear characteristics under different working conditions, and On the basis of previous parameter identification, use the neural network for refined modeling, that is $\{\begin{array}{c}{k}_{f,\mathrm{new}}=P_{kf}.Math.k_f,c_{f,\mathrm{new}}=P_{cf}.Math.c_f\\ k_{r,\mathrm{new}}=P_{kr}.Math.k_r,c_{r,\mathrm{new}}=P_{cr}.Math.c_r\end{array},$ the updated iterative parameters are, and the neural network obtains N.sub.out=[P.sub.kf, P.sub.kr, P.sub.cf, P.sub.cr] according to the time series signal N.sub.in=[z.sub.c, .sub.c, {umlaut over (z)}.sub.c, .sub.c, {dot over ()}.sub.c, {umlaut over ()}.sub.c, a, t, z.sub.tf, z.sub.tr].

5. The method for identifying time-varying characteristics of heavy-load vehicle suspension according to claim 1, characterized in that the deep learning network includes a three-layer LSTM network layer and a fully connected network layer, where: each LSTM layer has 128 hidden units, and the fully connected network layer outputs the suspension stiffness and damping correction coefficients.

Description

BRIEF DESCRIPTION OF DRAWINGS

[0018] FIG. 1 shows a schematic diagram of the vertical apparent dynamics identification method for mining trucks.

[0019] FIG. 2 shows a schematic diagram of the nonlinear dynamics model considering longitudinal-vertical coupling effects.

[0020] FIG. 3 shows a schematic diagram of the hybrid modeling architecture.

[0021] FIG. 4 shows a schematic diagram of the original signals collected by the sensors.

[0022] FIG. 5 shows a schematic diagram of the mining truck speed in each test section.

[0023] FIG. 6 shows a comparison of vertical acceleration at the CoG (Part 1).

[0024] FIG. 7 shows a comparison of vertical acceleration at the CoG (Part 5).

[0025] FIG. 8 shows a comparison of front suspension vertical velocity.

[0026] FIG. 9 shows a comparison of rear suspension vertical velocity.

[0027] FIG. 10 shows the variation pattern of front and rear suspension stiffness.

[0028] FIG. 11 shows the variation pattern of front and rear suspension damping.

DETAILED DESCRIPTION OF INVENTION

[0029] As shown in FIG. 1, the embodiment involves a method for identifying time-varying characteristics of heavy-load vehicle suspension. Sensors collect the sequence state of mining truck operations, and a deep learning network predicts parameter-related factors based on the sequence state. At the same time, a linear dynamics model considering longitudinal-vertical coupling effects roughly estimates the suspension stiffness and damping coefficients as key suspension parameters. A nonlinear dynamics model considering longitudinal-vertical coupling effects then predicts the state at the next time step based on the parameter-related factors and key suspension parameters.

[0030] As shown in FIG. 2, the nonlinear dynamics model considering longitudinal-vertical coupling effects has a defined positive direction. It includes a body with mass m.sub.c and moment of inertia Ic, and front and rear axles with unsprung masses m.sub.tf and m.sub.tr respectively. The suspension forces transmitted to the body from the front and rear axles are F.sub.f and F.sub.r, and the road excitations for the front and rear tires are z.sub.qf and z.sub.qr. The degrees of freedom (DoF) of the model include the vertical displacement of the center of gravity (CoG) z.sub.c, the pitch angle .sub.c, and the vertical displacements of the front and rear unsprung masses z.sub.tf and z.sub.tr. F.sub.c represents the inertial force acting on the CoG due to acceleration a and velocity v along the x-axis, and M.sub.c represents the moment of inertia. The suspension stiffness and damping coefficients are k and c, with subscripts f and r representing the front and rear suspensions, and subscript t representing the tire. The distances from the CoG to the front and rear axles are a.sub.c and b.sub.c respectively, and the vertical distance between the CoG and the pitch center PC is h.

[0031] According to D'Alembert's principle, the longitudinal-vertical dynamics model of the mining truck is

[00001]$\{\begin{array}{c}0=F_f+F_r-m_c{\stackrel{.Math.}{z}}_c\\ 0=F_r{b}_c-F_f{a}_c-h.Math.m_{c}a-I_c{\stackrel{.Math.}{}}_c\\ 0=k_{\mathrm{tf}}\left(z_{\mathrm{qf}}-z_{\mathrm{tf}}\right)-F_f-m_{\mathrm{tf}}{\stackrel{.Math.}{z}}_{\mathrm{tf}}\\ 0=k_{\mathrm{tr}}\left(z_{qr}-z_{\mathrm{tr}}\right)-F_r-m_{\mathrm{tr}}{\stackrel{.Math.}{z}}_{\mathrm{tr}}\end{array},$ [0032] where the force exerted by the body on the suspension is

[00002]$\{\begin{array}{c}{F}_f=k_f\left(z_{\mathrm{tf}}-z_c+a_c{}_c\right)+c_f\left({\stackrel{.}{z}}_{\mathrm{tf}}-{\stackrel{.}{z}}_c+a_c{\stackrel{.}{}}_c\right)\\ F_r=k_r\left(z_{tr}-z_c-b_c{}_c\right)+c_r\left({\stackrel{.}{z}}_{tr}-{\stackrel{.}{z}}_c-b_c{\stackrel{.}{}}_c\right)\end{array}.$

[0033] Preferably, the vertical acceleration an velocity at the CoG position are calculated based on the nonlinear mapping relationship of the IMU

[00003]$\{\begin{array}{c}{\stackrel{.}{z}}_c=A\left({}^i{u}^i\right)\\ {\stackrel{.Math.}{z}}_c=A\left({}^i{u}^i+{}^i{}^i{u}^i\right)\end{array},$

specifically including:

[0034] 1) According to the Euler theorem of rigid body kinematics, the motion of a point in the absolute coordinate system r is decomposed into the translational displacement R in the relative coordinate system plus the vector sum of rotation around the base point:

[00004]$r=R+u$$u=A{u}^i,u^i={\left[x^i,y^i,z^i\right]}^T,$

where: A is the coordinate rotation matrix, u.sub.i is the projection on the coordinate axis after rotation. OXYZ is the inertial reference frame (absolute reference frame) fixed on the ground; O.sub.iX.sub.mY.sub.mZ.sub.m is the coordinate system (relative coordinate system) fixed on the IMU, with the three-axis direction the same as the inertial reference frame; O.sub.iX.sub.iY.sub.iZ.sub.i represents the coordinate system definition of the IMU itself. P is the point to be measured, u is the vector from the coordinate axis reference point to the measured point.

[0035] 2) According to the actual situation, reasonably determine the rotation order of the coordinate axis and calculate the rotation matrix: According to vehicle kinematics knowledge, the roll angle and pitch angle are relatively small compared to the yaw angle. Therefore, the rotation order when selecting the coordinate axis is: first rotate around the Z axis (corresponding to the yaw angle .sub.c); then rotate around the Y axis (corresponding to the pitch angle .sub.c); finally rotate around the X axis (corresponding to the roll angle .sub.e). Obtain the rotation matrix expression A (.sub.c, .sub.c, .sub.c)=A.sub.A.sub.A.sub., where: the positive direction of rotation is the right-hand spiral direction, and the expressions of the three rotation matrices are

[00005]$A_{}=\left(\begin{array}{ccc}\cos {}_c& -\sin {}_c& 0\\ \sin {}_c& \cos {}_c& 0\\ 0& 0& 1\end{array}\right);$$A_{}=\left(\begin{array}{ccc}\cos {}_c& 0& \sin {}_c\\ 0& 1& 0\\ -\sin {}_c& 0& \cos {}_c\end{array}\right);$

3) Calculate the speed of

[00006]$A_{}=\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos {}_c& -\sin {}_c\\ 0& \sin {}_c& \cos {}_c\end{array}\right)$

point P, and take the derivative of both sides of the expressions of the three rotation matrices with respect to time. The rotation matrix is regular, and its derivative with respect to time is calculated through linear algebra to obtain {dot over (A)}={tilde over ()}A, where: the angular velocity matrix of each axis is an antisymmetric matrix

[00007]$\stackrel{~}{}=\left(\begin{array}{ccc}0& -{}_z& {}_y\\ {}_z& 0& -{}_x\\ -{}_y& {}_x& 0\end{array}\right),$

obtain the speed expression of the measured point

[00008]$\begin{array}{c}\stackrel{.}{r}=\stackrel{}{R}+\stackrel{}{A}{u}^i=\stackrel{}{R}+\stackrel{~}{}A{u}^i\\ =\stackrel{.}{R}+u=\stackrel{.}{R}+A\left({}^i{u}^i\right)\end{array}$

and the acceleration expression of the measured point (the angular acceleration of each axis)

[00009]$\begin{array}{c}\stackrel{.Math.}{r}=\stackrel{.Math.}{R}+u+u\\ =\stackrel{.Math.}{R}+A\left({}^i{u}^i+{}^i{}^i{u}^i\right)\end{array},$

that is, calculate the acceleration (including angular acceleration) of any point rigidly connected to the IMU installation position.

[0036] Thus, the motion of the vehicle body is obtained, and the longitudinal-vertical dynamics model of the mining truck is completed.

[0037] The key suspension parameters, namely the suspension stiffness and damping coefficients, are obtained in the following way: Take the front suspension as an example, under the input of force F.sub.f, obtain the vertical displacement of the sprung mass z.sub.f according to the transfer relationship from the suspension force to the suspension displacement in the dynamics system of the mining truck, and use the system identification method to calculate the suspension stiffness, and calculate the displacement of the unsprung position, specifically including:

[0038] 1) According to the Laplace transform relationship between force and displacement

[00010]$\{\begin{array}{c}{F}_f=\left(k+cs\right)\left(Z_f-Z_{\mathrm{tf}}\right)\\ F_f=\left(k_{\mathrm{tf}}+m_{\mathrm{tf}}{s}^2\right)Z_{\mathrm{tf}}\end{array},$

where: the above are all defined, and the script is the Laplace transform of the corresponding variable.

[0039] 2) By eliminating the vertical displacement of the unsprung mass, the transfer relationship from the suspension force to the suspension displacement in the dynamics system of the mining truck

[00011]$\frac{Z_f}{F_f}=\frac{m_{\mathrm{tf}}{s}^2+c_{f}s+\left(k_{\mathrm{tf}}+k_f\right)}{m_{\mathrm{tf}}{c}_f{s}^3+m_{\mathrm{tf}}{k}_f{s}^2+k_{\mathrm{tf}}{c}_{f}s+k_{\mathrm{tf}}{k}_f},$

where: s is the time domain differential operator. According to the measured real-time timing signal of the mining truck, use the continuous-time system simplified refined instrumental variable method (SRIVC) for system identification.sup.[22], and rewrite the system transfer equation into the standard form with limited continuous-time white noise e(t)

[00012]$Z_f=G\left(s\right)F_f+e\left(t\right)$$G\left(s\right)=\frac{B\left(s\right)}{A\left(s\right)}=\frac{b_2{s}^2+b_{1}s+b_0}{a_3{s}^3+a_2{s}^2+a_{1}s+1},$

where: the parameter vector to be identified is =[a.sub.1, a.sub.2, a.sub.3, b.sub.0, b.sub.1, b.sub.2].sup.T;

[0040] 3) Use the continuous-time state variable filter to perform low-pass filtering on the entire equation to obtain the pre-filtered time derivatives of the input and output required for identification. Since the filtered regression vector .sub.f is related to the noise vector, the identification result is asymptotically biased. To eliminate the influence of the noise vector, define the instrumental vector at step j

[00013]${}_f\left(t_k\right)={\frac{1}{A_j\left(s\right)}\left[-s{Z}_f,.Math.,-s^3{Z}_f,F_f,.Math.,s^2{F}_f\right]}_{t_k}^T$${}_f\left(t_k\right)={\frac{F_f\left(t_k\right)}{A_j\left(s\right)}\left[\frac{-s{B}_j\left(s\right)}{A_j\left(s\right)},.Math.,\frac{-s^3{B}_j\left(s\right)}{A_j\left(s\right)},1,.Math.,s^2\right]}^T,$

calculate the parameter estimate at step j+1

[00014]$Z_f^f\left(t_k\right)=\frac{1}{A_j\left(p\right)}{Z}_f\left(t_k\right)$${}_{j+1}={\left[\underset{k=1}{\overset{N}{.Math.}}{}_f\left(t_k\right){}_f^T\left(t_k\right)\right]}^{-1}\left[\underset{k=1}{\overset{N}{.Math.}}{}_f\left(t_k\right)Z_f^f\left(t_k\right)\right];$

continue iterating until the error between the identification result and the previous identification result is small enough, that is

[00015]$\frac{{.Math.{}_{j+1}-{}_j.Math.}_2}{{.Math.{}_j.Math.}_2}<;$

[0041] 4) According to the identified system transfer function and the original input-output relationship, use the normalized root mean square error (NRMSE) and the final prediction error (FPE) as evaluation indicators to evaluate the accuracy of the model.

[0042] In this embodiment, the high-order term is selected as the calculation method of parameter regression, and the calculation methods of cf=a3/b2, kf=a2/b2 are used. Considering that this test is conducted under the same initial conditions, the suspension stiffness parameters should also be the same. Therefore, the least squares method is used to calculate the suspension stiffness parameters under all working conditions, and the average value of the error evaluation indicators is given. Similarly, the rear suspension also uses the same calculation process. In addition, it also provides a basis for the establishment and training of subsequent network models, and is directly used for the calculation of physical models.

[0043] The state at the next time step is identified by combining data-driven and neural network methods to identify suspension parameters, specifically including:

[0044] 1) Take the front suspension as an example, perform backward differentiation on {umlaut over (z)}.sub.tf [0045] to obtain the discrete system expression of the suspension at time k

[00016]$0=k_{\mathrm{tf}}^k\left(z_{\mathrm{qf}}^k-z_{\mathrm{tf}}^k\right)-F_f^k-m_{\mathrm{tf}}\frac{{\stackrel{.}{z}}_{\mathrm{tf}}^k-{\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}}{t^k},$ where: the superscripts k and k1 represent the current state and the last state respectively. The sampling time difference of each state; since the displacement of the unsprung mass is obtained through parameter identification, its differentiation is not accurate, and it needs to be combined with other equations for identification. The speed of the unsprung mass at time k is

[00017]$F_f^k=k_{f,\mathrm{new}}\left(z_{\mathrm{tf}}^k-z_c^k+a_c{}_c^k\right)+c_{f,\mathrm{new}}\left({\stackrel{.}{z}}_{\mathrm{tf}}^k-{\stackrel{.}{z}}_c^k+a_c{\stackrel{.}{}}_c^k\right),$$F_f^{k-1}=k_{f,\mathrm{new}}\left(z_{\mathrm{tf}}^{k-1}-z_c^{k-1}+a_c{}_c^{k-1}\right)+c_{f,\mathrm{new}}\left({\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}-{\stackrel{.}{z}}_c^{k-1}+a_c{\stackrel{.}{}}_c^{k-1}\right),$ where: the subscript new represents the parameters modified by the neural network. When the suspension force at time k1 is introduced, the entire algorithm will iterate, causing the accumulation of errors each time. Introduce the forward differentiation of the body DoF to represent the speed of the unsprung mass at time k1

[00018]$F_f^{k-1}=\frac{m_c{\stackrel{.Math.}{z}}_c^{k-1}{b}_c-I_c{\stackrel{.Math.}{}}_c^{k-1}-h.Math.m_c{a}^{k-1}}{L_c},$${\stackrel{.Math.}{}}_c^{k-1}=\frac{{\stackrel{.}{}}_c^k-{\stackrel{.}{}}_c^{k-1}}{t^{k-1}},{\stackrel{.Math.}{z}}_c^{k-1}=\frac{{\stackrel{}{z}}_c^k-{\stackrel{}{z}}_c^{k-1}}{t^{k-1}},$ [0046] obtain the speed of the unsprung mass

[00019]${\stackrel{}{z}}_{\mathrm{tf}}^{k-1}$ and then obtain

[00020]$\{\begin{array}{c}{\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}=\frac{F_f^{k-1}-k_{f,\mathrm{new}}\left(z_{\mathrm{tf}}^{k-1}-z_c^{k-1}+a_c{}_c^{k-1}\right)}{c_{f,\mathrm{new}}}+{\stackrel{.}{z}}_c^{k-1}-a_c{\stackrel{.}{}}_c^{k-1}\\ {\stackrel{.}{z}}_{\mathrm{tf}}^k=\left[k_{\mathrm{tf}}^k\left(z_{\mathrm{af}}^k-z_{\mathrm{tf}}^k\right)-F_f^k\right]\frac{t^k}{m_{\mathrm{tf}}}+{\stackrel{.}{z}}_{\mathrm{tf}}^{k-1}\end{array},$ and solve the function to obtain the force transmitted by the front suspension at time k: Ff:

[00021]$f\left(F_f^k\right)=0;$

[0047] 2) For the rear suspension, the equation

[00022]$\{\begin{array}{c}{z}_{\mathrm{tf}}^{k-1}=\frac{F_r^{k-1}-k_{r,\mathrm{new}}\left(z_{\mathrm{tr}}^{k-1}-z_c^{k-1}-b_c{}_c^{k-1}\right)}{c_{r,\mathrm{new}}}+{\stackrel{.}{z}}_c^{k-1}+b_c{\stackrel{.}{}}_c^{k-1}\\ z_{\mathrm{tr}}^k=\left[k_{\mathrm{tr}}^k\left(z_{\mathrm{qf}}^k-z_{\mathrm{tr}}^k\right)-F_r^k\right]\frac{t^k}{m_{\mathrm{tr}}}+z_{\mathrm{tr}}^{k-1}\\ F_r^k=k_{r,\mathrm{new}}\left(z_{\mathrm{tr}}^k-z_c^k-b_c{}_c^k\right)+c_{r,\mathrm{new}}\left({\stackrel{.}{z}}_{\mathrm{tr}}^k-{\stackrel{.}{z}}_c^k-b_c{\stackrel{.}{}}_c^k\right)\end{array}$

can be obtained, and the F.sub.r at time k is calculated. Therefore, the motion of the vehicle body is

[00023]$\{\begin{array}{c}{\stackrel{.Math.}{}}_c^k=\frac{F_r^k{b}_c-F_f^k{a}_c-h.Math.m_c{a}^k}{I_c}\\ {\stackrel{.Math.}{z}}_c^k=\frac{F_f^k+F_r^k}{m_c}\end{array},$

and according to the forward differentiation method, the predicted state of the vehicle body is summarized as

[00024]$\{\begin{array}{c}{\stackrel{.}{}}_c^{k+1}={\stackrel{.Math.}{}}_c^k{t}^k+{\stackrel{.}{}}_c^k\\ {\stackrel{.}{z}}_c^{k+1}={\stackrel{.Math.}{z}}_c^k{t}^k+{\stackrel{.}{z}}_c^k\end{array},$

obtaining the iterative method of the mining truck body parameters in discrete time.

[0048] 3) The front and rear suspension parameters of the mining truck show nonlinear characteristics under different working conditions. On the basis of previous parameter identification, use the neural network for refined modeling, that is, the updated iterative parameters

[00025]$\{\begin{array}{c}{k}_{f,\mathrm{new}}=P_{\mathrm{kf}}.Math.k_f,c_{f,\mathrm{new}}=P_{\mathrm{cf}}.Math.c_f\\ k_{r,\mathrm{new}}=P_{\mathrm{kr}}.Math.k_r,c_{r,\mathrm{new}}=P_{\mathrm{cr}}.Math.c_r\end{array},$

and the neural network obtains N.sub.out=[P.sub.kf, P.sub.kr, P.sub.cf, P.sub.cr] according to the time series signal N.sub.in=[z.sub.c, .sub.c, {umlaut over (z)}.sub.c, .sub.c, {dot over ()}.sub.c, {umlaut over ()}.sub.c, a, t, z.sub.tf, z.sub.tr].

[0049] As shown in FIG. 3, the neural network uses a recurrent neural network (RNN), including: a three-layer long short-term memory (LSTM) network and a fully connected layer, where: the hidden layer size is 128, and the fully connected layer calculates the suspension correction parameters and outputs them.

[0050] The loss function of the neural network is the mean square error (MSE), and the compared physical signals are the vertical acceleration and pitch rate of the vehicle body sprung mass, specifically:

[00026]$MSE=\frac{1}{n}.Math.{\left(N_{out}^{test}-N_{out}^{target}\right)}_t^2.$

[0051] To avoid the problems of gradient disappearance or gradient explosion during the training of the neural network, the Xavier initialization method is used to initialize the RNN network. The linear output layer uses the normal distribution initialization method to enhance training stability and accelerate network convergence.

[0052] Through specific actual experiments, the TR100A mining truck produced by TEREX is selected as the experimental platform. The test site is located in the Ruomaoshan mining area in Wuhu City, Anhui Province, China. In the experiment, the truck position information and all kinematic information were collected by the IFS2100 high-robustness, high-precision positioning and attitude setting system designed by DAISCH. All data were collected at 100 Hz, and the past 2-step state is input into the hybrid neural network model and physical model for calculation.

[0053] The sensor is installed on the cabin floor, and its x-axis and z-axis positive directions are the same as the vehicle coordinate system. However, the y-axis positive direction is opposite to the vehicle coordinate system, and the positive direction needs to be considered during calculation.

[0054] The experiment is conducted on a flat paved road, and the mining truck is in an unloaded state. For simplicity, the ground excitation is ignored compared to the tire size during transportation. The experiment included longitudinal acceleration and deceleration tests, and the specific tests are shown in Table 1.

TABLE-US-00001 TABLE 1 Mining Truck Test Status Test Number Status Velocity Part 1 Acceleration 0-20 km/h Part 2 20-30 km/h Part 3 Deceleration 10-0 km/h Part 4 20-0 km/h Part 5 30-0 km/h

[0055] The mining truck size parameters, stiffness characteristics, inertia information, and sensor installation position required during the calculation process are shown in Table 2, where x, y, and z represent the coordinates of the CoG relative to the IMU installation position. Calculate the suspension stiffness parameters, and the average fitting accuracy (1-NRMSE) under all working conditions is 84.76%, and the average FPE is 0.024. The fitting result has high accuracy and is used as the reference input for subsequent model and network calculations.

TABLE-US-00002 TABLE 2 Mining Truck Parameters Parameter Value (Unit) Parameter Value (Unit) a.sub.c 2.33 (m) m.sub.tf 8150.45 (kg) b.sub.c 2.24 (m) m.sub.tr 18423.55 (kg) x 3.30 (m) k.sub.tf 5.25 107 (N/m) y 1.24 (m) k.sub.tr 7.96 107 (N/m) z 0.45 (m) k.sub.f 1.87 106 (N/m) h 1.12 (m) k.sub.r 1.93 106 (N/m) m.sub.c 42046/(kg) c.sub.f 1.59 104 (Ns/m) I.sub.c 234333 (kg .Math. m.sup.2) c.sub.r 2.00 104 (Ns/m)

[0056] After clarifying the parameters and working conditions of the mining truck, use the above equations to train the network and establish the model. The specific comparison results are given below. Before calculation and network training, all collected data are aligned on the same timeline and converted according to the theory proposed in this invention.

[0057] As shown in FIG. 4, it is the original value of the longitudinal acceleration, vertical acceleration, and pitch angular velocity measured by the sensor. As shown in FIG. 5, it is the truck speed in each part of the test set.

[0058] During the network training process, distinguish the training set and test set for all working conditions in Table 1. To conduct reasonable training and comparison under different working conditions, use the first 40% and the last 40% of the data of all working conditions as the training set, and the remaining 20% of the data as the test set.

[0059] To better verify the accuracy of this invention, select the physical model, data-driven model, and hybrid model for comparison. Among them, the physical model is the classic half-vehicle model, and the data-driven model uses the LSTM network. All training parameters participating in the training of this invention are consistent to verify the applicability of each model. The design of mining trucks generally prevents the cargo from falling due to excessive vertical acceleration during transportation, thereby avoiding transportation capacity loss and road interference. In addition, obtaining the suspension movement speed is also beneficial for subsequent suspension control and tire dynamic load calculation. Therefore, this invention focuses on these two indicators and compares the calculation results of different models.

[0060] Select Part 1 and Part 5 of the working conditions for display to verify the model effect during acceleration and deceleration. As shown in FIGS. 6 and 7, they are the comparison results of the vertical acceleration of the center of gravity in the two working condition test groups. In the physical model, the calculation of the unsprung position needs to solve the differential equation during the calculation process. This process will introduce errors and oscillations, causing severe oscillations in the final state calculation at the next moment. However, this method has a simple calculation process and can reflect the movement law of the mining truck. The accuracy of the data-driven model is higher than that of the pure physical model, but it is limited by the limited samples and the large number of iterations (about 10-15 epochs, depending on the sample size), which limits its generalization ability. This invention can effectively reflect the real motion state and converge within a limited number of iterations (about 2-3 epochs). From the comparison results, it can also be seen that the hybrid model still has strong robustness under some harsh working conditions. Overall, due to the combination of the advantages of the physical model and the data-driven model, the hybrid model performs better.

[0061] In Part 4 of the working conditions, the comparison results of the vertical speed of the front and rear suspension sprung positions calculated by the vertical speed of the center of gravity (obtained by integrating the acceleration) and the pitch angular velocity are shown in FIGS. 8 and 9. It can be seen that after integral processing, the accuracy of the physical model solution results is significantly improved, and it can represent the movement law of the vehicle body. However, like the data-driven model, the error is still large under specific or harsh working conditions. The hybrid model accurately reflects the sprung speed change pattern in most cases, but in a very few working conditions, the accuracy is not as good as the data-driven model. The reason for this phenomenon may be that the hybrid model only uses two consecutive states, and due to sensor noise, the physical model may amplify the predicted value, resulting in more significant prediction errors at individual points. But overall, the calculation accuracy of the hybrid model is higher than that of the other two models, and its prediction results are relatively more reliable.

[0062] To give a qualitative comparison result of the three models, this invention uses the MSE of all test sets as a measure. The smaller the value, the higher the accuracy of the model. The comparison results of the MSE of the vertical acceleration of the center of gravity a.sub.c, the front suspension sprung speed v.sub.f, and the rear suspension sprung speed v.sub.r are shown in Table 3.

TABLE-US-00003 TABLE 3 Comparison of MSE in the experiment MSEvalue Physicsmodel Data-drivenmodel Hybridmodel a.sub.c(m/s2)2 13.69 0.18 0.13 v.sub.f(10 2 m/s)2 2.36 0.76 0.32 v.sub.r(10 2 m/s2)2 2.23 0.78 0.31

[0063] It can be seen that this invention has good performance in predicting different motion parameters. It is worth mentioning that the hybrid model method of this invention only relies on the vehicle body sensor signal and does not require unsprung information. This advantage enables direct prediction of the dynamic characteristics of the mining truck without the need to build a complex mining truck sensor measurement system.

[0064] In addition, hydro-pneumatic suspension is particularly common in mining trucks due to its excellent durability. However, due to its design, the modeling of the nonlinear stiffness and damping characteristics of hydro-pneumatic suspension is challenging and difficult to directly apply to the calculation of truck dynamic characteristics. This invention provides the suspension parameters of mining trucks, which are crucial for the structural design and performance prediction of mining trucks. As shown in FIGS. 10 and 11, they are the variation patterns of the front and rear suspension stiffness and damping coefficients.

[0065] When the mining truck accelerates, that is, in Part 1 and Part 2 of the working conditions, the rear suspension is compressed, and the front suspension is tensioned. When the hydro-pneumatic suspension is compressed, the internal gas pressure increases, and the stiffness increases, and vice versa. Therefore, as shown in FIG. 10, the stiffness of the rear suspension should be higher than that of the front suspension. When the mining truck decelerates, the stiffness of the front suspension is significantly more significant. In addition, when observing the change of suspension stiffness during the switching process of the mining truck's acceleration and deceleration, it is not difficult to find that the front suspension has obvious changes. This phenomenon is because the cross-section and internal gas mass of the front suspension are relatively small, so under the same volume and pressure change rate (isothermal assumption), its stiffness change is obvious. Since the rear suspension mainly bears the load-bearing function, its cross-section is larger, so the suspension stiffness change is smaller during the working condition switching. Usually, the front and rear suspension stiffness changes should be opposite.

[0066] From the perspective of the damping coefficient, as the speed of the mining truck increases, its vertical movement will also intensify. At this time, the speed of the liquid flowing through the throttle hole in the hydro-pneumatic suspension will also tend to accelerate. According to fluid mechanics theory, the damping will increase in this case, as shown in FIG. 11. However, the influencing factors of damping are very complex, and a relatively reasonable explanation is given here.

[0067] Compared with the existing technology, this method can only calculate the time-varying stiffness and damping values of the front and rear suspensions of the mining truck when it is running on unstructured roads through the IMU.

[0068] The above specific implementation can be partially adjusted by those skilled in the art without departing from the principles and purposes of the invention. The protection scope of the invention is subject to the claims and is not limited by the above specific implementation. All implementation schemes within its scope are bound by the invention.