FOURIER-BASED ANALYTICAL FRAMEWORK
FOR NON-HOMOGENEOUS BALLAST-LESS
RAILWAY TRACKS

David Gonzalez Velasco, Denis Duhamel, Tien Hoang, Gilles Foret, Hervé Lenglin, Bertrand Findinier

► To cite this version:

David Gonzalez Velasco, Denis Duhamel, Tien Hoang, Gilles Foret, Hervé Lenglin, et al. FOURIER-BASED ANALYTICAL FRAMEWORK FOR NON-HOMOGENEOUS BALLAST-LESS RAILWAY TRACKS. COMPDYN 2025, Jun 2025, Rhodes (Grèce), Greece. hal-05459872

HAL Id: hal-05459872

https://hal.science/hal-05459872v1

Submitted on 15 Jan 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
FOURIER-BASED ANALYTICAL FRAMEWORK FOR
NON-HOMOGENEOUS BALLAST-LESS RAILWAY TRACKS

David Gonzalez Velasco<fnref n="1" />*, Denis Duhamel<fnref n="1" />, Tien Hoang<fnref n="1" />, Gilles Foret<fnref n="1" />, Hervé Lenglin<fnref n="2" />,
Bertrand Findinier<fnref n="2" />

e-mail: {david.gonzalez-velasco,denis.duhamel,tien.hoang,gilles.foret}@enpc.fr

BP69, 62904 Coquelles Cedex, France

e-mail: {herve.lenglin, bertrand.findinier}@eurotunnel.com

**Keywords:** Railway Tracks Behaviour, Dynamic Loads, Analytical Model, Periodic Structures.

**Abstract.** Understanding and predicting the behavior of structures under operational loads has long been a focal point for engineers and researchers. This is crucial to ensuring safety, meeting regulatory standards, and optimizing structural performance. Regarding railway tracks modeling, both numerical and analytical approaches have been extensively utilized, with analytical methods holding particular appeal due to their ability to capture global structural responses while requiring low computational effort.

In the context of railway tracks, typically seen as semi-infinite structures, periodicity assumptions are often adopted to reduce the problem to a finite system of equations. Existing analytical models generally consider static loading, treating the applied forces as constant and moving. However, dynamic loads have been shown to produce significantly more critical effects on structural behavior.

This work introduces an alternative analytical model capable of capturing the response of ballast-less railway tracks under both static and dynamic periodic loads. The proposed loading function can adopt various forms, provided it remains L-periodic, allowing for the study of a wide range of scenarios. The track structure is modeled as periodic with a defined spatial periodicity L. This L-periodic section is defined as a railway track segment constituted by a finite number of discrete and independent supports, allowing the model to account for non-homogeneous regions and local variations in structural properties. This modelling choice allows to converge into a solution by using Floquet's theorem and representing the problem's governing equations as Fourier series.

The proposed analytical framework highlights the amplified impact of dynamic loads on railway tracks. It enables the quantification of overloads induced by dynamic forces, laying the groundwork for developing strategies to mitigate their effects and reduce their severity. This approach not only enhances our understanding of dynamic railway track behavior but also provides insights for designing more resilient and efficient infrastructure.

<footnote n="1">¹Laboratoire Navier, UMR 8205, ENPC, Institut Polytechnique de Paris, Université Gustave Eiffel, CNRS 6 et 8 Av. Blaise Pascal, 77420 Champs-sur-Marne, France</footnote>

<footnote n="2">² Getlink</footnote>
# 1 INTRODUCTION

Modeling and predicting the behavior of structures is of great interest in the industrial world for designing systems that comply with government regulations to ensure safe operation. Additionally, such modeling ensures that functionality and design ultimately result in an optimal structure. The modeling of railway tracks falls into this category of structures warranting analysis. For those interested in modeling railway tracks, various approaches can be adopted, with numerical and analytical methods being the most prominent.

In what regards this paper, analytical approaches for studying beam-like structures began with Mead's work [7], where the characterization of wave propagation to infinity in periodically supported systems was studied, emphasizing the free propagation of undamped waves without any applied loads. Later, in [8], receptance functions were employed to derive an equation for wave propagation in one-dimensional infinite structures, such as Timoshenko beams, under periodic conditions. The loads were modeled as harmonic functions at specific points. In subsequent works [9], this topic was further explored, focusing on the flexural waves generated by applied excitation and reaction forces on supports, modeled as free-propagation waves.

Authors in [1] revisited the studies in [2, 10], proposing an analytical model of an Euler-Bernoulli beam periodically supported and excited by a harmonic force moving at a constant speed. The structure considered is continuous and homogeneous, and the problem is solved in the frequency domain using the Fourier transform. Similarly, [6, 11] addressed the problem to study phenomena such as damping in supports and critical speeds.

The problem where rails are represented as Euler-Bernoulli beams periodically supported over a three-dimensional viscoelastic layer was tackled in [12], simulating the interaction between rails, supports, and the ground. The dynamic response spectrum of the structure was found to be more significant at low frequencies. Analytical expressions were developed to calculate vertical displacement caused by moving loads, studying the effects of critical speeds, damping, and periodicity on the track's response.

In a similar analytical vein, [5] addressed the problem of an Euler-Bernoulli beam periodically supported with point loads moving at constant speed. This model introduced discrete supports and a foundation exhibiting linear or nonlinear behavior. The problem was treated in the Fourier domain, and iterative methods employing harmonic balance techniques were used to address the nonlinear behavior in supports.

The comparative of the dynamic response of railway tracks using two beam models, Euler-Bernoulli and Timoshenko, periodically supported, was done in [3]. The proposed model characterizes non-homogeneous sections of the structure, introducing the concept of periodic motifs that reduce the calculation of an infinite system of equations to a finite one. This was achieved using Floquet's theorem, subsequently solving the system of equations in the Fourier domain.

This work proposes an extension to the analytical models developed in [5] and [3]. Here, the railway track is represented as a periodically supported beam that can be segmented into distinct patterns to capture its non-homogeneous nature. In our approach, the train-generated load is defined by a spatially dependent function that satisfies specific periodic criteria. By incorporating this train-generated loading function, our model accommodates a spectrum of load types, from static to dynamic, thus enabling a more comprehensive evaluation of the dynamic effects on the track structure. Ultimately, this framework provides a versatile tool for estimating the impact of dynamic loads, offering new insights for the design and optimization of railway infrastructure.
# 2 RAILWAY TRACK ANALYTICAL MODEL DEVELOPMENT

## 2.1 MODELING CHOICES

In this work, our objective is to model the behavior of a ballastless railway track from an analytical perspective, treating the structure as a one-dimensional (1D) problem. As illustrated in Figure 1, the track is divided into repeating, identical patterns (or segments). Each pattern contains a finite number of point supports, which may or may not be spaced uniformly. The mechanical characteristics of these supports can also vary, allowing for the study of both homogeneous and non-homogeneous sections of track.

<figure><img src="image_2.png" /><figcaption>Figure 1: Periodic non-homogeneous railway track configuration subjected to train generated moving loads.</figcaption></figure>

The overall modeling approach separates the structure into three main components, the supports, the rail and the loads.

In what regards the supports, these are as independent, point-like systems that may exhibit linear or nonlinear behavior. Here, a specific common type of ballastless track support is considered, the Low Vibration Track (LVT) system. This can be idealized by a mass connected to two spring-damper arrangements (as shown in figure 2), in the lower and upper part of the mass. In such a system the rail is connected to a concrete block via an under-rail pad. The block itself (with mass *m*) rests on a second elastic stage (composed by an under-sleeper and a rubber boot), which in turn connects to the slab.

<figure><img src="image_3.png" /><figcaption>Figure 2: Low Vibration Track (LVT) support modelling choice for ballastless railway tracks with upper (1. Under-rail pad) and lower (2. Under-sleeper pad and rubber boot) elastic stages [4].</figcaption></figure>
The under-rail pad and the under-sleeper pad are each represented by their own stiffness components, $k_{r,p}(\omega)$ and $k_{f,p}(\omega)$, respectively, it is often convenient to define a single equivalent support stiffness $k_{s,p}(\omega)$ by,

$$ \frac{1}{k_{s,p}(\omega)} = \frac{1}{k_{r,p}(\omega)} + \frac{1}{k_{f,p}(\omega)} \qquad (1) $$

Where $k_{s,p}(\omega) = k_1 + i\omega\xi_1$ and $k_{f,p}(\omega) = k_2 + i\omega\xi_2 - m\omega^2$. This relationship assumes that the two spring-like elements act in series. It provides a concise way to incorporate the combined effect of both pads in each support.

On the other hand the train-generated forces acting on the system are introduced via the term $Q_T(x - NL, t - \frac{NL}{v})$, which is spatially and temporally defined with a $L$ periodicity. Details on this term are going to be presented further in this paper.

Thus, the total force acting on the system in the time domain is given by,

$$ F(x, t) = \sum_{N \in \mathbb{Z}} \sum_{p=0}^{m-1} R_p(t - \frac{NL}{v})\delta(x - pl - NL) - \sum_{N \in \mathbb{Z}} Q_T(x - NL, t - \frac{NL}{v}) \qquad (2) $$

Where the first term represents the reaction forces generated by the supports, and the second term corresponds to the forces due to train passage. As depicted in figure 1, $R_p$ denotes the reaction force at the $p$-th support, $L$ is the period length, $l$ is the spacing between consecutive supports, and the Dirac delta functions localize these forces in space.

In the frequency domain, the total force can be written as,

$$ \hat{F}(x, \omega) = \sum_{N \in \mathbb{Z}} \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{-i\omega \frac{NL}{v}} \delta(x - pl - NL) - \sum_{N \in \mathbb{Z}} \hat{Q}_T(x - NL, \omega) e^{-i\omega \frac{NL}{v}} \qquad (3) $$

About the rail, it is idealized as a beam-like element in 1D form for analytical treatment. We propose to use Timoshenko's beam theory, which is well-suited for analyzing a wide range of frequencies and yields results that more closely reflect reality, particularly when the train's speed is high and the spacing between supports is small [4]. The dynamic behavior of the beam is described by the following equations,

$$ \left\{ \begin{array}{l} \rho S \frac{\partial^2 w_r(x,t)}{\partial t^2} = \kappa S G \left( \frac{\partial^2 w_r(x,t)}{\partial x^2} - \frac{\partial \phi_r(x,t)}{\partial x} \right) + F(x,t) \\ \rho I \frac{\partial^2 \phi_r(x,t)}{\partial t^2} = E I \frac{\partial^2 \phi_r(x,t)}{\partial x^2} + \kappa S G \left( \frac{\partial w_r(x,t)}{\partial x} - \phi_r \right) \end{array} \right. \qquad (4) $$

Here, $E$ represents the Young's modulus, $I$ is the moment of inertia of the beam's cross-section, $w_r(x, t)$ denotes the vertical displacement of the beam, and $\phi_r(x, t)$ the rotation of the beam's cross-section. In addition, $\rho$ is the mass density, $S$ the cross-sectional area, and $\kappa$ the shear correction factor.

## 2.2 DYNAMIC GOVERNING EQUATIONS DERIVATION

### 2.2.1 FOURIER-BASED PERIODICALLY SUPPORTED BEAM EQUATIONS

The method presented here reduces the number of unknowns by assuming periodic conditions. Because the system is periodic, the response of a support located at $x = 0$ is identical
to that of a support located at $x = L$ (they share the same equivalent stiffness, $K_q = K_{q+Nm}$ for all $N \in Z$, and experience similar loading). In fact, the response of the support at $x = L$ is simply delayed by a time interval of $\frac{L}{v}$ compared to the support at $x = 0$. This periodicity is also imposed on the applied forces so that they are spatially and temporally periodic, ensuring a consistent formulation of the problem. Accordingly, the relationships between the reaction forces at supports in different locations and the applied forces are expressed as,

$$ \begin{cases} R_{Nm+p}(t) = R_p(t - \frac{NL}{v}), & \forall N \in Z \\ Q_T(x, t) = Q_T(x - NL, t - \frac{NL}{v}), & \forall N \in Z \end{cases} \qquad (5) $$

Here, $R_p(t)$, defined for $0 \le p < m$ at the coordinate $x = pl$, represents the reaction force of the p-th support.

In our analysis, we adopt the approach used by Hoang and Claudet, which relies on decomposing the periodic functions that describe the phenomenon into Fourier series. To solve the system of equations given in equation 4, Hoang et al. [5] propose applying the Fourier transform with respect to time $t$, followed by a Fourier series expansion with respect to the spatial variable $x$. The Fourier transform of Timoshenko's beam dynamic equations with respect to $t$ yields,

$$ \left\{ \begin{aligned} \kappa S G \frac{\partial \hat{\phi}_r(x, \omega)}{\partial x} &= \kappa S G \frac{\partial^2 \hat{w}_r(x, \omega)}{\partial x^2} + \rho S \omega^2 \hat{w}_r(x, \omega) + \hat{F}(x, \omega) \\ -\kappa S G \frac{\partial \hat{w}_r(x, \omega)}{\partial x} &= E I \frac{\partial^2 \hat{\phi}_r(x, \omega)}{\partial x^2} + (\rho I \omega^2 - \kappa S G) \hat{\phi}_r(x, \omega) \end{aligned} \right. \qquad (6) $$

Similarly, applying the same procedure to equation 2 and multiplying by a factor $e^{i\omega\frac{x}{v}}$, allows us to rewrite the force term in its frequency-domain representation (see equation 3) as,

$$ \hat{F}(x, \omega) e^{i\omega\frac{x}{v}} = \sum_{N \in Z} \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i\omega\frac{x-NL}{v}} \delta(x - pl - NL) - \sum_{N \in Z} \hat{Q}_T(x - NL, \omega) e^{i\omega\frac{x-NL}{v}} \qquad (7) $$

Observe that the product $\hat{F}(x, \omega) e^{i\omega\frac{x}{v}}$ is periodic with period $L$. Exploiting this periodicity, we introduce the new functions,

$$ \hat{w}_r(x, \omega) = \Psi(x, \omega) e^{-i\frac{\omega}{v}x} \quad \text{and} \quad \hat{\phi}_r(x, \omega) = \Phi(x, \omega) e^{-i\frac{\omega}{v}x} \qquad (8) $$

This change of variables allows us to recast equation 6 into a more tractable form. The modified equations yields,

$$ \kappa S G \left( \frac{\partial \Phi}{\partial x} - \frac{i\omega}{v} \Phi \right) = \kappa S G \left( \frac{\partial^2 \Psi}{\partial x^2} - 2 \frac{i\omega}{v} \frac{\partial \Psi}{\partial x} - \frac{\omega^2}{v^2} \Psi \right) + \rho S \omega^2 \Psi + e^{i\frac{\omega}{v}x} \hat{F}(x, \omega) \qquad (9) $$

$$ \kappa S G \left( \frac{i\omega}{v} \Psi - \frac{\partial \Psi}{\partial x} \right) = EI \left( \frac{\partial^2 \Phi}{\partial x^2} - 2 \frac{i\omega}{v} \frac{\partial \Phi}{\partial x} - \frac{\omega^2}{v^2} \Phi \right) - (\kappa S G - \rho I \omega^2) \Phi $$

By expanding $\Psi(x, \omega)$ and $\Phi(x, \omega)$ into Fourier series, we can solve the system in equation (9) by applying Floquet's theorem. The Fourier series expansions are given by,

$$ \Psi(x, \omega) = \sum_{n=-\infty}^{\infty} \mathbb{P}_n(\omega) e^{i2\pi n \frac{x}{L}} \quad \text{and} \quad \Phi(x, \omega) = \sum_{n=-\infty}^{\infty} \mathbb{Q}_n(\omega) e^{i2\pi n \frac{x}{L}} \qquad (10) $$
Where $\mathbb{P}_n(\omega)$ and $\mathbb{Q}_n(\omega)$ are the Fourier coefficients corresponding to $\Psi(x, \omega)$ and $\Phi(x, \omega)$, respectively.

Following a similar approach, we now consider the force function in the governing equations. To analyze it, we express the force function as a Fourier series. The total force function given in equation (3) comprises both the reaction forces and the train-generated forces. In terms of representing the reaction force contribution via a Fourier series, we need to extract its Fourier coefficients from the L-periodic components present in equation (7). In general, the Fourier coefficients for any L-periodic function $f(x)$ are defined as,

$$c_n = \frac{1}{L} \int_{0}^{L} f(x) e^{-i2\pi n \frac{x}{L}} dx \quad (11)$$

Using this definition to obtain the Fourier coefficient $\mathbb{R}_n(\omega)$ for the reaction forces, we have,

$$\begin{aligned} \frac{\hat{\mathbb{R}}_n(\omega)}{L} &= \frac{1}{L} \sum_{p=0}^{m-1} \int_0^L \hat{R}_p(\omega) e^{i\omega \frac{x-NL}{v}} \delta(x - x_p - NL) e^{-i2\pi n \frac{x}{L}} dx \\ &= \frac{1}{L} \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i\left(\frac{\omega}{v} - \frac{2\pi n}{L}\right)x_p} \end{aligned} \quad (12)$$

A separate procedure is employed for the train-generated force term, which will be developed and explained in detail later. For the moment, we denote the result of that Fourier expansion by $\mathbb{Q}_T(x, \omega)$.

From this point, we substitute equations (10) and (12) into equation (9) and expand the derivatives with respect to $x$. As a result, equation (9) becomes,

$$\begin{aligned} \kappa SG \left(\frac{i2\pi n}{L} - \frac{i\omega}{v}\right) \mathbb{Q}_n(\omega) &= -\kappa SG \left(\frac{\omega}{v} - \frac{2\pi n}{L}\right)^2 \mathbb{P}_n(\omega) + \rho S \omega^2 \mathbb{P}_n(\omega) + \frac{1}{L} \hat{\mathbb{R}}_n(\omega) - \mathbb{Q}_{T,n}(\omega) \\ \kappa SG \left(\frac{i\omega}{v} - \frac{i2\pi n}{L}\right) \mathbb{P}_n(\omega) &= -EI \left(\frac{\omega}{v} - \frac{2\pi n}{L}\right)^2 \mathbb{Q}_n(\omega) - (\kappa SG - \rho I \omega^2) \mathbb{Q}_n(\omega) \end{aligned} \quad (13)$$

By isolating $\mathbb{Q}_n(\omega)$ from the second equation of system (13) and substituting it back into the first, we obtain the Fourier coefficients for $\Phi(x, \omega)$ and $\Psi(x, \omega)$ as

$$\begin{aligned} \mathbb{P}_n(\omega) &= \tilde{\mathbb{P}}_n(\omega) \hat{\mathbb{R}}_n(\omega) - L \tilde{\mathbb{P}}_n(\omega) \hat{\mathbb{Q}}_{T,n}(\omega) \\ \mathbb{Q}_n(\omega) &= \tilde{\mathbb{Q}}_n(\omega) \hat{\mathbb{R}}_n(\omega) - L \tilde{\mathbb{Q}}_n(\omega) \hat{\mathbb{Q}}_{T,n}(\omega) \end{aligned} \quad (14)$$

The quantities $\tilde{\mathbb{P}}_n(\omega)$, $\tilde{\mathbb{Q}}_n(\omega)$, and $z_n(\omega)$ are computed according to the following formulas [5],

$$\left\{ \begin{aligned} \tilde{\mathbb{P}}_n(\omega) &= \frac{\kappa SG - \rho I \omega^2 + EI \left(\frac{\omega}{v} - \frac{2\pi n}{L}\right)^2}{L z_n(\omega)} \\ \tilde{\mathbb{Q}}_n(\omega) &= -\frac{i \kappa SG \left(\frac{\omega}{v} - \frac{2\pi n}{L}\right)}{L z_n(\omega)} \\ z_n(\omega) &= \kappa SG EI \left(\frac{\omega}{v} - \frac{2\pi n}{L}\right)^4 - \rho I S \omega^2 (\kappa G + E) \left(\frac{\omega}{v} - \frac{2\pi n}{L}\right)^2 + \rho S \omega^2 (\rho I \omega^2 - \kappa SG) \end{aligned} \right. \quad (15)$$

Finally, substituting equation (14) into equation (10) yields the final expressions for $\hat{w}_r(x, \omega)$ and $\hat{\phi}_r(x, \omega)$:
$$ 
\begin{aligned}
\Psi(x, \omega) &= \sum_{n=-\infty}^{\infty} \hat{\mathbb{R}}_n(\omega) \tilde{\mathbb{P}}_n(\omega) e^{i2\pi n \frac{x}{L}} - L \sum_{n=-\infty}^{\infty} \hat{\mathbb{Q}}_{T,n}(\omega) \tilde{\mathbb{P}}_n(\omega) e^{i2\pi n \frac{x}{L}} \\
\Phi(x, \omega) &= \sum_{n=-\infty}^{\infty} \hat{\mathbb{R}}_n(\omega) \tilde{\mathbb{Q}}_n(\omega) e^{i2\pi n \frac{x}{L}} - L \sum_{n=-\infty}^{\infty} \hat{\mathbb{Q}}_{T,n}(\omega) \tilde{\mathbb{Q}}_n(\omega) e^{i2\pi n \frac{x}{L}}
\end{aligned} 
\quad (16) $$

### 2.2.2 TRAIN-GENERATED FORCES EQUATIONS

Regarding the train-generated forces acting on the system, these are incorporated through the term $Q_T(x - NL, t - \frac{NL}{v})$ defined earlier. In accordance with the periodicity conditions established in equations 5, this term is given by,

$$ Q_T \left( x - NL, t - \frac{NL}{v} \right) = \sum_{j=1}^{K} Q_j(x - NL) \delta \left( x + D_j - vt \right) \quad (17) $$

Here, the function $Q_j(x - NL)$ enables the introduction of various forces as long as they adhere to the periodicity constraints. In essence, $Q_j(x - NL)$ is equivalent to $Q_j(x)$ defined on the interval $[0, L]$ and then extended periodically to capture the contribution of the forces generated by the train's passage. Typically, the forces on railway structures can be categorized into two types, static loads that represent the train's weight, and dynamic loads that travel along the track at a constant speed $v$. Each dynamic load is separated from the next by a distance $D_j$, corresponding to the spacing between consecutive train wheels.

In the frequency domain, equation 17 is transformed into,

$$ \hat{Q}_T(x - NL, \omega) = \frac{1}{v} \sum_{j=1}^{K} Q_j(x - NL) e^{-i\omega \frac{(x+D_j)}{v}} \quad (18) $$

For simplicity in the expressions that follow, we introduce the notation,

$$ \hat{Q}(x, \omega) = \sum_{N \in \mathbb{Z}} \sum_{j=1}^{K} Q_j(x - NL) e^{-i\omega \frac{D_j}{v}} \quad (19) $$

From this point onward, we focus on determining the Fourier coefficient for $Q_T(x - NL, \omega)$ using the same approach applied to $\Psi(x, \omega)$, $\Phi(x, \omega)$, and the support reaction forces. In particular, we employ the Fourier coefficient formula given in equation (11). Thus, the Fourier coefficient for the train-generated forces is obtained as,

$$ \hat{Q}_{T,n}(\omega) = \frac{1}{Lv} \int_{0}^{L} \hat{Q}(x, \omega) e^{-i2\pi n \frac{x}{L}} dx \quad (20) $$

### 2.2.3 RAIL DISPLACEMENT EXPRESSION FROM THE FORCES

To express equation 16 in a more compact and intuitive form, we introduce the functions,

$$ 
\begin{aligned}
\eta(x, \omega) &= \sum_{n=-\infty}^{\infty} \tilde{\mathbb{P}}_n(\omega) e^{i \frac{2\pi n x}{L}} \\
\gamma(x, \omega) &= \sum_{n=-\infty}^{\infty} \tilde{\mathbb{Q}}_n(\omega) e^{i \frac{2\pi n x}{L}}
\end{aligned} 
\quad (21) $$
Using the first relation in equation 21 in conjunction with equation 8, we can rewrite the
vertical displacement of the beam, $\hat{w}_r(x, \omega)$, as,

$$ \hat{w}_r(x, \omega) = \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x)} \eta(x - x_p, \omega) - L \sum_{n=-\infty}^{\infty} \hat{Q}_{T,n}(\omega) \tilde{P}_n(\omega) e^{i \left( \frac{2\pi n}{L} - \frac{\omega}{v} \right) x} \quad (22) $$

To obtain a general expression for the beam's vertical displacement at any spatial point $x$ and frequency $\omega$, we now focus on the train-generated force term. Replacing the Fourier expansion of the train-generated forces (as given in equation 20) into the second term of equation 16 yields,

$$ \begin{aligned} L \sum_{n=-\infty}^{\infty} \hat{Q}_{T,n}(\omega) \tilde{P}_n(\omega) e^{i \frac{2\pi n}{L} x} &= \frac{1}{v} \sum_{n=-\infty}^{\infty} \left( \int_0^L \hat{Q}(s, \omega) e^{-i \frac{2\pi n}{L} s} ds \right) \tilde{P}_n(\omega) e^{i \frac{2\pi n}{L} x} \\ &= \frac{1}{v} \int_0^L \hat{Q}(s, \omega) \left( \sum_{n=-\infty}^{\infty} \tilde{P}_n(\omega) e^{i \frac{2\pi n}{L} (x-s)} \right) ds \quad (23) \\ &= \frac{1}{v} \int_0^L \hat{Q}(s, \omega) \eta(x-s, \omega) ds \end{aligned} $$

We then define,

$$ \tilde{Q}(x, \omega) = \frac{1}{v} \int_{0}^{L} \hat{Q}(s, \omega) \eta(x - s, \omega) ds \quad (24) $$

So that equation 22 becomes,

$$ \hat{w}_r(x, \omega) = \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x)} \eta(x - x_p, \omega) - \tilde{Q}(x, \omega) e^{-i \frac{\omega}{v} x}, \quad (25) $$

valid for all $x \in [0, L]$

Furthermore, the function $\eta(x, \omega)$ originally defined in [4], is given explicitly by

$$ \left\{ \begin{aligned} \eta(x, \omega) &= \frac{e^{i \frac{\omega x}{v}}}{2EI(\lambda_1^2 + \lambda_2^2)} \left[ \frac{C_1 \sin \lambda_1 (L-x) + e^{-i \frac{\omega L}{v}} \sin \lambda_1 x}{\lambda_1 \cos L\lambda_1 - \cos \frac{\omega L}{v}} \right. \\ & \qquad \left. - \frac{C_2 \sinh \lambda_2 (L-x) + e^{-i \frac{\omega L}{v}} \sinh \lambda_2 x}{\lambda_2 \cosh L\lambda_2 - \cos \frac{\omega L}{v}} \right] \\ C_{1,2} &= 1 - \frac{\rho I \omega^2 \mp EI \lambda_{1,2}^2}{\kappa SG} \\ \lambda_{1,2}^2 &= \sqrt{\frac{\omega^4}{4} \left( \frac{\rho I}{EI} - \frac{\rho S}{\kappa SG} \right)^2 + \frac{\rho S \omega^2}{EI}} \pm \frac{\omega^2}{2} \left( \frac{\rho I}{EI} + \frac{\rho S}{\kappa SG} \right) \end{aligned} \right. $$

Equation 25 thus represents the vertical displacement of the beam (rail) as a result of the system's inherent periodicity and the dynamics described by the Timoshenko beam theory.

The term $\tilde{Q}(x, \omega)$ encapsulates the interaction between the beam's behavior and the external loading. Since this study addresses a direct problem, the imposed function $\hat{Q}(x, \omega)$ makes equation (25) solvable. Section 3 presents diverse case studies in which specific forms of $\hat{Q}(x, \omega)$ are used to evaluate the behavior of the railway track under static and dynamic loading.
To express the rail's vertical displacement at a particular support, let $x = x_q$, where $x_q$ denotes the coordinate of support $q$. Then $\hat{w}_q(\omega) = \hat{w}_r(x_q, \omega)$. Since $\eta(x, \omega)$ is L-periodic, we define $\eta_{q,p}(\omega) = \eta(x_q - x_p, \omega) = \eta(x_q - x_{p+m}, \omega)$. Thus, equation (25) becomes,

$$ \hat{w}_q(\omega) = \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x_q)} \eta_{q,p}(\omega) - \tilde{Q}(x_q, \omega) e^{-i \frac{\omega}{v} x_q} \quad (26) $$

Expressing equation (26) in matrix form yields,

$$ \underline{C_p \hat{R}} = \underline{\tilde{Q}} + \hat{w} \quad (27) $$

Equation (27) represents the relationship between the support reaction forces and the vertical displacement of the beam. However, neither the reaction forces $\hat{R}$ nor the displacements $\hat{w}$ are known a priori. Therefore, an additional relation is needed, this is provided by the dynamic stiffness equation of the support.

## 2.3 SUPPORT'S DIRECT RESPONSES COMPUTATION

The method proposed here uses a direct formulation to find support's reaction forces and vertical displacements. This approach is based on the support stiffness equation. In matrix form, one may write,

$$ \hat{w} = -\underline{D \hat{R}} \quad (28) $$

Where the unknown vector is the reaction force vector $\hat{R}$ and $D$ is a diagonal matrix containing the inverses of the dynamic stiffnesses of the supports,

$$ \underline{D} = \begin{pmatrix} k_{s,0}^{-1} & 0 & \cdots & 0 \\ 0 & k_{s,1}^{-1} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & k_{s,m-1}^{-1} \end{pmatrix} $$

Alternatively, the support stiffness can be expressed as,

$$ \underline{\hat{R}} = -\underline{K \hat{w}} \quad (29) $$

Where,

$$ \underline{K} = \begin{pmatrix} k_{s,0} & 0 & \cdots & 0 \\ 0 & k_{s,1} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & k_{s,m-1} \end{pmatrix} $$

Substituting equation (28) into equation (27) leads to,

$$ \underline{\hat{R}} = \underline{\underline{A}^{-1} \tilde{Q}} \quad (30) $$

with $\underline{A} = \underline{\underline{C_p + D}}$. On the other hand, substituting equation (29) into equation (27) gives,

$$ \hat{w} = -\underline{\underline{B}^{-1} \tilde{Q}} \quad (31) $$
With $\underline{\underline{\underline{\underline{B}}}} = \underline{\underline{\underline{\underline{C}}}}_p \underline{\underline{\underline{\underline{\underline{K}}}}} + \underline{\underline{\underline{\underline{\underline{I}}}}}$.

At each frequency, the $m \times m$ matrices A and B are computed and inverted. Equation (30) is preferable when analyzing the impact of stiffer supports in a motif, while equation (31) is better suited for assessing the effect of supports with low (or even zero) stiffness.

### 3 EXAMPLES

In the present section, a dynamic loading function is going to be imposed. To simplify expressions and make the analysis easier to perform, we are going to treat the case where the structure is going to be loaded by the passage of only one wheel ($D_j = 0$ and $K = 1$). The structure is going to be generally composed of 40 supports, and when referring to a non-homogenous track, we consider the case where the track contains a broken sleeper at half span of the structure.

The mechanical parameters describing the Timoshenko's beam behavior are presented in table 1.

<table><caption>Table 1: Parameters used for railway track response computation.</caption><thead><tr><th>Parameter</th><th>Notation</th><th>Value</th></tr></thead><tbody><tr><td>Rail mass</td><td>ρS</td><td>60 kg m<sup>-1</sup></td></tr><tr><td>Rail flexion stiffness</td><td>EI</td><td>6.38 MN m<sup>2</sup></td></tr><tr><td>Timoshenko's ratio</td><td>κ</td><td>0.4</td></tr><tr><td>Shear modulus</td><td>G</td><td>8.077 GPa</td></tr><tr><td>Train's speed</td><td>v</td><td>30 m s<sup>-1</sup></td></tr><tr><td>Train-generated force</td><td>Q<sub>d</sub></td><td>100 kN</td></tr><tr><td>Block's mass</td><td>m</td><td>100 kg</td></tr><tr><td>Supports spacing</td><td>l</td><td>0.6 m</td></tr><tr><td>Railpad stiffness</td><td>k<sub>r</sub></td><td>192 MN m<sup>-1</sup></td></tr><tr><td>Railpad damping</td><td>ξ<sub>r</sub></td><td>1.97 MN s m<sup>-1</sup></td></tr><tr><td>Under-sleeper stiffness</td><td>k<sub>f</sub></td><td>26.4 MN m<sup>-1</sup></td></tr><tr><td>Under-sleeper damping</td><td>ξ<sub>f</sub></td><td>0.17 MN s m<sup>-1</sup></td></tr></tbody></table>

#### 3.1 STEP FUNCTION

By introducing a dynamical force of a step function kind, $\hat{Q}(x)$ takes the following form,

$$ \hat{Q}(x) = \sum_{N \in \mathbb{Z}} Q_d \mathcal{H}(x - NL) \quad (32) $$

As the step function is defined as a function taking the value 1 within the interval $[a, b]$ ($a$ and $b$ being coordinates in the spatial domain $[0, L]$) and 0 elsewhere. $\hat{Q}_{T,n}$ takes the form,

$$ \hat{Q}_{T,n} = \frac{1}{Lv} \int_{0}^{L} \sum_{N \in \mathbb{Z}} Q_d \mathcal{H}(x - NL)e^{-i2\pi n \frac{x}{L}} dx = \begin{cases} Q_d \frac{i(e^{-i2\pi n \frac{b}{L}} - e^{-i2\pi n \frac{a}{L}})}{2\pi n}, & n \neq 0 \\ Q_d \frac{(b-a)}{L}, & n = 0 \end{cases} \quad (33) $$

Then, by setting $x = x_q$ with $x_q \in [0, L]$ and using the relation defined in equation 24, the analytical expression for the beam's displacement at a given frequency and coordinate is,
$$
\begin{align}
\hat{w}_q(\omega) &= \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x_q)} \eta_{q,p}(\omega) - \tilde{Q}(x_q, \omega) e^{-i \frac{\omega}{v} x_q} \nonumber \\
&= \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x_q)} \eta_{q,p}(\omega) - \frac{Q_d}{v} e^{-i \frac{\omega}{v} x_q} \sum_{n=-\infty}^{\infty} \int_a^b \tilde{P}_n e^{i \frac{2\pi n}{L} (x_q - s)} ds \tag{34} \\
&= \sum_{p=0}^{m-1} \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x_q)} \eta_{q,p}(\omega) - \frac{Q_d}{v} e^{-i \frac{\omega}{v} x_q} \int_{x_q-b}^{x_q-a} \eta(s, \omega) ds \nonumber
\end{align}
$$

From the displacement expression above, the reaction forces are determined by applying the dynamic stiffness equation of the supports. Figure 3 shows the numerical results for two cases, one for a homogeneous track, and another for a non-homogeneous track. In the non-homogeneous case, the support located at the mid-span of the spatial domain is broken (its stiffness is set to zero). In this analysis, the step function was defined over the interval between supports $x_{q=-10}$ and $x_{q=10}$. The presented results show the peak reaction forces experienced by the supports over the time period simulated, corresponding to the passage of one wheel.

<figure><img src="image_4.png" /><figcaption>Figure 3: Support's reaction forces for a step type load on a homogeneous track (left) and a non-homogeneous track (right).</figcaption></figure>

The results highlight two key observations. First, as expected, when the step function takes the value $Q_d$, the response corresponds to the static load case. Second, the presence of non-homogeneities leads to significant overloads on the supports. Specifically, at the borders of the defined step function, the overloads on supports can be about twice as high as those in the regular area (see the left graph in figure 3), this effect is attributed to dynamic loading. Additionally, overloads are observed due to spatial non-homogeneities, as evidenced by the mid-span location in the right-hand graph of figure 3.

Finally, the step function result can be validated as well analytically. When the step function is set to 1 over the entire spatial domain $[0, L]$, the resulting expression should match that obtained on [3], who analyzed the problem for a static moving load. As shown below,
$$
\begin{align}
\hat{w}_q(\omega) &= \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x_q)} \eta_{q,p}(\omega) - \frac{Q_d}{v} e^{-i \frac{\omega}{v} x_q} \sum_{n=-\infty}^{\infty} \int_{x_q-L}^{x_q} \tilde{\mathbb{P}}_n e^{i 2\pi n \frac{s}{L}} ds \nonumber \\
&= \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x_q)} \eta_{q,p}(\omega) - \frac{Q_d}{v} e^{-i \frac{\omega}{v} x_q} \sum_{n=-\infty}^{\infty} \delta_{0n} \tilde{\mathbb{P}}_n L \tag{35} \\
&= \hat{R}_p(\omega) e^{i \frac{\omega}{v} (x_p - x_q)} \eta_{q,p}(\omega) - \frac{Q_d}{v} \tilde{\mathbb{P}}_0(\omega) L e^{-i \frac{\omega}{v} x_q} \nonumber
\end{align}
$$

3.2 SINUSOIDAL FUNCTION

Let's consider a general sinusoidal type function of the following kind, defined on the spatial interval $[0 L]$,

$$
\hat{Q}(x) = Q_d \sin \frac{\pi(x-a)}{b-a} \mathbb{1}_{[a \ b]}(x) \qquad (36)
$$

Where $Q_d$ is the load amplitude and $\mathbb{1}_{[a \ b]}(x)$ is the indicator function that restricts the function to the interval $[a \ b]$.

We will study the case in which this sinusoidal function loads the system. In this example,
the function is defined such that it takes positive values during the excitation; specifically, we
set $a = 11.395 \text{ m}$ and $b = 11.405 \text{ m}$, meaning that the load is applied over a span of $1 \text{ cm}$.

<figure><img src="image_5.png" /><figcaption>Figure 4: Rail's vertical displacement over time and reaction forces for a sinusoidal type load function.</figcaption></figure>

Figure 4 illustrates the results from the numerical implementation of the sinusoidal load case. Upon reviewing these graphs, it is evident that the dynamic loads result in significantly higher reaction forces at the supports, up to twice the magnitude to the ones observed for the static case load. While these forces do propagate to neighboring supports, their greater impact is observed at the specific location of the applied load. This highlights the potential risk of dynamic loading for the directly affected support, as it experiences substantially larger forces compared to those of the static load.
# 4 CONCLUDING REMARKS

In this work, we developed an analytical model for the dynamic response of ballastless railway tracks using periodic conditions and Timoshenko beam theory. By dividing the track into repeating segments with supports that may exhibit varying properties, we were able to reduce the complex, infinite-dimensional problem to a finite system via Fourier series and Floquet's theorem.

By incorporating a loading function defined in both time and space, our model is highly versatile. It can accommodate any desired load, provided it meets the periodic conditions, thus extending the scope of the model to simulate a wide range of realistic loading scenarios. This flexibility allows for a comprehensive assessment of the structure's behavior under various imposed loads, making the model a powerful tool for both design and maintenance of railway infrastructure.

Moreover, our numerical analyses show that dynamic loads can produce reaction forces nearly twice as high as those from static loads. In particular, the support directly under the load experiences the highest force, even though these elevated forces spread to neighboring supports. This observation highlights that local disturbances can critically increase the stress on individual supports.

# REFERENCES

[1] P. Belotserkovskiy, On the oscillations of infinite periodic beams subjected to a moving concentrated force, Journal of Sound and Vibration 193 (1996) 705-712.

[2] P. Belotserkovskiy, Parametric oscillations of a railway vehicle wheel, connected with sleepers spacing, Proceedings of the Moscow Institute for Railway Transport Engineers 866 (1992) 15-23.

[3] B. Claudet, D. Duhamel, G. Foret, T. Hoang, J.L. Pochet, F. Sabatier, "A Comparison of Beam Models for the Dynamics of Railway Tracks on a Non-Uniform Viscoelastic Foundation", ICSV27, Ostrava, Czech Republic, Jul 2021.

[4] T. Hoang, Comportement mécanique et dégradation d'une voie ferrée en tunnel, PhD Thesis Université Paris-Est (2011).

[5] T. Hoang, D. Duhamel, G. Foret, H. Yin, G. Cumunel, Response of a periodically supported beam on a nonlinear foundation subjected to moving loads, Journal of Nonlinear Dynamics 86 (2016) 953-961.

[6] A. Mallik, S. Chandra, A. B. Singh, Steady-state response of an elastically supported infinite beam to a moving load, Journal of Sound and Vibration 291 (2006) 1148-1169.

[7] D. Mead, Free wave propagation in periodically supported, infinite beams, Journal of Sound and Vibration 11 (1970) 181-197.

[8] D. Mead, A new method of analyzing wave propagation in periodic structures; applications to periodic timoshenko beams and stiffened plates, Journal of Sound and Vibration 104 (1986) 9-27.

[9] D. Mead, Y. Yaman, The response of infinite periodic beams to point harmonic forces: A flexural wave analysis, Journal of Sound and Vibration 144 (1991) 507-529.

[10] D. Mead, Vibration response and wave propagation in periodic structures, Journal of Engineering for Industry 93 (1971) 783-792.

[11] X. Sheng, M. Li, C. Jones, D. Thompson, Using the fourier-series approach to study interactions between moving wheels and a periodically supported rail, Journal of Sound and Vibration 303 (2007) 873-894.
[12] A. Vostroukhov, A. Metrikine, Periodically supported beam on a visco-elastic layer as a model for dynamic analysis of a high-speed railway track, International Journal of Solids and Structures 40 (2003) 5723–5752.