Dynamics of the small-span railway bridge under moving loads

. The paper presents the results of dynamic analysis of the small-span railway bridge, subjected to an action of moving trains. Numerical simulations were performed using three different load models: series of moving forces, series of moving single-mass and double-mass oscillators. The parameters of the vehicle were taken from the existing EN57 train. The parameters of the bridge were taken from the existing steel span of 10,24 m long. In both cases, the dynamic parameters were identified based on free-response measurements using modal identification techniques. Vibrations of the midpoint of the bridge as well as the mass of the oscillator have been analyzed. Numerical results obtained for individual load models were compared with the results of in-situ tests performed under operating conditions.


Introduction
In 1829 the George Stephenson's 'Rocket' locomotive was presented. A one year later, in 1830, the first steamonly railway line was opened between Liverpool and Manchester. These historical events gave rise to the railway era in the world. In 1847 the Dee Bridge railway disaster took place with five fatalities. This tragic accident marked the beginning of the practical considerations of the bridge-vehicle interaction dynamics.
A moving vehicle induces vibrations of a bridge span, a bridge in turn affects the vibrations of a vehicle. Thus, we have a complex, mutually coupled dynamic system whose exact analysis is very complicated. The first attempts to describe the problem were presented by Krylov and Timoshenko [1,2]. They considered a single force model moving on the Euler-Bernoulli beam with a constant speed. Later on Saller and Jeffcott took into account a mass of the vehicle as well [3,4]. Many later studies refer to the above mentioned, basic dynamic models. The influence of the various parameters (velocity, mass, stiffness, damping) on the dynamic response of the bridge was investigated, e.g. [5][6][7].
A rapid development of railway transport took place in the second half of 20th century. The first high-speed train was introduced to traffic on Shinkansen railway in Japan in 1964. It shuttled between Tokyo and Osaka with the speed of 210 kph. Since then constant tendency to increase a traffic speed is observed. Currently, the conventional vehicles can move with the speed of up to 350 kph on many European and non-European railways. The higher train speed, the greater dynamic impact to a bridge. For this reason, bridge-vehicle interaction has become an important part of design and research work, particularly for high-speed railways. The simplest sprung vehicles were modelled in the form of single-mass, elastic or viscoelastic oscillators [8][9][10][11]. The spring force variation was assumed as the measure of the dynamic coupling between a bridge and a vehicle.
A typical train consists of the sequence of locomotives and carriages which in fact causes the cyclic excitations. It was simulated by developing the single load models into the series of loads [12][13][14]. Load spacing was corresponded to the distance of axes of bogies or wheelsets. A regular load spacing and a constant speed can cause the intrusive resonant effects [15,16].
The main requirements for high-speed railways refer to the safety, stability and riding comfort. The comfort of passengers depends on the car body accelerations. From the safety point of view the important issue is the relation between vertical and horizontal force in the wheel-rail contact area. The multi-body mechanical models of vehicles have been developed to accurately analyse the problem [17][18][19][20]. The influence of the railway track was also studied.
The advanced FEM software and numerical solution methods allow for the detailed modelling of structures. But the more complex model, the more parameters must be determined. For this reason, in the case of technical requirements, the dynamic effects of loads are simplified very often [21][22][23][24].
The paper deals with the dynamic analysis of the small-span railway bridge under moving train. Three different models of railway vehicle were investigated and compared: series of moving forces, series of moving single-mass and double-mass oscillators. The first model is commonly used in the design practises. The successive two models are simplified sprung mechanical systems, which allow for the bridge-vehicle interaction consideration. The FE model was developed according to the existing bridge span and the existing railway vehicle. The model was validated on the basis of modal parameters identification results (frequencies and damping ratios). Eigensystem Realization Algorithm (ERA) was used. Numerical results of the train passing the bridge were finally referred to the field measurements performed under operating conditions.

The Bridge
The bridge is situated along the E65 railway line Warsaw-Gdynia (km 321,153) and crosses Radunia River in Gdańsk (Fig. 1). The E65 railway line is one of the major transport route in Poland. It belongs to the IV European Transport Corridor and characterises of high traffic density. The passenger train can travel up to 160 kph.  The bridge deck consists of I-shaped crossbeams spaced every 2,56 m and longitudinal beams spaced in 1,80 m distance. The railway track consists of the rails fastened to the wood sleepers.
This type of structure was widely used in the case of the small-span steel bridges with an open deck. In 2014, within the modernization of E65 railway line, the steel spans were replaced by the reinforced concrete plate girders. Just before the reconstruction the field measurements of the steel bridge were performed. The response signals were collected for modal parameters identification and the FE model validation.

FE model
The numerical model was developed using the SOFiSTiK software (Fig. 3). Beam elements with 6dof's of each node were used. Additional masses related to the existing equipment (sidewalk, railway track) were also applied. The model was validated according to in-situ measurements and modal identification results of the actual span. The criterion was the equality between measured and theoretical modes.
The open railway track was included in the substitute manner. Only the rails were considered as the structural finite elements. The connection between rails and longitudinal beams was assumed as rigid. Spring elements of the "infinite" stiffness were used between corresponding nodes (Fig. 3).

Field measurements
The field measurements of the bridge were performed in order to collect: free-decay signals for modal parameters identification and model validation, force vibration signals for comparison with the results of the bridge-vehicle numerical simulations. Displacement and acceleration data of the external span were collected during the passage of trains (Fig. 4). In total, twenty seven measurements were conducted for the passenger and the freight trains. -A1v/L, A1v/P -vertical acceleration of the midspan (left and right girder respectively), -A1h/L -horizontal acceleration of the mid-span (left girder), -A2v/L, A2v/P -vertical acceleration in the ¼ span length (left and right girder respectively), -A2h/L -horizontal acceleration in the ¼ span length (left girder), -UG1/L, UG1/P -vertical displacement of the mid span (left and right girder respectively).

Modal parameters identification
Modal parameters (frequencies, modal damping and mode shapes) were identified based on the free-decay responses. The free responses were collected just after the passages of trains. Eigensystem Realization Algorithm (ERA) was used for system identification. Acceleration signals from all measuring points were applied. The Average Normalized Power Spectral Density function (ANPSD) was also calculated for better verification of the results.

Eigensystem Realization Algorithm (ERA)
In 1985 Juang and Pappa proposed an ERA algorithm for modal parameters identification of LTI systems [25]. ERA is a direct, time-domain system identification based on the evolution of Ho-Kalman minimum realization problem [26]. It allows to determine the state-space model based on the finite-time and noisy measuring data. Linear discrete-time state-space model with initial condition x 0 = x(0) = 0 is described by the system of equations: where A, B, C and D are the state (transition), input, output and transmission matrices, respectively. The state of the system at time t k = kt, (t -time step, k  N) is characterized by the state vector x k = x(t k ). The external control and the system response are described by the input vector u k = u(t k ) and the output vector y k = y(t k ), respectively. In practice system matrices are unknown and are the result of the identification process. The basic principle of ERA involves identification of system matrices A, B and C from the free-vibration responses. The ERA method starts from forming Hankel matrix based on the measured Markov parameters. Let the value y k,p be the response at k-th time step in p-th sensor. Markov parameters represent a combination of instantaneous signal values from each measuring point. The response in k-th time step from all measuring points can be put into the vector (Markov parameter) y k = col{y k,1 , y k,2 , y k,3 , …, y k,p }. Markov parameters form the following general Hankel matrix: . (2) The method parameters  and  decide about the size and shape of the Hankel matrix, i.e. number of samples used. In the special case, for k = 1, we obtain: y y y y y y y y y y y y H y y y y Estimation of matrices A, B and C is possible using the Singular Value Decomposition (SVD) of the matrix H(0): where S (p×β) is a rectangular diagonal matrix containing singular values of H(0) in decreasing order Columns of matrices U (p×p) and V (×) are the left and right singular vectors, respectively, corresponding to each individual value of  i . The order of the system is specified by choosing n first (largest) singular values. The plot of  i may be helpful in such selection [27]. Because of the presence of the measurement noise, singular values can be divided into large (significant) and close to zero values. Significant values are associated with the real modes while the others can be computational or be the result of noise and should be rejected. Therefore, the square submatrix S n (n×n) with significant singular values can be extracted from the matrix S. Corresponding selection of the first n columns of matrices U and V leads to submatrices U n (p×n) and V n (×n). It can be noted that the reduced matrices U n , S n , and V n approximate the fullsize Hankel matrix. Finally, the minimum realization of the discrete LTI system of order n is defined by the triple: where matrices A e , B e , C e are the estimated quantities.
Auxiliary matrices E r (×p) and E p T (p×p) contain an appropriate number of identity and zero matrices of order p: To extract frequencies, damping and mode shapes it is necessary to transform the system matrices (5) from physical into modal coordinates. This can be done using eigenvalue decomposition: where  is a diagonal matrix containing the discretetime complex conjugate eigenvalues and  contains the corresponding eigenvectors as columns. The matrix  contains the information about eigenfrequencies and modal damping. Before extracting modal parameters the matrix  must be transformed from discrete-time into continuous-time form:  c = f s ln(), (f s is the sampling frequency). The i-th eigenvalue and it's complex conjugate become [28]: where  i and  ni are the damping ratio and natural frequency of the i-th mode, respectively. From (7) where |·|, Re(·), Im(·) are the modulus, real and imaginary part of a complex number. Mode shape vectors can be obtained from the column vectors of the estimated output matrix written in modal coordinates: C e m = C e  Each coordinate becomes: where i = 1, 2, 3, …, p and j = 1, 2, 3, …, n. Because of the inevitable measurement noise, a clear border between significant and close to zero singular values not always exists. Besides the noise and other possible signal distortions introduce uncertainty to the results [29]. For this reason it is recommended to overestimate the system order (rank of the Hankel matrix) [30]. So in order to distinguish the true modes from computational ones the Modal Amplitude Coherence (MAC) of the observability and controllability matrices can be applied [31]. The MAC criterion checks the compatibility in time between each individual modal component response of the identified model (vectors q e i ) and the real system (vectors q i ) and is defined as normalized dot product between these two vectors: A close to unity MAC value means that both vectors coincide, so corresponding modal parameters characterize the identified model (they are classified as structural). The values of MAC ≥ 0,97 are acceptable in most cases [32].

Averaged Normalized Power Spectral Density
The ANPSD function for the set of p signals is defined as [33]: where f k is the k-th discrete frequency (k = 1, 2, …, N), PSD i is the auto-power spectrum (Power Spectral Density) of the i-th channel (i = 1, 2, …, p), N is the number of discrete frequencies. The PSD function for ith signal can be obtained as: where FFT i is the Fast Fourier Transform of the i-th signal and "*" means Hermitian conjugation.

Identification results
The representative free response signals after the passage of EN57 train are shown on Fig. 6 and Fig. 7.
Corresponding stabilization diagram of ERA is shown on Fig. 8 together with the ANPSD function plotted. It can be seen, that MAC criterion is satisfy only for one mode of the frequency of 22,18 Hz. The identified mode corresponds to the bending theoretical mode of the natural frequency of 21,91 Hz (FEM solution). Modal parameters of the identified mode are summarized on Fig. 9.

Model validation
The FE model was validated based on the results of the 1 st identified mode. Natural frequencies depend on the mass and stiffness distribution of the structure. In the commercial FEM software stiffness and mass matrices are generated on the basis of defined cross-sections, material properties and the applied discretization. Theoretical frequencies can be adjusted to the measured ones by the mass actualization. It is a simple and commonly used technique. A slightly more difficult is to determine the coefficients of damping model. The Rayleigh damping model is implemented in SOFiSTiK software. In the model, the damping matrix is introduced by a combination of the mass and stiffness matrices: where a and b are the coefficients to be determined. Determination of Rayleigh coefficients is possible based on frequencies and damping ratios of two identified modes. In practical applications, if only one mode is identified, the damping model reduces to individual components of the sum, i.e. C = aM (mass proportional damping) or C = bK (stiffness proportional damping). In the case of Radunia bridge the stiffness proportional damping model was applied with b calculated as: To verify the correctness of damping definition logarithmic decrement (LDD) of the theoretical response was calculated and compared with LDD of the identified model (ERA). The identified model response for the first mode (ERA) is shown on fig. 10 in comparison with measurements. The theoretical response was obtained by impulse excitation applied to the mid-span of the bridge model.
The results of model validation are summarised in Table 1. The achieved level of compliance was considered as sufficient.

The Train
The EN57 traction unit was adopted as the moving load (see Fig. 4). Despite the aged construction it is still one of the most popular regional trains in Poland. Basic RSR configuration consist of three carriages: middle "engine" carriage (S) and two external "camshaft" carriages (R). In Ref. [34] the singular load models were defined on the basis of half of EN57 R carriage. In current study the simplified mechanical models of vehicle are: series of moving forces, series of moving single-mass oscillators, series of moving double-mass oscillators (Fig. 11). Individual loads of each model were put in the wheelsets distance in length and the wheels distance across the track width. Therefore, each carriage was modelled by eight individual forces (oscillators). In the case of the single-mass oscillators total mass of the carriage was assumed to be sprung. The double-mass oscillators model is more accurate. The wheelset mass is divided from the total mass and assumed to be unsprang. Fig. 11. Load models of EN57 train: a) scheme of RSR configuration, b) series of moving forces, c) series of moving single-mass oscillators, d) series of moving double-mass oscillators Dynamic parameters of the train were identified in [35] with the reference to each individual load model. These parameters are given in Table 2. Designations are as follows: P S , P R -load value for S and R carriage, respectively, M S , k S , c S -mass, spring stiffness and damping coefficient of the single-mass oscillators for S carriage, M R , k R , c R -mass, spring stiffness and damping coefficient of the single-mass oscillators for R carriage, M S1 , M S2 , k S1 , c S1 -sprung mass, unsprung mass, spring stiffness and damping coefficient of the double-mass oscillators for S carriage, M R1 , M R2 , k R1 , c R1 -sprung mass, unsprung mass, spring stiffness and damping coefficient of the double-mass oscillators for R carriage. The mass of passengers was also added to the structural mass of carriage. Number of passengers were referred to the number of seating and a human weight was assumed as 80 kg.

Bridge-vehicle interaction
Numerical simulations were performed using Newmarkβ procedure. According to convergence analysis a time step was assumed as Δt = 0,002 s. This value was enough for the solution stability. The RSR configuration of EN57 train was applied in simulations. The train speed was 49 kph. A comparison between measurements and simulations was of the main interest. Basic assumptions are as follows: the mass of each carriage is evenly distributed, only the vertical motion of the vehicle is possible, the train moves with a constant speed, vibrations of the vehicle and the bridge are measured from the static equilibrium position. On Fig. 12 visualization of the numerical model is shown for the double-mass oscillator case. On Fig. 13 the total response of vertical displacement is compared. Figures 14-16 show the comparison between measured and numerical results of vertical acceleration for individual vehicle models.   The vibration results of vertical displacement show a good compliance between measurements and simulations. Some differences occur in the maximum response for the first and last bogie. A possible reason is the lack of symmetry of the R carriage, as was assumed in simulations. The displacement level for all vehicle models is similar, but looking at the character of vibrations one can see that the response obtained for the double-mass oscillators load model is the closest to measurements.
A significant difference between theoretical models is visible in the vertical acceleration response. The series of moving forces and the series of moving single-mass oscillators significantly lower the results. The results for the double-mass oscillators are in greater compatibility with measurements. A root-mean square value (RMS) for each acceleration response was also calculated (Fig.  17). The RMS value is the measure of signal energy and is often used as the useful indicator for more objective comparison between two signals. In the calculations the forced vibration data was taking into account (time interval 0 ÷ 5,3 s). On Fig. 18 the vertical displacements of the mass of the oscillator for both sprung models of vehicles are compared. The same is depicted on Fig. 19 for the acceleration responses. From the above one can see that the response level is similar for both models of vehicle. The vibration character reminds a free-decay response. The time when the oscillator leaves the beam is 0,75 sec.

Conclusions
The bridge-vehicle dynamic interaction is considered in the paper. Physical parameters of the model were taken from the existing bridge and the existing railway vehicle. Field measurements were conducted in operating conditions for response data collection and system identification. The knowledge of existing dynamic parameters allows for the correct model validation.
Proper definition of train load is an essential part of the dynamic analysis of railway bridges. The concentrated force formulation omits the inertial terms of moving mass. Moreover, the influence of the dynamic bridge-vehicle interaction cannot be taken into account. An accurate sprung model of railway vehicle should concern at least two mass discretization with one mass unsprung (wheelset's mass). The sprung-only load models are more appropriate for road bridges.
Passing over the small-span bridge works like an impulse excitation for the vehicle. However, the vehicle acceleration response is ten times smaller than the bridge span response (for the conditions considered in the paper). It should be noted that the track stiffness outside the bridge was assumed to be infinite. The actual track stiffness could change the results quantitatively.
The connection between rails and longitudinal beams was assumed to be stiff for simplicity. In order to define the actual track parameters the spring stiffness can be validated based on in-situ measurements.
The developed models of vehicles were applied to the small-span steel bridge of a simple structure. However there is no limitation of using them for more complex structures. The vertical track parameters and rail irregularities can also be taken into account.