Calibration and validation of the dynamic response of two slab track models using data from a full-scale test rig

For the development of accurate and reliable simulation models, the procedure of calibration and validation against measurement data is essential. In this paper, a finite element model and a waveguide finite element model of a slab track are calibrated and validated against hammer impact measurement data from a full-scale test rig. The finite element model is three-dimensional, where the rails are modelled as Rayleigh – Timoshenko beams and the concrete slab and support layer are modelled using linear shell elements. In the waveguide finite element model, a constant track cross-section described by two-dimensional finite elements is assumed, and the vibration in the direction perpendicular to the cross-section is described by propagating waves that are decaying expo- nentially. Measured frequency response functions (FRFs) are compared with the corresponding calculated FRFs from the two modelling approaches. The calibration is conducted in two steps using (i) a parameter study and (ii) a genetic algorithm. For multiple excitation positions and sensor locations, both calibrated models capture the trend of the Single-Input Multiple-Output measurements with rather small deviations compared to the overall dynamic range. This implies that both models can successfully represent the dynamic response of the test rig and can be considered as validated. Conceptualization, Resources, Supervision. Wanming Zhai: Conceptualization, Resources, Supervi-sion, Project administration.


Introduction
In high-speed railway applications, the usage of slab track has increased in recent decades [1,2]. Compared to ballasted track, slab track has several advantages including reduced maintenance costs and higher track stability. The main drawbacks with slab track are increased construction cost and higher noise levels. To build optimised slab track structures, the railway industry is dependent on accurate and reliable track models. Depending on the application, different types of models have been developed [3][4][5]. Traditionally, finite element models based on beam elements have been used, cf. [6][7][8][9][10][11][12][13][14][15]. For applications related to noise radiation modelling, the waveguide finite element method has been applied to model the rail dynamics up to high frequencies, e.g. 4 kHz in Ref. [16] and 80 kHz in Ref. [17]. Recently, several threedimensional (3D) slab track models have been developed, cf. [18][19][20][21][22][23][24][25][26]. The main benefit of using a 3D model instead of a two-dimensional (2D) one is that the influence of non-symmetric wheel load excitation and support conditions can be studied, whereas the main drawback is the increased computational cost. However, accounting for the periodicity of the track structure and solving the dynamics in the wavenumber domain, cf. [27][28][29][30], can decrease this computational cost.
For the development of accurate and reliable simulation models, the procedure of calibration and validation against measurement data is essential. In the literature, different measurement strategies have been described depending on what part of the slab track that was analysed, and whether the measurements were carried out in the field or in a laboratory. Often an impact hammer including an integrated force transducer is used to excite the structure in a Single-Input Multiple-Output (SIMO) test. The response is typically measured using piezoelectric accelerometers [31,32]. A major benefit of laboratory tests as compared to field tests is that track properties can be more controlled. However, dynamic track models have also been calibrated based on measurements in the field [33][34][35].
In a study by Cox et al. [36], a dynamic characterisation of different floating slab track systems and direct fixation fastening systems was performed by measuring receptances (displacement over force) in a fullscale test rig. In the test, 12 m long rails were used and boxes of sand were placed at the ends to reduce boundary effects. In 2017, another full-scale test rig was presented by Wang et al. [37]. In their study, the dynamic performance of the China Railway Track System (CRTS) series was analysed by conducting so-called wheel-drop tests. Zangeneh et al. [38] analysed the dynamic response of portal frame railway bridges and used a model updating approach to calibrate their finite element model. Sainz-Aja et al. [39] calibrated a slab track model based on finite elements by using measurement data from a full-scale test rig. The measurements included variations of load amplitude and frequency, with a focus on low frequencies up to 5.6 Hz. Tarifa et al. [40] analysed the fatigue life of reinforced concrete slabs. Their study consisted of several tests including three-point bending tests of full-scale slabs. When a simulation model is calibrated against measurement data, the calibration is typically conducted using an optimisation algorithm. Depending on the application, different optimisation strategies can be used. If the objective function is convex and has continuous derivatives, classical gradient-based methods can be used. Gradients that can not be calculated and expressed in closed form can be approximated by finite difference methods. An alternative to gradient-based algorithms is to use a Response Surface Methodology (RSM). As an example, Shevtsov et al. [41] used an RSM to optimise railway wheels. Although gradient-based algorithms and RSM can be powerful optimisation tools, they fail in several engineering applications [42]. Typical reasons are high computational cost, multiple objective functions and complexity of the objective function(s). A remedy can then be to apply a stochastic optimisation method, e.g. a genetic algorithm (GA). GAs have been used frequently in the railway context in different calibration and optimisation processes [35,[43][44][45][46].
In this paper, results from SIMO hammer impact measurements on a section of the CRTS III design are presented. The measurements were performed in the full-scale test rig described by Wang et al. [37] and Zhai et al. [47]. From the measured receptances, the dynamic behaviour of the CRTS III is analysed. Furthermore, a 3D finite element model and a waveguide finite element model are calibrated and validated with the measurement data. The calibration consists of two steps including (i) a parameter study and (ii) a genetic algorithm. From the comparisons between the models and the measurements, the applicability and accuracy of the presented calibrated and validated slab track models are shown. Finally, the benefits and drawbacks of the two different models are discussed.

Test set-up
In Section 2.1, predefined parameter values and the geometry of the analysed track section are presented. Details about how the measurements were executed are described in Section 2.2.

Track geometry and parameters
The measurements were performed in the State Key Laboratory of Traction Power at the Southwest Jiaotong University in Chengdu, P.R. China. The full-scale test rig has a total length of about 55 m and includes five sections with different types of slab track and one section with ballasted track. An overview of the test rig is shown in Fig. 1(a). This work focuses on the section with the China Railway Track System III (CRTS III), which has a total length of 16.5 m. The CRTS III section is located in the central part of the test rig, while the adjacent track sections are a CRTS II section and a floating slab track section.
The geometry of the track is presented in Figs. 1(b) and 2. The considered track section consists of three pre-fabricated concrete slabs. Each concrete slab has a length of 5.5m, which corresponds to eight rail seat distances. An adjustment layer made of self-compacting concrete (SCC) connects the slabs with the concrete support layer. The support layer rests on a soil layer, which is designed to have homogeneous properties for several metres of depth. The used rail profile, CN60, resembles closely the standard UIC60 rail profile. The fasteners are of type WJ-8, which is a common fastening system on Chinese ballastless highspeed lines. The parameter values of the studied track are given in Table 1.
Deviations in geometry and boundary conditions between the test rig and the models need to be considered in the model validation. One such deviation is that the rails are discontinuous within the investigated section of the test rig. The rail segments are joined using suspended rail joints (fishplates) with six bolts, see Figs. 2 and 3. The locations of the rail discontinuities are about three rail seat distances away from either end of the CRTS III section. A thin concrete floor resting on the soil is cast on both sides of the concrete support layer (not accounted for in the models). Furthermore, the exact properties of the connection between each pair of concrete slabs are unclear.

Measurement execution
The aim of the measurements was to obtain transfer functions describing the vibration of the slab for an excitation on the rail. A sledgehammer with a steel tip and a weight of 8 kg was used to generate an excitation impulse on the rail. The hammer was manually guided to a marked excitation position. A typical amplitude spectrum of the voltage signal of the hammer is shown in Fig. 5. As visible, the steel tip of the hammer produced a fairly flat response up to approximately 7 kHz; the narrow peaks are discussed below. It would have been desirable to use a hammer with a slightly softer tip in order to excite only the frequency range that was used for the calibration, but this was not available. An alternative of using a hammer for excitation is to use a so-called wheel- drop test, cf. Ref. [37]. In this paper, the usage of wheel-drop tests was not considered since it is difficult to control and determine the input force.
The acceleration was measured at five locations on the slab and support layer as indicated in Fig. 4. The sensors were numbered (1)(2)(3)(4)(5) with squares indicating tri-axial sensors, whereas circles indicate uniaxial, vertical sensors. Sensor 3 is placed below the rail, next to the rail seat, see Fig. 3(b). The excitation positions are shown in grey. Each excitation was carried out both vertically on the rail head and laterally on the inside of the rail head. In order to mount the sensors on the track, each sensor was screwed to a steel plate with dimensions 3 cm × 3 cm × 4 mm, which was glued on the concrete using a two-component glue, see Fig. 3(b). The influence of the mass of the steel plate on the measured transfer function is considered negligible in the frequency range below 2 kHz. The measurement equipment is listed in Table 2. All sensors measured the vertical acceleration and the tri-axial sensors were used to additionally measure the lateral acceleration.
At each excitation position, the measurements were executed four times. The force signal of the hammer was used to trigger the measurement with a block size of 4 s at a sampling frequency of f s = 16 kHz.
The recorded time signals were transformed to the frequency domain. Harmonic peaks were found in the spectra as shown in Fig. 5. With the peaks being spaced by 512 Hz, and exactly one frequency bin wide, they were considered non-physical and a moving median operation with a width of 1 Hz was applied to smooth the spectra. The data acquisition system is assumed to be responsible for this phenomenon, as these peaks occurred in all channels and measurement positions throughout the measurement campaign. From the measured accelerations and the force signal, transfer functions were calculated. For each transfer path, the magnitudes of the four measured transfer functions were averaged to produce one magnitude spectrum.
This method of measuring the frequency response has limitations      Fig. 6 shows the variance between four individual measurements and the mean value. Using four measurements per position proved to be an acceptable compromise between accuracy, repeatability, and time efficiency. The mean receptance used in the following is indicated in black. Note the increasing span between the highest and lowest measured receptance towards the low end of the frequency range. This variance is reflected in the coherence of the measurements. Fig. 7 shows coherences of some measured transfer functions as an indication of the signal quality and linearity. These coherences are representative of most measurements. It is clear that the quality of the measured signal below 60 Hz is low. Therefore, frequencies below 60 Hz were excluded from the further processing.

Modelling and tuning of parameters
The measured dynamic responses are compared to simulations with two modelling approaches, namely a discretely coupled waveguide finite element model (WFE) and a finite element (FE) model. In the following, the models are introduced and compared to the measurements. For the calibration and validation, the sensors and excitation positions shown in Fig. 4 were used. Four of the excitation positions (I c , III c , I f and III f ) were used for the calibration, while the other two excitation positions (II c and II f ) were used for the validation. For both the WFE and FE models, the calibration was performed using a parameter study and a genetic algorithm. This procedure is described in Section 3.1. The description of each model and the simulation results are given in Sections 3.2 and 3.3.

Calibration strategy
In the procedure for calibrating the models to the measurement data, the first step is to decide what parameters to include in the calibration process, see Section 3.1.1. In this paper, the calibration consists of two steps. In the first step, the stiffness of the rail pad is calibrated by a parameter study, see Section 3.1.2. In Section 3.1.3, the second step of the calibration is described, where the damping of the rail pad and the stiffness and damping of the soil are determined using a Genetic Algorithm (GA). An alternative to using the adopted two-step calibration process would be to include all four design parameters in the GA. The performance of such a calibration process has been investigated and it was found that the mode around 150 Hz was calibrated with higher accuracy when the adopted two-step calibration process was used. Since the result of the SIMO measurements is a large number of receptances that are used in the calibration and validation, the information needs to be compressed for visual presentation. In Section 3.1.4, a description is given of how the results are visualised in the upcoming parts of the paper.

Parameter selection
Several issues need to be considered when selecting what track parameters to include in the calibration. The most important ones are (i) the influence of the parameter on the considered track response(s) in the frequency range of interest and (ii) the uncertainty of the parameter values. There is no need to calibrate a parameter that already by its design can be specified with high accuracy. In this paper, the models will  be compared to the measurement data in terms of receptances.
In general, the stiffness and damping of the support and resilient layers in the track meet both criteria in terms of significantly affecting the receptance and including an inherent uncertainty. In this paper, frequencies from 60 Hz up to 1500 Hz are studied. This frequency range is selected since the coherence for these frequencies is high and since both models can be expected to give accurate results in the selected frequency range. From a sensitivity analysis, it was concluded that the stiffness and damping of the soil and the rail pad have a significant effect on the receptances in the studied frequency range. Therefore, these parameters are included in the calibration process. From the sensitivity analysis, it was also found that the influence of the stiffness and damping of the SCC layer is negligible since this layer is significantly stiffer than the rail pad and the soil. The properties describing the SCC layer are, hence, not included in the calibration process. Finally, the sensitivity analysis did also show that the rail parameters, as well as the density, dimensions, Young's modulus and Poisson's ratio of the concrete parts, have a significant effect on the receptances. These properties are, however, not included in the calibration process since their parameter values can be assumed to be specified with relatively high accuracy.

Parameter study
In Section 3.1.1, it was concluded that four properties will be included in the calibration process: (a) the stiffness of the rail pad, (b) the damping of the rail pad, (c) the stiffness of the soil and (d) the damping of the soil. Note that the number of parameters that will be tuned is different between the finite element (FE) model and the waveguide finite element (WFE) model since both vertical and lateral parameters are included in the WFE model, whereas only vertical parameters are included in the FE model.
Consistently for all the measurements with excitation and response in the vertical direction, a peak in the response can be seen around 140-160 Hz. This peak corresponds to the cut-on frequency of a vertical rail bending wave and is strongly affected by the rail pad stiffness. Considering the dispersion relation of the waves in the rail, the vibration at this cut-on frequency is identical with the vertical rigid body motion of the rail. A description of typical rail cross-sectional mode-shapes is presented by Thompson [48]. For a given design of ballasted track, the cut-on frequency for this vertical rail bending wave was found to occur at around 200 Hz [31]. The other parameters used in the calibration, i.e. the damping of the rail pad and stiffness and damping of the soil, do not significantly affect the frequency of this wave. Hence, in order to reach a good agreement between simulations and measurements around 140-160 Hz, the rail pad stiffness has to be tuned.
The exact resonance frequency of this mode-shape varies slightly between the different measured receptances. Based on an average over all measured receptances included in the calibration, the resonance frequency is 147 Hz. To determine the rail pad stiffness that meets this resonance frequency, parameter studies were conducted for both the FE model and the WFE model for vertical rail pad stiffnesses spanning from 10-50 kN/mm. For each of the considered rail pad stiffnesses, all different receptances in the vertical direction to be used in the calibration were calculated and the average resonance frequency of the mode was determined. The influence of the rail pad stiffness on the frequency of the mode is shown in Fig. 8. From the figure, the best match between the simulation and measurement using the FE model is achieved if the rail pad stiffness is k rp = 34 kN/mm. A marginally higher rail pad stiffness of 37 kN/mm was found for the WFE model.
No specification of the lateral stiffness for this rail pad has been found in the literature. In a preliminary study (not shown here), a genetic algorithm was used with a larger set of parameters, including the vertical and lateral rail pad stiffnesses. From this study, it was found that the vertical and lateral rail pad stiffnesses typically converged to similar values. Thus, in the WFE model, the same stiffness is used for both the vertical and lateral stiffnesses. Note that this does not conform with specifications for other types of rail pads, as the lateral stiffness is often specified as lower than the vertical one [49].
The calibrated values of the vertical rail pad stiffness (34 kN/mm and 37 kN/mm) are slightly higher than the static stiffness value provided by Wang et al. [37] (25 kN/mm). However, the stiffness of the rail pad depends on the frequency of excitation, the magnitude of the preload, temperature, strain amplitude and strain history [50]. In particular, it has been noted from other measurements that the ratio between the static and dynamic rail pad stiffness may vary by a factor of 2-8 [51]. With this in mind, the tuned dynamic rail pad stiffnesses 34 kN/mm and 37 kN/mm seem to be reasonable.

Genetic algorithm
The damping of the rail pad and the stiffness and damping of the soil affect the magnitude of the receptance over a wide frequency range. An efficient approach for more complex calibration problems with several interacting parameters is to use a genetic algorithm (GA). In this paper, the selected objective function to be minimised is inspired by the objective function used by Andersson and Abrahamsson [52]. For each considered receptance, the logarithmic difference e between the receptance magnitudes of the model H X i and the measurement H A i is calculated as where n p is the number of receptances used in the calibration and n ω is the number of considered frequency lines. In this paper, 138 frequency lines were considered, spanning from 60Hz to 1500Hz. By collecting all values of e i,k in a vector e, the objective function can be written as where Q is a non-negative definite weighting matrix with dimensions n p n ω × n p n ω . In this paper, since the parameters used in the GA have a significant influence on the receptances over a wide frequency range, the weighting matrix is defined as the identity matrix. In the calibration of the design parameters, only the magnitudes of the receptances are considered. Although the procedure applied here can easily be extended to account for the errors in both magnitude and phase, see Ref. [52], it was found that considering only the error in magnitude resulted in the best agreement between the two models. In this paper, both models are calibrated using the standard genetic algorithm implemented in Matlab with a population size of 50 for each of 30 generations (iterations). In the algorithm, the initial population is generated by randomly selecting values within the predefined bounds of the design variables that are given in Tables 3 and 4. From the initial population, the algorithm creates a sequence of new populations. The main idea is that a member in a previous population has a higher probability of being a member in the next population if its setting leads to a low value of the objective function. In practice, this is achieved by using the concepts of elitism, mutation and crossover, see Ref. [42].
To determine the proper population size and number of generations, a convergence study was conducted to make sure that the minimum of the objective function was found. In Fig. 9, the convergence of the models can be seen. In the figure, both the smallest value and the mean value of the objective function are shown as a function of generation number. The values of the objective function have been normalised with the smallest value of the GA in the final generation. From the figure, it can be seen that the algorithm has converged, which implies that terminating the algorithm after 30 generations is a valid option. Note that the convergence shown in the figure is obtained by analysing one simulation for each modelling approach. However, it has been verified that the convergence rate of the other simulations used to calculate the mean value and standard deviation of the design variables presented in Sections 3.2.2 and 3.3.2 are similar.
For each member in the genetic algorithm, a new track model has to be generated. For the finite element model, this implies that the Python script that generates the track model in Abaqus has to be called multiple times from Matlab. By using parallel computations during each generation, the computational cost of the GA is reduced significantly. In Fig. 10, a flowchart of the GA for the finite element model is shown. For the waveguide finite element model, the same steps were used, but Abaqus was not involved. Instead, an in-house Python software was used for the calculation of the track receptance.
As briefly discussed above, an alternative of using the GA is to use a gradient-based algorithm. The main advantage of using a gradient-based algorithm is that the optimal solution can be found faster. In this paper, where the number of design variables is relatively small (3-4 depending on which model is used), both a GA and a gradient-based algorithm would probably work. However, since it cannot be proven that the objective function is convex, a minimum calculated with a gradientbased algorithm is not necessarily a global minimum. Since the computational demands when using the GA in this application are manageable, the GA was considered as the best optimisation method for the calibration.

Result visualisation
The quality of the match between the measurements and the tuned models is analysed by comparing the individual receptances. However, due to the extensive number of receptances (20 for the FE model and 32 for the WFE model), this information needs to be compressed for visual presentation. In this paper, the similarity between the measured and simulated receptances is determined by calculating their difference in dB, denoted ΔH dB . This similarity measure is calculated for each exci- Observe that the receptance is expressed as 10log 10 instead of the usual 20log 10 . The difference spectra are then averaged for each one-third octave band and displayed in a surface plot using the colour coding shown in Fig. 11. Here, the brown colour indicates an over-estimation by the model, while the blue colour means an under-estimation.

Discretely coupled waveguide finite element model
The waveguide finite element (WFE) method can be used to model structures that are sufficiently long in one dimension and have a constant cross-section along this dimension [16]. Such a structure acts as a waveguide with propagating, exponentially decaying waves. The WFE method uses the assumption of an infinitely long waveguide. Thus, when discretising the three-dimensional structure, the longer dimension is described by wave functions, and only the cross-section needs to be discretised with two-dimensional finite elements. This process vastly   reduces the required number of degrees of freedom (DOFs) in the numerical model. The model has been developed for the prediction of noise radiation due to wheel-rail contact loading and vibrations in the vertical and lateral directions [53,54].

Model description
Here, the rail and the slab on concrete support are modelled as two separate WFE-models, similar to Ref. [53]. Their cross-sections are meshed using 9-node, quadrilateral isoparametric elements. The dimensions of the model correspond to the dimensions given in Fig. 1(b). Each WFE model is formulated as presented e.g. in Ref. [16], relating the forces and displacements in the discretised domain. Fig. 12 shows the nodes of the 2D meshes of the rail and the concrete and soil parts of the modelled track cross-section. Thick lines indicate a fixed boundary condition for the nodes on the boundary. In early-stage investigations, elements on the outer top edge of the soil layer were used to represent the thin floor layer present adjacent to the test rig. However, due to a low influence on the simulation results, these elements represented soil material in the final optimisation.
Assembling the element matrices in the global matrix system gives the expression with the mass matrix M and the nodal displacements and forces Ũ and F , respectively. The matrices K i are generalised stiffness matrices and κ is the wavenumber. The time dependency e jωt is used with the circular frequency ω. The equation has been transformed to the wavenumber domain using the Fourier transform.
The free response F = 0 → is solved by prescribing a frequency ω and solving the quadratic eigenvalue problem. This generates complexconjugate pairs of eigenvalues corresponding to wavenumbers κ n , representing propagating, decaying waves, and corresponding left and right eigenvectors Ũ nL and Ũ nR . The total displacement response of the structure to the forced excitation is calculated by a superposition of these waves. At the origin of the coordinate direction along the structure (index 0), this superposition is expressed in the wavenumber domain as with and D(κ n ) = − 2κ n K 2 − jK 1 .
The force vector F 0 contains the nodal excitation forces over the length of the rail expressed in the wavenumber domain. The coupling of the rail with the other components of the track is formulated in the spatial domain. The calculation includes four steps. First, a transfer function matrix is calculated separately for each waveguide using the WFE approach described above. Transfer functions between all coupling points are evaluated for each waveguide. The two waveguides are coupled at n l locations, i.e. at all rail seats. These do not need to be uniformly spaced [30]. In each of those locations, the structures can be coupled in multiple (n n ) degrees of freedom, e.g. in multiple nodes of the FE-mesh along the width of the rail and in both vertical and lateral directions. This leads to a total number of n t = n l n n connections between the structures.
The nodal displacement u at any point i on any of the structures ξ (here ξ ∈ [1, 2]) is a superposition of the response due to an excitation force F e,ξ on the structure and the response due to the reaction forces F in the coupling points. This can be written as with the index g (g = 1,2,…,n t ) indicating the coupling point. The cross receptance H ξ i * describes the displacement at the location i for a unit force input F * . These transfer functions are evaluated individually for the rail and the remaining parts of the track according to Eq. 5. The evaluation position i in Eq. 8 can be iteratively placed at each coupling point. This generates n t equations, of which each one has the displacement u i and the connection forces F as unknowns. Assembling these n t equations into one system of equations generates a symmetric transfer function matrix H ξ of size n 2 t . Secondly, a model for the connection itself is introduced. The rail pad is modelled using linear springs, in which damping is included by assuming a complex stiffness. Each spring is represented by its receptance H p such that at each coupling point g. The vertical and lateral directions are modelled as decoupled.
In the third step, the transfer function matrices from step one and the connection matrix are assembled and formulated in a system of equations in order to calculate the forces acting on each rail pad. Inserting Eq. 8 for each structure into Eq. 9 and assembling the receptances H p to a single diagonal matrix H p of size n 2 t produces with the vector of connection forces F . Note that the excitation term on the right hand side contains only the transfer functions of the structures which are excited.
Finally, this linear system of equations is solved for F . To connect this to the wavenumber domain calculation above, the reaction forces are expressed as a wavenumber spectrum at the origin, which can be introduced as an excitation in Eq. 5. This enables the calculation of the response due to the connection forces. A harmonic excitation of the coupled system can be implemented by superposing the free response of the excited structure with the response due to the reaction forces.
In this paper, a point load on a node of the mesh is assumed as indicated in Fig. 12 (b). The three nodes across the foot of the rail in which the rail and the slab are coupled are indicated in the same figure.
The bodies are coupled in both lateral and vertical directions (n n = 6). The number of rail seats is chosen to be n l = 18, corresponding to the 18 rail seats in the considered track section.
In total, 64 transfer functions were included in the calibration of parameters, corresponding to eight excitation positions and eight response positions. The considered excitations are at I c , III c , I f and III f in both the vertical and lateral directions, see Fig. 4. The five sensors shown in Fig. 4, three of which are measuring in both the vertical and lateral directions, are included. The calibration was carried out following the description in Section 3.1. As discussed in Section 3.1.2, the tuned (vertical and lateral) stiffness of the rail pad is 37 kN/mm. The upper and lower bounds for the optimisation variables (i.e. the rail pad and soil damping values) presented in Table 3 were chosen based on engineering judgement such that values commonly referenced in literature were enclosed, cf. Ref. [37]. Due to the stochastic nature of the optimisation methodology there is a variance in the optimisation results. Six optimisations were conducted and the mean and standard deviation of each variable is presented in Table 3.
The vertical soil stiffness k sg , described by Wang et al. [37], is formulated as a stiffness per unit area ([k sg ] = MPa/m). A rough comparison to the optimisation result can be made by assuming a linear elastic material model. The Young's modulus E, the stress σ and strain ε are related by with a force F that is acting on the material with area A and height l 0 . Rearranging to solve for E gives an approximation of the Young's modulus. In this model, where a soil thickness of l 0 = 0.33m was chosen, the value in [37] would correspond to62.3 MPa, which is in the same range as the calculated 86.5 MPa. Note that damping is introduced via a complex stiffness and so the damping coefficients do not have a unit. In a pre-study to the optimisation, it was found that the discontinuities in the rail have a significant effect on the measured response in the slab due to the modal behaviour of the rail section. These additional modes can not be directly modelled with a WFE-approach. Instead, a discontinuity in the support was included at the locations of the rail ends. There, the rail pad stiffness and damping of one rail seat closest to either end of the rail section was increased by an empirically determined factor three to introduce a discontinuity in the longitudinal direction and to better approximate the measured responses. The presence of this discontinuity also means that including neighbouring track sections in the model does not influence the result considerably and thus, these are not included in the model.

Optimisation results
As described in Section 3.1.4, the resulting match between the simulated and measured receptances is visualised by a colour scheme in one-third octave bands. In the following, a selection of these matches is presented. Fig. 13 shows the transfer functions for three vertical excitations on both rails at different distances to the sensors. Note that the upper three rows in each sub-figure correspond to lateral response channels. For the 50 Hz one-third octave band and below, it is observed that the model generally underestimates the response. In the 64 Hz to 125 Hz one-third octave bands, the model tends to slightly overestimate the response for both lateral and vertical channels. Above that, the model tends to slightly underestimate the response. Fig. 14 shows the similarity measure for the same excitation points, but for a lateral excitation. The darker colouring indicates larger differences between the model and the measurements. A fairly close match is observed for the lateral channels at frequencies above 80 Hz. For lateral excitation, it is noted that the vertical displacements are to varying degrees underestimated by the model.
In general, it can be noted that larger differences are observed for low frequencies. This is assumed to be a direct result of the low measurement accuracy in this frequency range as described in Section 2.2. In each of Figs. 13 and 14, the similarity measures are shown for three excitation positions. However, only positions I and III were included in the GA. The central figure, representing excitations at location II, shows matches of similar quality as the other two figures. As the calibrated model is able to predict receptances outside of the calibration process with the same accuracy, the model can be considered to be validated.
A notable observation is that a closer match is achieved when the excitation direction aligns with the measurement direction. This is further investigated in Fig. 15, which presents the receptances for one excitation position and one sensor position, with the excitation and measurement in lateral and vertical directions. In accordance with the findings from Figs. 13 and 14, the single-directional receptances (vertical to vertical or lateral to lateral) match more closely than the crossdirectional receptances. It can be assumed that these differences are partly due to the simplifications in the model. The simplifications include the assumption of a continuous slab and representing the complex behaviour of the rail pad by six linear, independently acting springs.
Finally, the obtained model is used to study the influence of the discontinuous rail. As described in Section 3.2.1, the rail pad stiffness of the outer rail pads is increased by a factor three to introduce a discontinuity in the rail support. Fig. 16 visualises this effect. It is observed that the prediction of the model produced by the GA follows the general trend of the measurement. However, when introducing the discontinuity, a visible modal pattern appears starting from about 80 Hz. This modal pattern matches some of the resonances and anti-resonances in the measured response. For this receptance, the updated model seems to resemble more closely with the measurement, especially around 300 Hz.

Finite element model
The second modelling approach that is considered in this paper is a finite element (FE) model. The model has been developed for the design of slab track structures considering the strength and fatigue life [55]. All three spatial dimensions are represented using finite elements, but only wheel-rail contact loading in the vertical direction is accounted for. Note that the FE-model can be extended to include also lateral dynamics, but this is outside of the scope of the paper. In Section 3.3.1, the model is described and in Section 3.3.2, the results of the calibration and validation are given.

Model description
The parametrised FE model has been developed in Abaqus using Python scripts and is described in detail in Ref. [55]. In Fig. 17, a schematic cross-section of the 3D FE model is shown. The dimensions of the model are according to the CRTS III system and are given in Figs. 1 (b) and 2. The rails are modelled as Rayleigh-Timoshenko beams. Regarding the modelling of the concrete slabs and support layer, it has been verified that the receptance characteristics are similar when using either shell or solid elements. In the calibration of the model, shell elements were used leading to a lower computational cost. A linear shell element (denoted S4 in Abaqus [56]) was employed and it has been verified that quadratic shell elements give similar results and are not required for this application. For a detailed comparison of shell versus solid elements and linear versus quadratic elements, see Ref. [55]. Furthermore, the influence of different mesh densities has been investigated. Depending on what frequency range to be studied, the required mesh density varies where a finer mesh is required at higher frequencies.
From the investigation, it was concluded that using a mesh density which has an average element side length of 14 cm leads to accurate results.
The rail pad, soil and SCC layer are modelled as distributed sets of non-interacting springs and viscous dampers. Since the rail pad distributes the force from the rail to the slab over a certain area, a set of springs and viscous dampers acting in parallel is used at each rail seat. This set distributes the load over an area corresponding to the side length of the rail pad in the longitudinal direction and the width of the rail in the lateral direction.
To determine the length of the track model, the trade-off between accuracy and computational cost needs to be considered. In the test rig, the CRTS III section has a length of 24 rail seat distances (corresponding to 16.5 m). Before the CRTS III section, there is a section of CRTS II and after the CRTS III section there is a section of a floating slab track, see Section 2 and Ref. [37]. In order to reduce the influence of boundary effects in the centre part of the model, track corresponding to six rail seat distances of both the CRTS II and the floating slab track were included in the model giving a total track model length of 24.7 m. Since the track sections of the CRTS II and the floating slab track are included in the model only to reduce boundary effects, they were for simplicity modelled in the same way as the CRTS III section.
As discussed previously, there are several fishplates connecting the different rail sections in the test rig. In the model, the fishplates are modelled as beam elements with a rectangular cross-section in a similar way as in Ref. [57]. At the joint, the rail is cut and a fishplate is added on either side of the rail. In the test rig, the fishplates were mounted to the rail with six bolts. In the model, the fishplates are connected to the rail using a tie constraint at each bolt location. The tie constraint imposes that no relative motion can occur between the rail and fishplates at the locations of the bolts. By comparing the receptances of the track model that includes the fishplates with another model where the rail is continuous, it was verified that additional modes are obtained when the rail is cut and fishplates are included (in particular in the high-frequency range above 500 Hz). In addition, it was found that each of the individual contributions of introducing a discontinuous rail and attaching fishplates to the rail affect the receptance significantly in the highfrequency range.
The calibration and validation of the model are performed using receptances of the track model calculated in the frequency domain. To calculate the receptances, the equations of motion are established as where u is a vector containing the displacement at all DOFs, F is the external load vector and M, C and K are the system matrices. Since only steady-state harmonic forces are considered when the receptances are calculated, Eq. 14 can be written as where E(ω) = − ω 2 M +jωC +K is the frequency-dependent dynamic stiffness matrix. To calculate the receptances to be compared with the measurements from the test rig, a vertical harmonic force F is applied at a node in the finite element model that corresponds to a hammer impact position in the measurements. By solving Eq. 15, all displacements are calculated. Let x i denote a vertical displacement at sensor i. For this test set-up, by extracting the displacements x i from u, the receptances are determined by

Optimisation results
As described in Section 3.1, the calibration is performed in two steps. In the first step, the stiffness of the rail pad was tuned using a parameter study and in the second step, the damping of the rail pad and the stiffness and damping of the soil were calibrated using a genetic algorithm (GA). From the parameter study, it was found that a vertical rail pad stiffness of 34 kN/mm should be used. In Table 4, the lower and upper bounds and the optimised results from the GA are shown. Since the GA is a stochastic optimisation methodology, there is some variability in the optimisation results. To investigate the extent of this variability, six simulations were used and the mean value and standard deviation of each design variable are presented in Table 4. To make sure that the selected bounds were reasonable, it was verified that the bounds give stiffness to damping ratios that are in the same order of magnitude as the ratios indicated by Nielsen [58]. Furthermore, it is confirmed that the optimised values of the design variables are not close to any of the bounds.
In Section 2, it was described that four excitation locations and five sensors were used in the calibration. For the finite element model, where no lateral dynamics is considered, this means that there are in total 20 receptances that can be used in the calibration. In Fig. 18, one example is shown when comparing the simulation and measurement results. In the figure, also the similarity measure for each one-third octave band is shown.
In Fig. 19, the similarity measure is shown for all of the 20 receptances that were used in the calibration. Each horizontal line corresponds to one receptance as a function of frequency in one-third octave bands. The results from Fig. 18 is shown as the top row in Fig. 19(c). Considering that the dynamic range of the receptance magnitudes is up to 100 dB, it is concluded that there is reasonably good agreement between the simulation and measurement results in Fig. 19. There are, however, some differences between the measurements and simulations. In particular, it can be seen that the model tends to over-estimate the  receptances at high frequencies (> 500 Hz). At these high frequencies, the influences of the fishplates and the boundary conditions between the slabs have a significant effect.
Using the calibrated parameter values of the track model, a validation of the track model has been conducted. The validation was done by comparing the measurement data with simulations when considering the excitation positions that were not used in the calibration (see II c and II f in Fig. 4). In Fig. 20, one example of a cross receptance is shown. Since this validation figure is similar to the calibration figure (Fig. 18) in terms of the similarity measure, the finite element model can be considered as validated. For completeness, Fig. 21 shows the validation version of Fig. 19. The fact that the appearances of Figs. 19 and 21 are similar strongly indicates that the model is validated.

Discussion
This paper presents two approaches for the modelling of the dynamic response of a slab track. Although both modelling approaches are based on finite elements, there are several fundamentally different assumptions. Here, a comparison of the methods in terms of their assumptions, advantages, limitations and numerical efficiency is presented.
As both models are based on the finite element method, numerical efficiency is an important factor leading to the need for different simplifications in the models. In the finite element (FE) model, the computational demands were reduced by modelling the SCC-layer and the soil as non-interacting springs and dampers. Modelling the concrete slab and support layer as linear shell elements instead of solid elements decreased the computational cost of the FE model even more. In the waveguide finite element (WFE) model, the length of the track is represented only by the assumption of propagating waves along its constant cross-section. In addition, since the geometries of the rail and the remaining parts of the track are evaluated independently, the global stiffness matrices of the WFE model are comparatively small. The independent calculation of the free response of the rail furthermore implies that this rail response can be pre-calculated since no rail parameters are altered during the calibration.
These simplifications have several implications regarding the advantages and limitations of each of the models. Due to the threedimensional mesh of the slab track when using the FE model, the boundaries between the panels and between track sections can be taken into account. Furthermore, it is possible to include a model of the rail joints within the structure. As shown in Section 3.2, this discontinuity has a large effect on the measured response. The most prominent advantage of the WFE model is that it can be used at higher frequencies compared to the FE model. The FE model, which uses a Rayleigh-Timoshenko beam description of the rails, gives accurate results up to about 1.5 kHz, whereas a WFE model of the rails has previously been used even up to 80 kHz [17]. In this paper, frequencies above 1.5 kHz were not considered, since the focus was on the response of the slab. The upper limit was chosen due to the increasing dynamic decoupling of the slab from the vibrations in the rail with increasing frequency. Above 1.5 kHz, this led to a low coherence in the measurements and, more practically, it implies that these frequencies are not as relevant when it comes to e.g. sound radiation from the slab or load on the foundation. Finally, as the accuracy of the hammer impact position on the rail head to some extent determines which rail modes are being excited, it is advantageous to be able to specify the input location across the rail head in the WFE model, especially at higher frequencies. However, this feature was not used in the calibration due to the inherent uncertainty of the impact position since the hammer was manually guided by a person.  The main reason to reduce the computational cost of the models is the employed Genetic Algorithm (GA). For each iteration in the GA, the receptances need to be evaluated for each member in the population. The computational demands of calculating the receptances depend on how many excitation positions and frequency lines that are considered. Here, 138 frequency lines were considered in the calibration, with four excitation positions in the FE model and eight excitation positions in the WFE model. With these settings and one CPU-core, the time to generate updated system matrices for the FE model and calculate the receptances was about 50 min. The corresponding time for the WFE model was only around one minute. Note that if it is relevant to only tune the rail pad stiffness in the calibration, the receptances of both the rail and the remaining parts of the track in the WFE model can be pre-calculated and the calibration can be reduced to a few seconds. When using the FE model, the computational cost is reduced by using 20 CPU-cores in parallel, see Fig. 10, which means that the average time to generate new system matrices and calculate the receptances for one member in one population is around 2.5 min. Finally, this means that running the GA for It is apparent that there exists a trade-off between model complexity (and thus, accuracy) and numerical efficiency. In this paper, the presented models were consciously chosen in order to provide an appropriate compromise. However, improvements could be made to both models to either increase their numerical efficiency or their accuracy. The FE model can be extended to include features such as lateral dynamics and a solid mesh for the SCC layer, see Ref. [21]. Furthermore, a rather small set of parameters has been included in the GA. In a prestudy to this paper, it was (as expected) found that when including a larger set of track parameters in the GA, the model converges towards a slightly closer match than presented here. Another factor that could benefit both models in terms of accuracy is to introduce a more elaborate model of the rail pad, for example by implementing a frequencydependent stiffness and damping.
Regarding the discrepancies when comparing the models with the measurements, there are several uncertainties about the test rig and the measurement that can be discussed. First of all, the boundary conditions between the panel slabs are unknown. In addition, the influence of the floor on both sides of the support layer has been neglected in the models as there was no further information about its properties. The influence of the floor might be especially relevant for the lateral dynamics studied in the WFE model. As mentioned previously, a large uncertainty is the discontinuous rails, leading to uncertain boundary conditions and additional modes due to reflections at the rail ends. This problem was elaborated on in the study performed by Cox et al. [36], where boxes of sand were used to reduce the boundary effects. Further, the impact position of the hammer influences the excitation of the rail. In this case, there is a possibility for inaccuracies in the order of 1 cm radius around the desired impact point. This inaccuracy is also relevant in terms of the impact angle, which might deviate from the desired purely vertical or lateral direction.
Overall, both calibrated models capture the trend of the measurements for multiple excitation positions and sensor locations within a rather small margin compared to the overall dynamic range. This implies that both models can successfully represent the dynamic response of the considered slab track.

Conclusion
In this paper, two models of the dynamic response of a slab track have been calibrated and validated using SIMO measurements in a full-  scale test rig. The measurements consisted of hammer impact measurements from which multiple cross-receptances were extracted. The calibration was divided into two steps. In the first step, the stiffness of the rail pad was tuned based on a parameter study, and in the second step, the damping of the rail pad and the stiffness and damping of the soil were calibrated using a genetic algorithm.
Both models capture the trend of the SIMO measurements for multiple excitation positions and sensor locations within a rather small margin compared to the overall dynamic range. This implies that both models can successfully represent the dynamic response of the test rig. Regarding the differences between the simulations and measurements, there are several uncertainties. Concerning the measurements, this includes the boundary conditions at the fishplates and between the slab panels, the accuracy of the excitation position and the influence of the adjacent structure on either side of the track. In addition, there are simplifications in the models that need to be considered including the assumption of using non-interacting springs and viscous dampers for several layers of the track in the finite element model and the assumption that the rail and remaining parts of the track are infinite waveguides in the waveguide finite element model. By using the calibrated and validated models obtained in this paper, a range of investigations and studies can be conducted. In particular, the finite element model will be used to assess the risk of crack initiation in the slab panel, while the waveguide finite element model will be used to model the sound radiation from slab track.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.