Operational Modal Analysis for Vibration Control Following Moving Window Locality Preserving Projections for Linear Slow-Time-Varying Structures

Modal parameters can reflect the dynamic characteristics of the structure and can be used to control vibration. To identify the operational modal parameters of linear slow-time-varying structures only from non-stationary vibration response signals, a method based on moving window locality preserving projections (MWLPP) algorithm is proposed. Based on the theory of “time freeze”, the method selects a fixed length window and takes the displacement response signal in each window as a stationary random sequence. The locality preserving projections algorithm is used to identify the transient modal frequency and modal shape of the structure at this window. The low-dimensional embedding of the displacement response data set calculated by locality preserving projections (LPP) corresponds to the modal coordinate response matrix, and the transformation matrix corresponds to the modal shape matrix. The simulation results of the mass slow-time-varying three degree of freedom (DOF) and the density slow-time-varying cantilever beam show that the new method can effectively identify the modal shape and modal natural frequency of the linear slow-time-varying only from the non-stationary vibration response signal, and the performance is better than the moving window principal component analysis (MWPCA).


Introduction
Vibration is an inherent property of structures and bad vibrations will cause damage to structures [1]. Operational modal analysis (OMA) identifies structural modal parameters (modal natural frequency, modal shape, damping ratio) from the output displacement response signal of the structure [2,3], which can be used for damage identification [4], structural design [5], structural health monitoring of aircraft wings [6], and building performance assessment [7].
Manifold Learning [8] has become a research focus in information science since it was first proposed in 2000. Assuming that the data are uniformly sampled from a lowdimensional manifold in a high-dimensional Euclidean space, the purpose of manifold learning is to find the low-dimensional embedding in the high-dimensional space and find the corresponding embedded transformation matrix [9]. Based on this concept, the dimensionality reduction process of manifold learning needs to keep the data after dimensionality reduction satisfying the geometric constraint relation related to high-dimensional space manifold. At present, the more widely applied algorithms of manifold learning include principal component analysis (PCA) [10], isometric feature mapping (Isomap) [11], Laplacian eigenmaps (LE) [12], locally linear embedding (LLE) [13] and so on. Many researchers have applied manifold learning to operational modal analysis. Wang et al. used principal component analysis for operational modal parameters identification to solve the problem of false modal in other identification methods [14]. Bai et al. used locally linear embedding to identify operational modal parameters of complex three-dimensional continuous structures [15]. Wang et al. used isometric feature mapping to identify the operational modal parameters of three-dimensional structures [16]. The method mentioned above is suitable for linear time-invariant structures.
However, operational modal parameters of many structures are time-varying [17]. A linear structure whose system parameters (i.e., mass, stiffness, or damping) change over time is called a linear time-varying (LTV) structure. For example, during the flight launch of missile or rocket, the mass of the whole structure is gradually reduced due to the constant consumption of fuel [18]. With the activity of human flow, the structure characteristics of large human gathering places change, such as stadiums [19]. Therefore, time-varying systems are different from constant system. Time-varying structures need to identify time-varying parameters, which can monitor the state of the structures. At present, the time-domain and frequency domain methods are used to identify the operational modal parameters of time-varying structures [20,21]. Zhou et al. analyzed and summarized the advantages and disadvantages of these methods in detail [22]. The time-varying structures can be divided into fast time-varying structures and slow-time-varying structures [23] and Ramnath pointed out the system whose variation in coefficient is much less than variation in solution is called slow time-varying system [24]. For example, when a train passes a bridge quickly, the train works together with the bridge to form a unified dynamic system. The mass distribution and stiffness distribution of the system change rapidly with time, which forms a fast time-varying structural dynamics problem [25]. The mechanical arms in industrial manufacture can be regarded as slow time-varying structures [26]. The vibration response signals of linear slow-time-varying structures cannot be obtained completely at one time, but need to be obtained by continuous sampling over time. Therefore, the moving window method based on the theory of "short time-invariant" can be well applied to the identification of operational modal parameters of linear slow-time-varying structures. At present, the moving window method has been applied to some algorithms to identify the operational modal parameters of linear slow-time-varying structures. Huang et al. proposed moving window EASI algorithm identify the operational modal parameters of the linear slow-time-varying system [27]. Wang et al. proposed a moving window second order blind identification method for identifying operational modal parameters of linear slow-time-varying structures, and the performance of this method is better than the moving window independent component analysis [28]. Guan et al. combined moving window with principal component analysis to effectively identify the operational modal parameters of the slow-time-varying system [29,30]. Huang et al. propose an operational modal analysis (OMA) method that uses eigenvector recursive PCA with a forgetting factor to identify the transient natural frequencies and transient modal shapes [31].
In addition, there are many time-domain methods for operational modal analysis. Some methods involve choosing a mathematical model to idealize the structural dynamic responses, including autoregressive moving average (ARMA) [32] and autoregressive (AR) model updating [33]. Based on singular value decomposition and QR factorization, Barros-Rodriguez et al. proposed a new method and applied it to the analysis of F-18 flutter flight test data. The method is capable of identifying the frequency and damping of the critical aircraft modes, those responsible for the flutter phenomenon [34]. Chen et al. developed a novel method for moving force identification (MFI) called preconditioned least square QR factorization (PLSQR) method which seeks to reduce the impact of identification errors caused by unknown noise [35]. Some methods, such as time-frequency analysis method, can use the response signal of a vibration sensor to identify multi-modal parameters, but generally can only identify the natural frequency and damping ratio of multi-modal, but cannot identify modal shapes. Locality preserving projections (LPP) is a linear dimensionality reduction method [36], which makes it fast and suitable for practical application. In addition, LPP has several nonlinear technology data representation characteristics, which makes it more accurate in preserving data characteristics. Therefore, a method based on moving window locality preserving projections (MWLPP) is proposed to identify the operational modal parameters of linear slow-time-varying structures. LPP first calculates the transformation matrix of mapping the high-dimensional data to the low-dimensional space, and then computes the low-dimensional embedding. The primary contributions of the article can be summarized as follows: (1) An operational modal parameter identification method based on LPP is proposed.
The main idea is to find out the one-to-one correspondence between the coordinate response matrix and the low-dimensional embedded data, and the one-to-one correspondence between the modal shape matrix and the transformation matrix. The operational modal parameter identification problem can be transformed into the manifold dimension reduction problem of structural vibration response data. (2) The LPP algorithm and moving window method are combined to identify the operational modal parameters of linear slow-time-varying structures. (3) By comparing the operational modal parameter identification method based on moving window principal component analysis method, MWLPP has higher accuracy and effectively reduces the modal missing. (4) To verify the identification ability of operational modal parameters of linear slow time-varying structures based on MWLPP, a mass slow-time-varying 3-DOF (degree of freedom) structure and a density slow-time-varying cantilever beam structure were designed. The MWPCA method was used to identify the operational modal parameters of non-stationary vibration response simulation data for comparison.
The remainder of this paper is organized as follows. In Section 2, The LPP algorithm is introduced to identify the operational modal parameters of linear time-invariant structures. The moving window combined with LPP algorithm is introduced to identify the operational modal parameters of linear slow-time-varying structures in Section 3. Section 4 presents the simulation verification results. Finally, we make a conclusion in Section 5.

Problem of OMA of Linear Time-Invariant Structures
According to the dynamic theory of the structure, the dynamic equation of n degrees of freedom (DOF) linear time-invariant vibration structure in the physical coordinate system is: M ..
where M ∈ R n×n , C ∈ R n×n and K ∈ R n×n are the mass matrix, damping matrix, and stiffness matrix of the structure. T is number of sampled data points. X(t) ∈ R n×T , . X(t) ∈ R n×T and .. X(t) ∈ R n×T are the time-domain sampling matrix of the displacement response signal, velocity response signal and acceleration response signal of the structure. F(t) ∈ R n×T is the time-domain sampling matrix of the external excitation.
T ∈ R n×T of n DOF small damping structure in modal coordinates is: where d is modal truncation and ranges from 1 to n. Φ = [ ∈ R n×d is the modal shape matrix constituted by d order modal shape → φ i ∈ R n×1 (i = 1, 2, · · · , d).
T ∈ R d×T is a modal response matrix constituted Appl. Sci. 2021, 11, 791 4 of 26 by d order modal response → q i (t) ∈ R 1×T (i = 1, 2, · · · , d). The main idea of OMA is to identify the modal shape matrix Φ ∈ R n×d and modal response matrix Q(t) ∈ R d×T of the structure only from the vibration displacement response signal X(t) ∈ R n×T . Finally, the natural frequency f and damping ratio ξ are identified by using the single degree of freedom (DOF) technique [37] from the modal response matrix Q(t) ∈ R d×T . For example, the natural frequency f was identified by fast Fourier transform (FFT), and the damping ratio ξ was identified by the random decrement technique (RDT) and the Hilbert transform (HT) from the modal response matrix Q(t) ∈ R d×T . The idea was proposed in [38,39].
When the value of each natural frequency is different, The main modal shape → φ i satisfies normalized orthogonality and the modal response → q i (t) of each other is independent.

OMA of Linear Time-Invariant Structures Based on LPP
LPP is a dimensionality reduction algorithm that preserves the geometric relations within the data set. A dataset of T real valued vectors in R d×d space, which has the same geometric properties as X(t). There is a transformation matrix A that makes . LPP calculates the d-dimensional embedding process of n-dimensional data as follows.
(1) Construct an adjacent-graph G of T real valued vectors in R n×n space: K nearest neighbor points of the node → x i (t)(i = 1, 2, · · · i · · · j · · · , T) are obtained by using the K-neighbor algorithm. If (2) Calculate the weight of the edge: let the matrix W ∈ R T×T represents the weight matrix, and w ij is the weight of the edge ( . The weights of the connected edges are calculated by heat kernel If the two nodes are not connected, the weight is 0. (3) Calculate transformation matrix A: the eigenvectors and eigenvalues of the generalized eigenvector problem are calculated as follows.
where A ∈ R n×d is the transformation matrix, A T X(t) is the d-dimensional embedding S(t) ∈ R d×T of X(t) ∈ R n×T . Therefore, X(t) ∈ R n×T has the following form: Comparing Equations (2) and (10), we can make the conclusion that the low-dimensional embedded S(t) ∈ R d×T corresponds to the modal response Q(t) ∈ R d×T , and the transformation matrix A ∈ R n×d corresponds to the modal shape Φ ∈ R n×d . Figure 1 illustrates OMA of time-invariant structure based on LPP algorithm. Compared with State matrix eigen-decomposition, LPP algorithm builds a graph incorporating neighborhood information of the data set. Using the notion of the Laplacian of the graph, LPP algorithm computes a transformation matrix which maps the data points to a subspace. This linear transformation optimally preserves local neighborhood information in a certain sense [36]. This feature helps retain information about the response signal.

Problem of OMA of Linear Time-Varying Structure
The operational modal parameters of time-varying structures vary with time. According to the dynamic theory of the structure, the motion equation of a linear time-varying structure of n degrees of freedom in the physical coordinate system in  are the mass matrix, damping matrix, and stiffness matrix of the time-varying structure.
( )  Compared with State matrix eigen-decomposition, LPP algorithm builds a graph incorporating neighborhood information of the data set. Using the notion of the Laplacian of the graph, LPP algorithm computes a transformation matrix which maps the data points to a subspace. This linear transformation optimally preserves local neighborhood information in a certain sense [36]. This feature helps retain information about the response signal. ..
where M(t) ∈ R n×n , C(t) ∈ R n×n and K(t) ∈ R n×n are the mass matrix, damping matrix, and stiffness matrix of the time-varying structure. X(t) ∈ R n×T , . X(t) ∈ R n×T and .. X(t) ∈ R n×T are the time-domain sampling matrix of the displacement response signal, velocity response signal and acceleration response signal of the time-varying structure. F(t) ∈ R n×T is the time-domain sampling matrix of the external excitation. According to the "time-freezing" theory [40], the mass, damping ratio and stiffness of a linear time-varying structure can be regarded as time-invariant for a short time τ ∈ [t begin , t end ]. Therefore, Equation (11) can be expressed as a set N consisting of a finite number of linear time-invariant structures N (τ)(τ ∈ [t begin , t end ]) within the complete time ..
The vibration response data of linear time-varying structures can be divided into a finite number of time-invariant parts by selecting moving window with fixed length T ∈ R d×L is the modal response matrix formed by the modal response vector → q i j (t) ∈ R 1×L (j = 1, 2, · · · , d) of the structure in i-th (τ ∈ t begin , t end ) window. When the order of each natural frequency is different, The main modal shape → φ i j satisfies normalized orthogonality and the modal response

OMA of Linear Time-Varying Structure Based on MWLPP
Based on the "time freeze theory" and fixed window length moving windows (MW), the non-stationary signals in each window are regarded as stationary signals. The vibration response signals in the window are identified by using the linear time-invariant OMA method. After OMA of the window is completed, the window moves to the right to delete some old data and add new data to form the vibration response data of the next window. A window corresponds to a moment, so as to identify the operational modal parameters of that moment. Finally, the modal parameter identification results of all windows (moments) are connected to form a continuous result. Moving windows has fixed window length L and moving step λ. Figure 2 shows the process of moving the moving window.  Next moving window data of fixed window length L

Begin
Vibration response data of time-varying structures measured by n sensors ， Set the length of the moving window L ，step ， Moving window number i = 1.
The response data in this window with LPP algorithm dimension reduction Window i : displacement response data of window length L.
Move to the next window ( ) Therefore, moving windows are used to track the time-varying properties of the structure. Comparing Equations (2) and (13), it can be seen that the decomposition of vibration response signals in the fixed window length window is consistent with the decomposition of vibration response signals in the linear time-invariant structure. The low- Figure 3 illustrates the moving window principle processing with time. Figure 4 illustrates OMA of linear time-varying structure based on MWLPP algorithm.  Next moving window data of fixed window length L

Begin
Vibration response data of time-varying structures measured by n sensors The response data in this window with LPP algorithm dimension reduction Window i : displacement response data of window length L.
Move to the next window ( ) Natural frequency

The Applicable Scope of the Method
In this paper, the application scope of MWLPP for OMA of linear slow-time-varying structure is, (1) The method is only suitable for linear slow-time-varying structures with small damping. If the damping is too high, the modals will be complex. Reference [16] show that the damping ratio reaches 10%, and the operational modal parameters can also be identified by manifold learning. However, the lower the damping ratio, the better the identification of modal parameters effect. Only for linear slow-time-varying structures, based on the "time-freezing" theory, Equation (12) can be expressed as Equation (13). (2) The number of vibration response sensors n should be greater than or equal to the d-order modal identified by the method. According to Equation (2), the displacement response signal in modal coordinates is approximately represented by d-order modal (d ≤ n). In addition, MWLPP can identify modal natural frequency, modal shape, and damping ratio of one mode with only 1 sensor. (3) The excitation vector to the structure should be approximately stationary Gaussian white noise. (4) The method can identify time-varying transient modal shapes, modal frequencies, and modal damping ratios from non-stationary vibration response signals. However, for linear slow-time-varying structures, mass reduction or motion will generate additional damping [41,42]. Therefore, the damping ratio identified by MWLPP cannot be directly compared with real values.

Simulation Verification
To verify the identification ability of operational modal parameters of linear slow timevarying structures based on MWLPP, a mass slow-time-varying 3-DOF structure and a density slow-time-varying cantilever beam structure were designed. The moving window principal component analysis (MWPCA) method was used to identify the operational modal parameters of non-stationary vibration response simulation data for comparison.

An Introduction to Simulation Systems and Data Sets
A mass slow-time-varying 3-DOF structure and a density slow-time-varying cantilever beam are designed in Matlab/Simulink. More simulation details are available in [28,30]. In the simulation, all the natural frequencies in finite element analysis (FEA) are real natural frequency, and the modal shapes in FEA are real modal shapes. In addition, 10% white Gaussian noise is added to the vibration displacement response signal of the slow-timevarying 3-DOF structure, and the modal parameters of the structure are identified.

A Mass Slow-Time-Varying 3-DOF Structure
The mass slow-time-varying 3-DOF structure is designed in Matlab/Simulink. The model of the structure is shown in Figure 5. The dynamic equation of the structure is shown in Equation (16). ..
..  1  1  12  2  1  12  2  1   2  2  2  2  3  3  2  2  2  3  3  2   3  3  3  3  3  3  3 The structural parameters are as follows: The damping ratio of the three objects is   Table  1, the highest natural frequency of the third order mode is relatively small, much less than 20 Hz. In general, the higher the sampling frequency is set, the more accurate the modal parameter identification results will be. However, at the same time, the vibration response semaphore and calculation amount of sampling will also be larger. Therefore, according to Nyquist theorem, the sampling frequency is set to 40 Hz, and the sampling interval of simulation is 0.025 s. The displacement response signal data set ( ) of the structure is obtained by Runge-Kutta algorithm in Matlab/Simulink module. The data set The structural parameters are as follows: The damping ratio of the three objects is c 1 = c 2 = c 3 = 0.01(Ns/m); the stiffness is k 1 = k 2 = k 3 = 1000(N/m); the initial displacement is 0. The quality is m 2 = m 3 = 1kg. The mass m 1 are time-varying, and the change rule is shown in the following Equation (17).
F 1 (t) is white Gaussian noise, and is applied to object m 1 . The simulation consists of two parts. In t ≤ 50 s, the mass of object m 1 is constant value. In 50 s < t ≤ 2000 s, the mass of object m 1 changes with time, and the system becomes a linear time-varying structure. The 3-DOF simulation structure has only third order mode at most. From Table  1, the highest natural frequency of the third order mode is relatively small, much less than 20 Hz. In general, the higher the sampling frequency is set, the more accurate the modal parameter identification results will be. However, at the same time, the vibration response semaphore and calculation amount of sampling will also be larger. Therefore, according to Nyquist theorem, the sampling frequency is set to 40 Hz, and the sampling interval of simulation is 0.025 s. The displacement response signal data set X(t) ∈ R n×T of the structure is obtained by Runge-Kutta algorithm in Matlab/Simulink module. The data set X(t) ∈ R n×T of three degrees of freedom structure is shown in Figure 6, and the sampling duration is 2000 s. Therefore, the final data set is X(t) ∈ R 3 * 80,000 . structure. The 3-DOF simulation structure has only third order mode at most. From Tab 1, the highest natural frequency of the third order mode is relatively small, much less th 20 Hz. In general, the higher the sampling frequency is set, the more accurate the mod parameter identification results will be. However, at the same time, the vibration respon semaphore and calculation amount of sampling will also be larger. Therefore, accordi to Nyquist theorem, the sampling frequency is set to 40 Hz, and the sampling interval simulation is 0.025 s. The displacement response signal data set ( ) of the stru ture is obtained by Runge-Kutta algorithm in Matlab/Simulink module. The data s ( ) of three degrees of freedom structure is shown in Figure 6, and the sampli duration is 2000 s. Therefore, the final data set is 3*80,000  The real natural frequency of the structure is shown in Figure 7. The real natural frequencies at t = 50.025 s, t = 1200 s and t = 1974.375 s are shown in Table 1   The real natural frequency of the structure is shown in Figure 7. The real natural frequencies at t = 50.025 s, t = 1200 s and t = 1974.375 s are shown in Table 1 below. The natural frequencies identified by MWLPP at t = 50.025 s, t = 1200 s and t = 1974.375 s are shown in Figure 8.

A Density Slow-Time-Varying Cantilever Beam Structure
The density slow-time-varying cantilever beam structure is more complicated, as shown in Figure 9. One end of the cantilever beam is a fixed support, the other end is a free end, and many actual engineering structure can be simplified as a cantilever beam. Therefore, the modal parameters identified in the cantilever beam data set are also valuable.

A Density Slow-Time-Varying Cantilever Beam Structure
The density slow-time-varying cantilever beam structure is more complica shown in Figure 9. One end of the cantilever beam is a fixed support, the other e free end, and many actual engineering structure can be simplified as a cantilever Therefore, the modal parameters identified in the cantilever beam data set are also ble.
Assuming that the shear deformation of the cantilever beam is ignored, the can beam structure is evenly divided into 40 elements using the finite element method. same time, the axial displacement of the beam is not considered, but only the vertic placement and rotation angle.  ... Assuming that the shear deformation of the cantilever beam is ignored, the cantilever beam structure is evenly divided into 40 elements using the finite element method. At the same time, the axial displacement of the beam is not considered, but only the vertical displacement and rotation angle.
In the finite element model of cantilever beam, the cantilever beam is divided into 40 units and the its equation of motion such as Equation (1) is established. In addition, mass matrix M e , damping matrix C e and stiffness matrix K e of each window can be expressed as: where L represents window length, β M and β K are proportional coefficients. Then the mass matrix M e , damping matrix C e and stiffness matrix K e are assembled into the total mass matrix MM, damping matrix CC and stiffness matrix KK of the system: Therefore, the modal parameter can be calculated by finite element method, and the mode frequency, mode shape and modal damping ratio can be expressed as: where MM r , and CC r are the r-th order of modal mass matrix, modal stiffness matrix and modal damping matrix respectively. f r , → φ r and ξ r are the r-th order of natural frequency, modal shape, and damping ratio. N is the number of sensor.
To achieve the slow-time-varying condition, the density of the cantilever beam is changed with time, as shown in Equation (25).
where rate of density change is 0.08, the total time is 4 s. Because the beam structure is a continuum rather than a multi-degree of freedom system, there are theoretically infinite modes. Even after discretization of a 40 degrees of freedom structure, there are theoretically 40 modes, the value of 700 Hz would not be sufficient then. In addition, the modal of the structure is generally concerned with the low frequency, and the analysis is only within 2000 Hz at most. Moreover, high-frequency noise is not considered in this article, so the sampling frequency of the structure is set at 10,000 Hz. The initial condition of the system is 0, and the density of the cantilever beam system remains unchanged in the first 0.5 s, in order to prevent the vibration system from being affected by random excitation in the initial stage of vibration. We took the data after 0.5 s for the simulation After 0.5 s, white noise excitation acts on the free end of the cantilever beam, and Newmark − β method is adopted to collect the vibration response signal of displacement of each node on the beam. The time step of Newmark − β method is 1/10, 000s, the parameter of β is 0.5, and the parameter of γ is 0.25, and the damping coefficient is β M = 4 × 10 −4 , β K = 1 × 10 −7 . The numerical simulation and algorithm writing of this paper are completed by using Matlab language and software. Finally, the displacement response X(t) ∈ R 40×40,000 is obtained, and the 1st, 20th and 40th elements are plot in Figure 10. As long as the stability of the algorithm is guaranteed, the parameters could be set according to the requirements, unless the algorithm is very sensitive to noise numerically and generates bias. Of course, the higher sampling frequency, the Newmark − β method for solving the vibration response signal of simulation data set are more accurate, the effect of modal identification are more accurate. The real natural frequency variation of the first three modes of the structure is shown in Figure 11. The real natural frequencies of t = 0.5 s, t = 2 s and t = 3.795 s are shown in Table 2.  The real natural frequency variation of the first three modes of the structure is shown in Figure 11. The real natural frequencies of t = 0.5 s, t = 2 s and t = 3.795 s are shown in Table 2. The real natural frequency variation of the first three modes of the structure is sho in Figure 11. The real natural frequencies of t = 0.5 s, t = 2 s and t = 3.795 s are shown Table 2.   Figure 11. The real natural frequency of the density slow-time-varying cantilever beam structure.

The Evaluation Indexes
In the simulation, the average model assurance criterion (MAC i−avg ) will be used to identify the average accuracy of modal shapes. The specific definition of MAC i−avg is, where → φ ij is the modal shape of j-th window of the i-th order identified by algorithm; → ψ ij is real modal shape of j-th window of the i-th order identified; m is the number of windows. It can be seen from Equation (26) that the value range of MAC i−avg is 0 ≤ MAC i−avg ≤ 1. The closer the MAC i−avg value is to 1, the higher the average identification accuracy of the modal shape of i-th order.
Average error rate (δ i−avg ) will be used to evaluate the natural frequency identified by algorithm. The specific definition of δ i−avg is, where f ij is the natural frequency of j-th window of the i-th order identified by algorithm; f ij−real is real natural frequency of j-th window of the i-th order; m is the number of windows. It can be seen from Equation (27) that the value range of δ i−avg is δ i−avg ≥ 0. The closer the δ i−avg value is to 0, the higher the average identification accuracy of the natural frequency of i-th order. It is important to choose an appropriate limited memory length L of moving window. Frequency resolution f and average frequency variation of the i-th modal in a window f L (i) was used to select window length. Window length L is proportional to sample frequency f s and frequency resolution f , f cannot be too small, because it cannot reflect the change of frequency, and f cannot be too large, because the change of f L (i) cannot be identified, at the same time, f L (i) cannot be larger too much than f , otherwise, LTV structure cannot be regarded a LTI structure in a window. Therefore, different structures may choose different window lengths. The idea was proposed in [28,30].
where variables f s , f end (i), f begin (i), t end , t begin are sample frequency, end-frequency of the i-th modal, begin-frequency of the i-th modal, end-time, begin-time of the whole data.
The model assurance criterion (MAC) will be used to identify the accuracy of modal shapes. The specific definition of MAC is, where → φ is the modal shape identified by algorithm; → ψ is real modal shape. The value range of MAC is 0 ≤ MAC ≤ 1. The closer the MAC value is to 1, the higher the identification accuracy of the modal shape.

The Parameter Settings
In the 3-DOF structure, the window length L is 1024, sample frequency f s is 40 Hz and sample time is 2000 s. Frequency resolution is f = 0.039 Hz, the first average frequency variation is f L (1) = 9.25 × 10 −4 Hz, the second average frequency variation is f L (2) = 0.0130 Hz, the third average frequency variation is f L (3) = 0.0426 Hz. Therefore, the length L = 1024 satisfies average frequency variation f L (i) cannot be lager too much than frequency resolution f .
In the cantilever beam structure, the window length L is 2048, sample frequency f s is 10,000 Hz and sample time is 4 s. Frequency resolution is f = 4.88 Hz, the first average frequency variation is f L (1) = 0.15 Hz, the second average frequency variation is f L (2) = 1.09 Hz, the third average frequency variation is f L (3) = 3.06 Hz. Therefore, the length L = 2048 satisfies average frequency variation f L (i) cannot be lager too much than frequency resolution f . The parameter K in LPP algorithm is 40, and the lowdimensional embedded dimension is 3.

A mass Slow-Time-Varying 3-DOF Structure
In the linear slow-time-varying structure, the operational modal parameters of the system change at any time. The modal shapes at all times are difficult to describe. Therefore, in the mass slow-time-varying 3-DOF structure, we select four moments including t = 100 s, t = 650 s, t = 1500 s and t = 1974.375 s. Table 3 shows the MAC value at the instantaneous moment. Figure 12 shows the modal shape of the instantaneous moment. Figure 13 shows the time-varying natural frequencies identified by the MWLPP. Figure 14 shows the timevarying MAC values of modal shapes identified by the MWLPP. Table 4 shows the average error δ i−avg of MWLPP and MWPCA in identifying natural frequencies. Figure 15 shows the change with time of damping ratio of slow-time-varying 3-DOF structure identified by MWLPP. In the identification process of operational modal parameters of structures, modal parameters of some windows (moments) are not identified, i.e., modal parameters are missing. Those windows that modal parameters are not identified are described as un-identified windows. Table 5 shows the number of un-identified windows between MWLPP and MWPCA. Table 6 shows the MAC i−avg of modal shapes identified by MWLPP and MWPCA.  structures, modal parameters of some windows (moments) are not identified, i.e., modal parameters are missing. Those windows that modal parameters are not identified are described as un-identified windows. Table 5 shows the number of un-identified windows between MWLPP and MWPCA. Table 6 shows the MAC i avg − of modal shapes identified by MWLPP and MWPCA.

A Density Slow-Time-Varying Cantilever Beam Structure
In the density slow-time-varying cantilever beam structure, we select four moments including t = 0.75 s, t = 1.75 s, t = 2.75 s and t = 3.75 s. Table 7 shows the MAC value at the instantaneous moment. Figure 16 shows the modal shapes of the instantaneous moment. Figure 17 shows the time-varying natural frequencies identified by the MWLPP. Figure  18 shows the time-varying MAC values of modal shapes identified by the MWLPP. Table  8 shows the average error i avg − δ of MWLPP and MWPCA in identifying natural frequen-   In the density slow-time-varying cantilever beam structure, we select four moments including t = 0.75 s, t = 1.75 s, t = 2.75 s and t = 3.75 s. Table 7 shows the MAC value at the  In the density slow-time-varying cantilever beam structure, we select four moments including t = 0.75 s, t = 1.75 s, t = 2.75 s and t = 3.75 s. Table 7 shows the MAC value at the instantaneous moment. Figure 16 shows the modal shapes of the instantaneous moment. Figure 17 shows the time-varying natural frequencies identified by the MWLPP. Figure 18 shows the time-varying MAC values of modal shapes identified by the MWLPP. Table 8 shows the average error δ i−avg of MWLPP and MWPCA in identifying natural frequencies. Table 9 shows the number of un-identified windows between MWLPP and MWPCA. Table 10 shows the MAC i−avg of modal shapes identified by MWLPP and MW-PCA.

A Mass Slow-Time-Varying 3-DOF Structure with 10% White Gaussian Noise
Ten percent white Gaussian noise is added to the vibration displacement response signal of the slow-time-varying 3-DOF structure, and the modal parameters of the structure are identified.
In the mass slow-time-varying 3-DOF structure with 10% white Gaussian noise, we also select four moments including t = 100 s, t = 650 s, t = 1500 s and t = 1974.375 s. Table 11 shows the MAC value at the instantaneous moment. Figure 19 shows the modal shape of the instantaneous moment. Figure 20 shows the time-varying natural frequencies identified by the MWLPP. Figure 21 shows the time-varying MAC values of modal shapes identified by the MWLPP. Table 12 shows the average error δ i−avg of MWLPP in identifying natural frequencies. In the identification process of operational modal parameters of structures, modal parameters of some windows (moments) are not identified, i.e., modal parameters are missing. Those windows that modal parameters are not identified are described as un-identified windows. Table 13 shows the number of un-identified windows by MWLPP. Table 14 shows the MAC i−avg of modal shapes identified by MWLPP.

Analysis of Simulation Results
(  Tables 5 and 9, the total number of un-identified windows in MWLPP method is lower than that in MWPCA whether the structure is a linear mass slow-timevarying 3-DOF structure or a complex linear density slow-time-varying cantilever beam structure. Comparing Tables 6 and 10, the MAC value of modal shape identified by MWLPP is higher than that MWPCA. Therefore, the MWLPP method is superior to the MWPCA method in identifying the operational modal parameters of linear slow time-varying structures. (5) From Figure 17 and Table 8, the first natural frequency seems contain fluctuations because the ordinate scale of the first mode natural frequency is small. According to Equation (29), in the cantilever beam structure, the first average frequency variation is f L (1) = 0.15 Hz, the second average frequency variation is f L (2) = 1.09 Hz, the third average frequency variation is f L (3) = 3.06 Hz. The smaller the average frequency variation is, the smaller the variation range of the natural frequency value is. Due to the small variation range of the first order natural frequency value, the fluctuation phenomenon of the first order natural frequency is obvious in Figure 15 [28].
In addition, according to Equation (27), the denominator of the first order is much smaller than that of the second and the third order when calculating the average error rate δ i−avg of natural frequency because the natural frequency value of the first order is smaller than that of the second and the third order. Therefore, in Table 8, the average error rate of natural frequency calculated in the first order is greater than that in the second and third order. In fact, the result of modal shape identification shows that the first order identification is indeed the best. (6) Figure 15 shows that the damping ratio identified by MWLPP algorithm fluctuates greatly because theoretical analysis and numerical simulations indicate that a decreasing or moving mass and density will generate additional damping in the LTV structures [41,42]. In addition, compared with theoretical values, the identification accuracy of damping ratios has a certain error. It is generally predictable because the identification of the damping ratio itself is a difficult problem in the field of structural dynamics, and easily affected by the adopted identification algorithm. Therefore, the time-varying transient mode damping ratio identified by MWLPP is not suitable to compare with the mode damping ratio calculated by finite element methods. time-domain method used in this article cannot use the average technique in the frequency domain to identify the time-domain modes from the non-stationary random vibration response signals. Therefore, the algorithm is unstable, and some moments cannot identify the modal parameters, or the identified modal natural frequency and modal shape are not good. In addition, references [16,43] indicate that the change of damping will affect the performance of the algorithm in identifying modal parameters. From the Figure 15, damping ratio of time-varying structure constantly changing, which leads to low precision of modal parameter identification in some time. In addition, References [28,30] indicate that the window length L will also affect the identification accuracy of modal parameters. From Equation (28), the window length L is proportional to the sample frequency f s and frequency resolution f . Frequency resolution f cannot be too small, because it cannot reflect the change of frequency, and cannot be too large, because the average frequency variation f L (i) cannot be identified. At the same time, f L (i) cannot be larger too much than f , otherwise, time-varying structure cannot be regarded a time-invariant structure in a window. (8) According to the noise simulation results of 3-DOF structure, the MWLPP algorithm turns out to be robust to noise. This is probably LPP algorithm computes a transformation matrix which maps the data points to a subspace. This linear transformation optimally preserves local neighborhood information in a certain sense [36]. This feature helps retain information about the response signal.

Conclusions
In this paper, a new method based on MWLPP for identifying operational modal parameters of linear slow-time-varying structures is presented. The low-dimensional embedding calculated by LPP algorithm corresponds to the modal coordinate response matrix, the transformation matrix corresponds to the modal shape matrix, and the operational modal parameters of the linear slow-time-varying structure are identified by the moving window method. Compared with the moving window principal component analysis method, MWLPP has higher accuracy and less modal missing in identifying the operational modal parameters of linear slow-time-varying structures.
However, fixed window length is an important parameter. How to determine and change the fixed length of moving window adaptively by non-stationary vibration response signal has not been solved completely. It is of great significance to apply this method to practical engineering structures. In addition, a further study is to establish a practical experiment to verify the effectiveness of the method.

Conflicts of Interest:
The authors declare that there is no conflict of interest regarding the publication of this paper and regarding the funding that they have received.