Design method for suspension parameters of high-speed railway vehicles considering uncertainties

Abstract

A design method for suspension parameters is proposed to ensure the desired performance of a railway vehicle despite uncertainties. The method used in this study consists of two loops: a sub-loop, which uses parameter uncertainties to determine the minimum value of target performance, and a main loop, which optimises suspension parameters such that the minimal performance achieves Pareto optimality. Thus, by selecting a suitable Pareto solution that ensures the minimum performance is equal to or higher than the required performance, the robustness of the designed vehicle against the parameter uncertainties can be achieved. The lowest low-frequency damping ratio in the service speed range was selected as one of the design targets. The other target was the critical speed. Two wheel-rail contact characteristics were applied, considering the worst conditions for each design target. The distribution of the parameter uncertainties that minimised the vehicle performance showed distinct characteristics between the longitudinal and lateral suspension parameters, which provides useful information for the maintenance of railway vehicles. The performance analysis for 100,000 random samples shows that the vehicle performed better than the targets with a high probability despite the parameter uncertainties. Further, the simulation results using multi-body dynamics simulation software validated the viability of the proposed method.

Introduction

In railway vehicles, suspension systems significantly affect running performance, such as stability, safety and ride comfort. As a result, suspension design becomes more important in high-speed railway vehicles to maintain high safety standards and quality service even at high speeds. Furthermore, the related performances of railway vehicles often conflict, deriving an adequate balance between them is crucial.

Optimisation techniques effectively solve problems with multiple input variables and one or several output variables. As a result, studies have been conducted on the use of optimisation in the design of suspension parameters of railway vehicles. He and McPhee ¹ introduced a two-loop method using a genetic algorithm and sequential quadratic programming to optimise the lateral stability of railway vehicles. In addition, they proposed a method for multibody systems using a hybrid multidisciplinary optimisation method and reported the trade-off relationship between lateral stability and curving performance. ² Johnsson et al. ³ introduced the Pareto optimisation method to enhance the safety and comfort of high-speed railway vehicles. The objective functions of safety and comfort were calculated from the response of a three-car model for lateral track irregularities and revealed a trade-off between safety and comfort based on the suspension parameters. The sensitivity of the optimal solution was also investigated for train speed, worn wheel profile, load and friction coefficients. It was found that speed and wornness significantly influenced safety and comfort. They further expanded their research to include the rate of wear of the wheel and rail and suggested some adaptive control strategies based on the optimisation results for the vehicle speed. ⁴ Mousavi Bideleh et al. ⁵ presented a multistep method for improving the computational efficiency of Pareto optimisation. Five scenarios with different track radii and vehicle speeds were applied to calculate the objective functions, and the advantages and disadvantages of the asymmetric suspension system were reported. Yao et al. ⁶ proposed optimal suspension parameters for hunting stability robustness, defined as five indices (linear stability indices for three vehicle speeds, standard deviations for suspension parameters and equivalent conicity). The K-Means clustering method was adopted to categorise Pareto fronts into six groups to meet the multiple objectives and matching rules for key suspension properties were suggested upon analysing each group characteristics. Chen et al. ⁷ suggested a multi-objective optimisation of high-speed train suspension parameters for hunting stability. The multi-objective function comprised stability indices for low and high wheel-rail contact conicity, and the robustness index defined as the standard deviation of the stability indices with the equivalent conicity distributed in 0.05–0.3. They also investigated the stability and wear Pareto optimisation of suspension parameters for a high-speed passenger train. ⁸ Stability indices were defined as the lateral ride comfort at low equivalent conicity and the bogie frame acceleration at high conicity. The contact spreading index and wear number were proposed for wheel wear performance. Li et al. ⁹ studied suspension parameter matching of high-speed locomotives for stability and comfort Pareto optimisation. They proposed three types of parameter combinations and investigated the performance characteristics for the vehicle speed and rail cant variation. Recently, to solve conflicting problems in vehicle performance, optimal suspension design methods with variable passive dampers such as the frequency-dependent yaw damper, have been introduced.^10,11 In addition, active suspension systems^12–14 have been developed for decades to improve vehicle performance under various environments. However, only a few systems were commercialised for railway vehicles.

Most components used in suspension systems differ from the designed values due to manufacturing tolerances, service periods and external environments. Uncertainties in design parameters can cause deviations in the optimal performance, and furthermore, in severe cases, make original optimal performance unsatisfactory to their requirements or constraints. Hence, the robustness of the optimal design against uncertainties in design parameters is important, but a related study has rarely been introduced. The robustness can be secured in two main ways. The first method is evaluating the robustness of the optimal results using sensitivity analysis after the design process is complete. Mousavi Bideleh ¹⁵ investigated the robustness of Pareto optimal values for the uncertainty of suspension parameters modelled as a lognormal distribution. The maximum entropy and M-DRM concept improved the computational efficiency with fairly appropriate accuracy. Gao et al. ¹⁶ proposed the Hermit orthogonal polynomial expansion method to analyse the sensitivity of the critical speed with uncertain parameters. Zhou et al. ¹⁷ carried out the global sensitivity analysis of the stochastic response of the railway vehicle with uncertain parameters based on the polynomial chaos expansion method and Sobol indices. Park ¹⁸ investigated the influence of the tolerance of suspension parameters on the optimal performance through parametric studies based on yaw and lateral stiffness and damping of the bogie frame. The disadvantage of the first method is that it requires large-scale additional computational work after the optimisation is completed. In addition, if the optimised parameters cannot meet the desired robustness, the optimisation or robustness evaluation should be repeated by changing their constraint conditions. Simple methods,^19,20 ignoring the interaction between suspension parameters and changing only one parameter at a time, can be applied with a low computational cost; however, they are less quantitatively accurate in evaluating the robustness. The other method is to apply robustness indicators to the objective functions of the optimisation process. It has the advantage that no additional work is required to evaluate the robustness. However, few studies have been conducted, and they have disadvantages. Martowicz et al. ²¹ proposed a response surface method to describe the relationship between the target performance and an uncertain parameter. Many simulations to acquire the surrogate model should precede before the optimal design. Standard deviation^6,7 of target performance for parameter uncertainties needs the distribution model of parameter uncertainties, and the modelling accuracy for uncertainty distribution affects the optimisation results. Thus, a design method that is more computationally efficient and always ensures the robustness for parameter uncertainties is required.

This study proposes a new design method, considering the uncertainty of suspension parameters. The proposed method requires no additional work and modelling or assumptions regarding the distribution of parameter uncertainties apart from the range of uncertainties. It consists of a sub-loop to find the minimal performance as the parameter uncertainties vary and a main loop to optimise the suspension parameters until the minimal performance reaches Pareto optimality, that is a condition of optimal trade-off among conflicting performance indicators. Thus, the railway vehicle always performs better than required, despite the parameter uncertainty, if the proper suspension parameters can be selected among the Pareto solutions.

The validity of the proposed method was verified through the design process and the results of suspension parameters for two target performances of a linearised high-speed railway vehicle model. The distribution of parameter uncertainties that minimises vehicle performance was examined using Pareto solutions. The Pareto set that satisfied the minimum required performance of the railway vehicle was chosen as the final optimal suspension parameters, and for 100,000 random samples assumed as their uncertainties, the variation characteristics of the vehicle performances were investigated. Finally, the performances of the reference vehicle model and the designed vehicle model were compared through simulation results using multi-body dynamics simulation software.

Design procedure

Linearised vehicle model

A 17-degrees-of-freedom railway vehicle model was used to design the suspension parameters and analyse the vehicle characteristics. This model was established in a previous study^18,22 and was slightly modified in the current study: the parameter of primary suspension was separated into a rod bush ( $k_{\mathrm{rx}},k_{\mathrm{ry}}$ ) and a coil spring ( $k_{\mathrm{cx}},k_{\mathrm{cy}}$ ); however, only the coil spring ( $k_{\mathrm{cz}}$ ) was applied for the vertical stiffness. In addition, it was assumed that the lateral dampers and the yaw dampers of the secondary suspension follow the Maxwell damper model. Other dampers were considered typical viscous dampers because the cutting frequencies of the dampers ²³ were sufficiently high compared to the modal frequencies of each component. The wheel-rail contact forces were calculated based on Kalker’s linear theory. Figure 1 shows the vehicle model comprising a carbody, two bogie frames and four wheelsets. The carbody and bogie frames were allowed to move in the lateral, roll and yaw directions, and each wheelset had lateral and yaw motions. Thus, the degrees-of-freedom of the vehicle model is 17. The equations of motion of the vehicle model on a straight track can be expressed as shown in equations (1)–(10):

$$\begin{array}{c}\multicolumn{1}{c}{{{\stackrel{··}{y}}_w}_j+2(k_{\mathrm{cy}}+k_{\mathrm{ry}})({y_w}_j-{y_b}_i)\mp 2(l_{\mathrm{wb}}{k}_{\mathrm{cy}}+(l_{\mathrm{wb}}-l_{\mathrm{kx}})k_{\mathrm{ry}}){\psi }_{b_i}}\\ -2l_{\mathrm{wz}}(k_{\mathrm{cy}}+k_{\mathrm{ry}}){{\varphi }_b}_i=F_j(V,{y_w}_j,{{\stackrel{·}{y}}_w}_j,{{\psi }_w}_j{{\stackrel{·}{\psi }}_w}_j),\hfill \end{array}$$

(1)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathrm{I}}_{\mathrm{wz}}{{\stackrel{··}{\psi }}_w}_j+I_{\mathrm{wy}}\left(\frac{r_0}{V}\right)\left(\frac{\lambda }{a_1}\right){{\stackrel{·}{y}}_w}_j+2l_{\mathrm{wa}}^2(k_{\mathrm{cx}}+k_{\mathrm{rx}})({{\psi }_w}_j-{{\psi }_b}_i)}\\ =M_j(V,{y_w}_j,{{\stackrel{·}{y}}_w}_j,{{\psi }_w}_j{{\stackrel{·}{\psi }}_w}_j)\hfill \end{array}$$

(2)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathrm{m}}_{\mathrm{b}}{{\stackrel{··}{y}}_b}_i+2c_{\mathrm{sy}}{{\stackrel{·}{x}}_{\mathrm{py}}}_i+2(k_{\mathrm{cy}}+k_{\mathrm{ry}})(2{y_b}_i+2l_{\mathrm{wz}}{{\varphi }_b}_i-{y_w}_j-{y_w}_{(j+1)})}\\ +2k_{\mathrm{sy}}({y_b}_i-y_v-l_{\mathrm{bz}}{{\varphi }_b}_i-2l_{\mathrm{cz}}{\varphi }_v\mp l_{\mathrm{bb}}{\psi }_v)=0,\hfill \end{array}$$

(3)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathrm{I}}_{\mathrm{bz}}{{\stackrel{··}{\psi }}_b}_i+(4l_{\mathrm{wb}}^2{k}_{\mathrm{cy}}+4{(l_{\mathrm{wb}}-l_{\mathrm{kx}})}^2{k}_{\mathrm{ry}}+4l_{\mathrm{wa}}^2(k_{\mathrm{cx}}+k_{\mathrm{rx}})){{\psi }_b}_i}\\ +2l_{\mathrm{ba}}^2{k}_{\mathrm{sx}}({{\psi }_b}_i-{\psi }_v)+2l_{\mathrm{bcsxy}}{c}_{\mathrm{sx}}{{\stackrel{·}{x}}_{\mathrm{px}}}_i\hfill \\ +2l_{\mathrm{wb}}{k}_{\mathrm{py}}({y_w}_{(j+1)}-{y_w}_j)-2l_{\mathrm{wa}}^2{k}_{\mathrm{px}}({{\psi }_w}_j+{{\psi }_w}_{(j+i)})\hfill \\ +2l_{\mathrm{bcsyx}}{c}_{\mathrm{sy}}{{\stackrel{·}{x}}_{\mathrm{py}}}_i=0,\hfill \end{array}$$

(4)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathrm{I}}_{\mathrm{bx}}{{\stackrel{··}{\varphi }}_b}_i+4l_{\mathrm{wa}}^2{k}_{\mathrm{pz}}{{\varphi }_b}_i+4l_{\mathrm{wa}}^2{c}_{\mathrm{pz}}{{\stackrel{·}{\varphi }}_b}_i-2l_{\mathrm{wz}}(k_{\mathrm{cy}}+k_{\mathrm{ry}})({y_w}_j+{y_w}_{(j+1)})}\\ +4l_{\mathrm{wz}}(k_{\mathrm{cy}}+k_{\mathrm{ry}}){y_b}_i+4l_{\mathrm{wz}}^2(k_{\mathrm{cy}}+k_{\mathrm{ry}}){{\varphi }_b}_i\hfill \\ +2l_{\mathrm{ba}}^2{k}_{\mathrm{sz}}({{\varphi }_b}_i-{\varphi }_v)+2l_{\mathrm{cszy}}^2{c}_{\mathrm{sz}}({{\stackrel{·}{\varphi }}_b}_i-{\stackrel{·}{\varphi }}_v)\hfill \\ -2l_{\mathrm{bz}}{k}_{\mathrm{sy}}({y_b}_i-y_v)-2l_{\mathrm{bcsyz}}{c}_{\mathrm{sy}}{{\stackrel{·}{x}}_{\mathrm{py}}}_i+2l_{\mathrm{cz}}{l}_{\mathrm{bz}}{k}_{\mathrm{sy}}{\varphi }_v\hfill \\ +2l_{\mathrm{ccsyz}}{l}_{\mathrm{bcsyz}}{\stackrel{·}{\varphi }}_v+2l_{\mathrm{bz}}^2{k}_{\mathrm{sy}}{{\varphi }_b}_i+2l_{\mathrm{bcsyz}}^2{{\stackrel{·}{\varphi }}_b}_i\hfill \\ +k_{\mathrm{sroll}}({{\varphi }_b}_i-{\varphi }_v)\pm 2l_{\mathrm{bz}}{l}_{\mathrm{bb}}{k}_{\mathrm{sy}}{\psi }_v=0,\hfill \end{array}$$

(5)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathrm{m}}_{\mathrm{v}}{\stackrel{··}{y}}_v-2c_{\mathrm{sy}}({{\stackrel{·}{x}}_{\mathrm{py}}}_i+{{\stackrel{·}{x}}_{\mathrm{py}}}_{(i+1)})+2k_{\mathrm{sy}}(2y_v-{y_b}_i-{y_b}_{(i+1)})}\\ +4l_{\mathrm{cz}}{k}_{\mathrm{sy}}{\varphi }_v+2l_{\mathrm{bz}}{k}_{\mathrm{sy}}({{\varphi }_b}_i+{{\varphi }_b}_{(i+1)})=0,\hfill \end{array}$$

(6)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathrm{I}}_{\mathrm{vx}}{\stackrel{··}{\varphi }}_v+2l_{\mathrm{cszy}}^2{c}_{\mathrm{sz}}(2{\stackrel{·}{\varphi }}_v-{{\stackrel{·}{\varphi }}_b}_i-{{\stackrel{·}{\varphi }}_b}_{(i+1)})}\\ +2l_{\mathrm{ba}}^2{k}_{\mathrm{sz}}(2{\varphi }_v-{{\varphi }_b}_i-{{\varphi }_b}_{(i+1)})\hfill \\ +2l_{\mathrm{cz}}{k}_{\mathrm{sy}}(2y_v-{y_b}_i-{y_b}_{(i+1)})\hfill \\ +4l_{\mathrm{cz}}^2{k}_{\mathrm{sy}}{\varphi }_v\hfill \\ +2l_{\mathrm{ccsyz}}{c}_{\mathrm{sy}}({{\stackrel{·}{x}}_{\mathrm{py}}}_i+{{\stackrel{·}{x}}_{\mathrm{py}}}_{(i+1)})\hfill \\ +2l_{\mathrm{cz}}{l}_{\mathrm{bz}}{k}_{\mathrm{sy}}({\varphi }_b+{{\varphi }_b}_{(i+1)})\hfill \\ +k_{\mathrm{sroll}}(2{\varphi }_v-{{\varphi }_b}_i-{{\varphi }_b}_{(i+1)})=0,\hfill \end{array}$$

(7)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathrm{I}}_{\mathrm{vz}}{\stackrel{··}{\psi }}_v+2l_{\mathrm{ba}}^2{k}_{\mathrm{sx}}(2{\psi }_v-{\psi }_b-{{\psi }_b}_{(i+1)})-2l_{\mathrm{bb}}{k}_{\mathrm{sy}}(y_b-{y_b}_{(i+1)})}\\ +2l_{\mathrm{bz}}{l}_{\mathrm{bb}}{k}_{\mathrm{sy}}({\varphi }_b-{{\varphi }_b}_{(i+1)})+4l_{\mathrm{bb}}^2{k}_{\mathrm{sy}}{\varphi }_v\hfill \\ +2l_{\mathrm{bb}}{c}_{\mathrm{sy}}({{\stackrel{·}{x}}_{\mathrm{py}}}_i-{{\stackrel{·}{x}}_{\mathrm{py}}}_{(i+1)})-2l_{\mathrm{bcsxy}}{c}_{\mathrm{sx}}({{\stackrel{·}{x}}_{\mathrm{px}}}_i+{{\stackrel{·}{x}}_{\mathrm{px}}}_{(i+1)})\hfill \\ -2l_{\mathrm{bcsyx}}{c}_{\mathrm{sy}}({{\stackrel{·}{x}}_{\mathrm{py}}}_i+{{\stackrel{·}{x}}_{\mathrm{py}}}_{(i+1)})=0,\hfill \end{array}$$

(8)

$$c_{\mathrm{sx}}{{\stackrel{·}{x}}_{\mathrm{px}}}_i+k_{\mathrm{dx}}({x_{\mathrm{px}}}_i-l_{\mathrm{bcsxy}}({{\psi }_b}_i-{\psi }_v))=0,$$

(9)

$$\begin{array}{l}{c}_{sy}{\stackrel{.}{x}}_{pyi}+k_{dy}(x_{pyi}-\left(y_{b_i}-y_v\right)+l_{bcsyz}{\varphi }_{b_i}\\ +l_{ccsyz}{\varphi }_v\pm l_{bb}{\psi }_v-l_{bcsyx}\left({\psi }_{b_i}-{\psi }_v\right))=0,\end{array}$$

(10)

for equations (1) and (2): $i=1(j=1,2)\mathrm{and}2(j=3,4)$ ; for equations (3)–(5): $i=1(j=1)\mathrm{and}2(j=3)$ ; for equations (6)–(8): $i=1$ ; for equations (9) and (10): $i=1,2$ .

Here, ${y_w}_j$ and ${{\psi }_w}_j$ are the lateral and yaw displacements of the jth wheelset, ${y_b}_i$ , ${{\varphi }_b}_i$ and ${{\psi }_b}_i$ are the lateral, roll and yaw displacements of the ith bogie frame, respectively. Also, $y_v$ , ${\varphi }_v$ , ${\psi }_v$ are the lateral, roll and yaw displacements of a carbody, ${x_{\mathrm{px}}}_i$ and ${x_{\mathrm{py}}}_i$ are the piston displacements of the ith yaw and lateral dampers, respectively. $F_j$ and $M_j$ are the linear creep force and creep moment of the jth wheelset and $V$ is the vehicle speed.

Figure 1.

High-speed rail way vehicle model: (a) rear view and (b) top view.

Equations (9) and (10) are the dynamic equations of yaw and lateral dampers. The yaw damper model only acts against the relative yaw motion between the carbody and bogie frames, whereas the lateral damper model affects the lateral, roll and yaw motion between the carbody and bogie frames. The other main parameters are listed in the Appendix.

Design variables, objective functions and design process

The low-frequency damping ratio, closely related to ride comfort, and the critical speed, a measure of hunting stability, were chosen as the objective functions. The goal is to maximise their lowest values, determined by the uncertainty of suspension parameters. The design variables were composed of the longitudinal and lateral stiffness of the rod bush, series stiffness and damping coefficient of the lateral and yaw dampers, respectively. Equations (11) and (12) express the design variables and boundary conditions, respectively.

$$\mathbf{d}=[k_{\mathrm{rx}},k_{\mathrm{ry}},k_{\mathrm{dx}},k_{\mathrm{dy}},c_{\mathrm{sx}},c_{\mathrm{sy}}],$$

(11)

$${\mathbf{L}}_b\le \mathbf{d}\le {\mathbf{U}}_b.$$

(12)

Assuming that the uncertainty range of each suspension parameter is the maximum allowable tolerance of the rubber components and hydraulic dampers for railway vehicle suspensions specified by European norms,^24,25 it can be expressed as:

$${\mathbf{L}}_u=[0.75,0.75,0.8,0.8,0.85,0.85],$$

(13)

$${\mathbf{U}}_u=[1.25,1.25,1.2,1.2,1.15,1.15],$$

(14)

where ${\mathbf{L}}_u$ and ${\mathbf{U}}_u$ are the sets of lower and upper limits for the normalised uncertainties, respectively. The series stiffness of the lateral and yaw dampers is modelled as a series connection between the stiffness of the rubber joints at both ends of the damper and the stiffness of the damper itself, and the uncertainty range of this series stiffness is estimated to be ±20% from numerical calculations for several cases.

The low-frequency damping ratio and critical speed can be calculated from the eigenvalue analysis of the linearised vehicle model. The state-space equation of the vehicle model, expressed as equations (1)–(10), can be written as:

$$\stackrel{·}{\mathbf{x}}=\mathbf{Ax},$$

(15)

$$\mathbf{A}=\left[\begin{array}{cc}0\hfill & \mathbf{I}\hfill \\ \multicolumn{1}{c}{-{\mathbf{M}}^{-1}\mathbf{K}}& -{\mathbf{M}}^{-1}\mathbf{C}\hfill \end{array}\right],$$

(16)

$$\mathbf{x}={[y_{\mathrm{wj}},y_{\mathrm{bi}},y_v,{\psi }_{\mathrm{wj}},{\psi }_{\mathrm{bi}},{\psi }_v,{\varphi }_{\mathrm{bi}},{\varphi }_v,{x_{\mathrm{px}}}_i,{x_{\mathrm{py}}}_i]}^T,$$

(17)

$$(\mathrm{\sigma }\mathbf{I}-\mathbf{A})\mathbf{v}=0,$$

(18)

where $\mathbf{A}$ is the system matrix, and $\mathbf{x}$ is the state vector. $\mathbf{M}$ is a diagonal matrix of the inertia terms, and $\mathbf{C}$ and $\mathbf{K}$ are the assembled matrices for the damping and stiffness coefficients, respectively. ¹⁸ The damped natural frequency ( $f_d$ ) and damping ratio ( ${\zeta }_d$ ) can be determined from the eigenvalue ( $\sigma$ ) of $\mathbf{A}$ :

$$f_d=\mathrm{Im}(\sigma )/2\pi ,$$

(19)

$${\zeta }_d=-\mathrm{Re}(\sigma )/\left|\sigma \right|.$$

(20)

The eigenvalues and eigenvectors of the yaw and lateral dampers were calculated along with those of the other variables; however, they were ignored in the eigenmodes of the vehicle. Thus, the linearised vehicle model has 17 eigenmodes; the modes with $f_d$ of 3 Hz or less were defined as low-frequency modes in this study. In addition, because the low-frequency mode with the minimum damping ratio dominates the low-frequency eigenbehaviour of the vehicle, its damping ratio represents the low-frequency damping ratio of the vehicle. Thus, the low-frequency damping ratio can be expressed as equation (21). Here, a separating frequency of 3 Hz was determined considering the kinematic hunting frequency of the bogie calculated by an analytical approximation.^22,26

$${\zeta }_{\mathrm{ld}}={{\zeta }_d|}_{f_d\le 3\mathrm{Hz}}$$

(21)

As the $\mathbf{C}$ and $\mathbf{K}$ matrices in equation (16) include the linear creep force and creep moment dependent on the vehicle speed, $f_d$ and ${\zeta }_d$ also change according to the vehicle speed, a vehicle speed should be given to calculate ${\zeta }_{\mathrm{ld}}$ . In this study, a target speed range was first determined instead of specifying a single or several vehicle speeds, and the lowest ${\zeta }_{\mathrm{ld}}$ in this speed range was adopted. This method can improve ${\zeta }_{\mathrm{ld}}$ over the entire vehicle speed range because it finds the new lowest ${\zeta }_{\mathrm{ld}}$ to improve at every iteration of the design process; thus, when the current lowest ${\zeta }_{\mathrm{ld}}$ has been improved and is no longer the lowest value, a new one is found and improved from the next iteration. By setting ${\zeta }_m$ as the minimised ${\zeta }_{\mathrm{ld}}$ owing to the parameter uncertainties within ${\mathbf{L}}_u$ and ${\mathbf{U}}_u$ , the robustness of ${\zeta }_{\mathrm{ld}}$ of the railway vehicle against the given uncertainty range can be secured by designing suspension parameters such that ${\zeta }_m$ becomes equal to or greater than the required value. Thus, the objective function for ${\zeta }_{\mathrm{ld}}$ is expressed as equation (22)

$${\zeta }_m=\underset{\begin{array}{c}{v}_1\le v\le v_2\hfill \\ \multicolumn{1}{c}{{\mathbf{L}}_u\le {\mathbf{u}}_1\le {\mathbf{U}}_u}\end{array}}{min}{\zeta }_{\mathrm{ld}}(v,{\mathbf{u}}_1\odot \mathbf{d}),$$

(22)

where ${\mathbf{u}}_1$ is the set of normalised parameter uncertainties to minimise ${\zeta }_{\mathrm{ld}}$ , and operator ⊙ denotes element-wise multiplying. Therefore, ${\mathbf{u}}_1\odot \mathbf{d}$ is the set of suspension parameters changed by their uncertainties.

In this study, the critical speed $v_c$ for bogie hunting was defined as the vehicle speed at which even one of the high-frequency damping ratios ( ${\zeta }_{\mathrm{hd}}$ ) becomes zero. Equation (23) expresses the objective function for the critical speed:

$$v_m=\underset{{\mathbf{L}}_u\le {\mathbf{u}}_2\le {\mathbf{U}}_u}{min}{v}_c({\mathbf{u}}_2\odot \mathbf{d}),$$

(23)

where ${\mathbf{u}}_2$ is the set of normalised parameter uncertainties to minimise $v_c$ . Because the wheel-rail contact characteristics also significantly influence the target performance, different wheel-rail contact conditions were applied to each objective function so that the suspension parameters were designed under the worst conditions. Table 1 lists the wheel-rail contact parameters applied to each objective function.

Table 1.

Wheel/rail contact parameters of the linearised vehicle model.

Property	W1 condition	W2 condition
	Low-frequency damping ( ${\zeta }_{\mathrm{ld}}$ )	Critical speed ( $v_c$ )
Equivalent conicity ( $\lambda$ ) (at 3 mm wheel shift)	0.0408	0.3382
$f_{11}$	1.0639 × 107	1.4420 × 107
$f_{12}$	2.5741 × 104	3.2573 × 104
$f_{22}$	103.0437	692.2468
$f_{33}$	1.1940 × 107	1.8106e × 107

Figure 2 presents the schematic of the proposed design method. As shown in Figure 2(a), the workflow comprises the main loop that optimises the suspension parameter $\mathbf{d}$ for acquiring the Pareto optimality of ${\zeta }_m$ and $v_m$ and the sub-loop that receives $\mathbf{d}$ from the main loop, calculates ${\zeta }_m$ and $v_m$ for ${\mathbf{u}}_1$ , ${\mathbf{u}}_2$ and output the results to the main loop. The multi-objective genetic algorithm (MOGA) was applied for the main loop. In the sub-loop, as shown in Figure 2(b), equations (22) and (23) are redefined as problems in finding the values of seven variables (six normalised parameter uncertainties and one vehicle speed) that minimise ${\zeta }_{\mathrm{ld}}$ and the absolute value of ${\zeta }_{\mathrm{hd}}$ . The interior-point algorithm ²⁷ (IPA) with high computational efficiency was applied; however, multiple IPA method with several initial values was used to increase reliability. Table 2 present a summary of the design problems, respectively. The system equation calculations, IPA and MOGA were implemented using MATLAB by MathWorks^®. MOGA from MATLAB was constructed based on NSGA-II.^28,29 As shown in Figure 2(b), a next-generation (modified) population was created in the MOGA through selection, crossover and mutation, based on biological evolution. In this study, tournament, intermediate and feasible adaptation methods ²⁹ were applied for selection, crossover and mutation, respectively. In each iteration of the two algorithms, the feasibility or fitness was evaluated and compared with the stop conditions. The stop condition for each algorithm is provided in Table 2.

Figure 2.

Schematic of the proposed design method: (a) workflow and (b) optimisation algorithm.

Table 2.

Design problems.

Level (Algorithm)	Problems	Conditions
Main (MOGA)	$\underset{{\mathbf{L}}_b\le \mathbf{d}\le {\mathbf{U}}_b}{max}\{{\zeta }_m\left(\mathbf{d}\right),v_m(\mathbf{d})\}$	Population size: 300 Function tolerance: $1\times {10}^{-4}$ Stall generation: 100 ${\mathbf{L}}_b\le \mathbf{d}\le {\mathbf{U}}_b$ ${\mathbf{L}}_b=[10,4,3,1.9,1,0.02]\times {10}^6$ ${\mathbf{U}}_b=[35,8,28,6,2.8,0.07]\times {10}^6$
Sub (IPA)	${\zeta }_m=\underset{\mathbf{u},v}{min}{\zeta }_{\mathrm{ld}}(v_{\mathrm{ld}},\mathbf{u}\odot \mathbf{d})$	Function tolerance: $1\times {10}^{-3}$ Step tolerance: $5\times {10}^{-3}$ Optimality tolerance: $5\times {10}^{-5}$ ${\mathbf{L}}_u\le \mathbf{u}\le {\mathbf{U}}_u$ ${\mathbf{L}}_u=[0.75,0.75,0.8,0.8,0.85,0.85]$ ${\mathbf{U}}_u=[1.25,1.25,1.2,1.2,1.15,1.15]$ $36\le v_{\mathrm{ld}}\le 360\mathrm{km}/\mathrm{h}$
	$\mathrm{g}(v_m,\cdots )=\underset{\mathbf{u},v}{min}\left|{\zeta }_{\mathrm{hd}}(v_{\mathrm{hd}},\mathbf{u}\odot \mathbf{d})\right|$

Design results

Figure 3 shows the Pareto solutions obtained using the proposed method. Because the objective functions are defined by equations (22) and (23), the Pareto fronts, indicated by blue circles in Figure 3(a), show minimal performances owing to the parameter uncertainties; however, they were maximised as much as possible without degrading each other by MOGA of the outer-loop. Red squares represent the nominal performances when the optimised suspension parameters without parameter uncertainties are applied, which were acquired from additional calculations when the optimisation process was completed. Because the minimum required performance was determined to be over 10% for the lowest ${\zeta }_{\mathrm{ld}}$ under the W1 condition and over 300 km/h for $v_c$ under the W2, a Pareto front with the lowest ${\zeta }_{\mathrm{ld}}$ of 10.3% and $v_c$ of 304 km/h was chosen as the final design result of the railway vehicle. Its nominal values are 13.9% and 427 km/h, respectively. Figure 3(b) to (g) present the Pareto sets from MOGA. The blue circles and red squares indicate the nominal values of the lowest ${\zeta }_{\mathrm{ld}}$ and $v_c$ , respectively. The dashed line denotes the optimised suspension parameters for the selected Pareto front. $k_{\mathrm{rx}}$ and $k_{\mathrm{ry}}$ were located close to their lower limits, whereas $k_{\mathrm{dy}}$ and $c_{\mathrm{sx}}$ were close to their upper limits. The nominal performance curves based on $k_{\mathrm{dx}}$ and $c_{\mathrm{sy}}$ exhibited similar shapes.

Figure 3.

Design results: (a) Pareto fronts (blue circles), nominal performances (red squares) and minimum required performance (dashed line), and Pareto set for (b) $k_{\mathrm{rx}}$ , (c) $k_{\mathrm{ry}}$ , (d) $k_{\mathrm{dx}}$ , (e) $k_{\mathrm{dy}}$ , (f) $c_{\mathrm{sx}}$ and (g) $c_{\mathrm{sy}}$ ; from (b) to (g), dashed lines, blue circles and red squares indicate the optimal suspension parameters corresponding to the selected Pareto front, the nominal values of the lowest damping ratio and the critical speed, respectively.

Figure 4(a) displays ${\zeta }_{\mathrm{ld}}$ and its modal frequency, $f_{\mathrm{ld}}$ , according to the vehicle speeds. In the reference vehicle model (M_O), ${\zeta }_{\mathrm{ld}}$ decreased rapidly to approximately 4.7% at the mid-speed section and increased gradually as the vehicle speed increased. $f_{\mathrm{ld}}$ increased with vehicle speed, indicating that the dominant mode is the bogie lateral movement coupled with the carbody sway motion. The designed vehicle model (M_N) with optimal suspension parameters exhibited a solid damping characteristic over the entire speed range and had a constant ${\zeta }_{\mathrm{ld}}$ at speeds above approximately 120 km/h. $f_{\mathrm{ld}}$ was also approximately constant, regardless of the vehicle speed, and the dominant mode was the carbody’s lateral or lower sway motion. For a vehicle model (M_D) with suspension parameters changed from optimal values by the uncertainties such that the lowest ${\zeta }_{\mathrm{ld}}$ is minimised, ${\zeta }_{\mathrm{ld}}$ decreased gradually with vehicle speed and $f_{\mathrm{ld}}$ increased. Thus, the dominant mode is the coupled motion of the two bogies and the carbody, which is similar to M_O. However, this oscillation mode can be successfully suppressed because the lowest ${\zeta }_{\mathrm{ld}}$ is approximately 10%. For a vehicle model (M_C) with suspension parameters degraded from optimal values such that $v_c$ is minimised, ${\zeta }_{\mathrm{ld}}$ is higher than that of M_N for the entire vehicle speed, and the dominant mode is the same as that of M_N. Figure 4(b) shows ${\zeta }_{\mathrm{hd}}$ and its modal frequency, $f_{\mathrm{hd}}$ , for W1. The M_O had a higher ${\zeta }_{\mathrm{hd}}$ in the mid-speed section than the other models; however, all models had similar values in the high-speed section. The $v_c$ of all the models under the W1 condition were higher than 550 km/h because their ${\zeta }_{\mathrm{hd}}$ values were larger than zero. The dominant mode of M_O changed from bogie yaw mode to bogie lateral mode at a vehicle speed of approximately 60 km/h and changed again to wheelset yaw mode at a vehicle speed of approximately 160 km/h. For M_N, M_D and M_C, although $f_{\mathrm{hd}}$ slightly increased with vehicle speed due to a slight change in mode shape, the dominant mode was bogie yaw mode regardless of the vehicle speed. Figure 4(c) shows ${\zeta }_{\mathrm{ld}}$ and $f_{\mathrm{ld}}$ for W2. The ${\zeta }_{\mathrm{ld}}$ values of M_O and M_D increased remarkably and were almost constant over the entire vehicle speed range. Their dominant mode also changed from coupled motion to a lateral or lower sway motion of the carbody. The ${\zeta }_{\mathrm{ld}}$ of M_N and M_C increased slightly because the dominant mode was already a lateral or lower sway motion of the carbody. Figure 4(d) shows ${\zeta }_{\mathrm{hd}}$ and $f_{\mathrm{hd}}$ under the W2 condition. The ${\zeta }_{\mathrm{hd}}$ of the M_C was zero at a vehicle speed of approximately 300 km/h. Thus, the $v_c$ of this model under W2 was approximately 300 km/h. The other models had $v_c$ above 400 km/h under the W2 condition. The M_O and M_D had similar $v_c$ of approximately 410 km/h, and the M_N had $v_c$ of approximately 430 km/h. The unstable mode of each model was very similar, with wheelset lateral and yaw modes accompanied by bogie lateral modes.

Figure 4.

Damping ratio and modal frequency of the high-speed vehicle model: (a) low-frequency and (b) high-frequency damping ratio under W1 condition, (c) low-frequency and (d) high-frequency damping ratio under W2 condition.

Since each of the 105 Pareto fronts represents the minimised performance by the parameter uncertainties although optimised by MOGA, it is possible to consider the corresponding ${\mathbf{u}}_1$ and ${\mathbf{u}}_2$ as the collateral results. Their distribution is exhibited in Figure 5. As shown in Figure 5(a), (c) and (e), the two performances have a strong contradictory relationship with uncertainties of the longitudinal suspension parameters. For the lateral suspension parameters in Figure 5(b), (d) and (f), when the parameters have their lowest values within the uncertainty limits, both performances are minimised with a high probability. However, for $k_{\mathrm{ry}}$ , both ends of the uncertainty have a similar frequency (probability) to minimise the critical speed. Thus, an interesting implication for the maintenance process is that the focus of attention should be the longitudinal suspension elements whose parameters have excessively deviated from their nominal values in both directions and the lateral elements whose parameters have significantly changed in the decreasing direction.

Figure 5.

Histogram of suspension uncertainty minimising the vehicle performances: (a) $k_{\mathrm{rx}}$ , (b) $k_{\mathrm{ry}}$ , (c) $k_{\mathrm{dx}}$ , (d) $k_{\mathrm{dy}}$ , (e) $c_{\mathrm{sx}}$ and (f) $c_{\mathrm{sy}}$ .

Validation

Variation of vehicle performances for arbitrary parameter uncertainties

Investigating the performance distribution for arbitrary parameter uncertainties to evaluate the validity of the design results is necessary, as the proposed method only considers the minimised performance of the vehicle. This study used 100,000 sets of normalised random uncertainties, modelled to be uniformly distributed between ${\mathbf{L}}_u$ and ${\mathbf{U}}_u$ . According to Figures 3(a) and 4(a), the lowest ${\zeta }_{\mathrm{ld}}$ of the M_N can be reduced to 10.3% owing to the parameter uncertainty; however, as shown in Figure 6(a), the probability of the lowest ${\zeta }_{\mathrm{ld}}$ being higher than 12% is 99.86%. The probability of the lowest ${\zeta }_{\mathrm{ld}}$ being between 12.25% and 12.5% had the highest probability, and the probability of the interval including the nominal lowest ${\zeta }_{\mathrm{ld}}$ was approximately 6.5%. Figure 6(b) presents the distribution of $v_c$ with a 5% damping criterion for the W1 condition, as 5% damping is usually used for evaluating $v_c$ as a conservative criterion in railway vehicle engineering. ³⁰ The maximum probability appears in the interval between 620 and 640 km/h, where the vehicle model has $v_c$ of more than 500 km/h with a probability of 94.84%. Figure 6(c) shows the distribution of the lowest ${\zeta }_{\mathrm{ld}}$ for W2. As stated earlier, in the W2 condition, the lowest ${\zeta }_{\mathrm{ld}}$ increased slightly compared with the W1 condition. All results for 100,000 uncertainties show the lowest ${\zeta }_{\mathrm{ld}}$ to be higher than 12%, and the interval between 12.5% and 12.75% had the highest probability. The $v_c$ for the parameter uncertainties under the W2 condition are distributed as shown in Figure 6(d). The damping criterion was set to zero percent to compare with the design results, the same as used in the design process. The probability of $v_c$ reaching the minimum value evaluated from the selected Pareto front was estimated to be extremely low. $v_c$ was predicted to be higher than 350 km/h with a probability of 99.29% and 400 km/h with a probability of 74.99%.

Figure 6.

Distribution of vehicle performances for uniformly random uncertainties in suspension parameters: (a) the lowest ${\zeta }_{\mathrm{ld}}$ , (b) $v_c$ for 5% damping under W1 condition, (c) the lowest ${\zeta }_{\mathrm{ld}}$ and (d) $v_c$ under W2 condition.

Consequently, in most arbitrary uncertainties, the two target performances of the vehicle model are higher than the minimum required, and the probability that they are lower than the minimum required is zero and equal to the minimum required is extremely low. It implies that the proposed method worked well as expected, and the minimum performance of the vehicle model is always secured regardless of the uncertainties of suspension parameters.

Simulation

As the suspension parameters were designed using a linearised vehicle model, several simulations were performed with VI-Rail – a commercial software for simulating the railway vehicle dynamics of VI-grade GmbH^©, to estimate the validity of the achieved performances. The vehicle model in VI-Rail had nonlinear characteristics such as nonlinear damping coefficient, nonlinear secondary air spring and nonlinear wheel/rail contact. The critical speed was evaluated for the M_O, M_N and M_C. In combination with the 60E1 rail, a wheel profile that satisfied the W2 condition was used. Figure 7(a) and (b) display the high-speed railway vehicle modelled in VI-Rail and the equivalent conicities of the wheel-rail profiles applied to satisfy conditions W1 and W2. Figure 7(c) to (e) show the lateral displacements of the leading wheelset of each vehicle model for the ramped track irregularity in the lateral direction. There is good agreement with the predicted values from the linearised vehicle model (Figure 4(d)), which implies that the combination of wheel-rail profiles used in VI-Rail has the characteristic of supercritical Hopf bifurcation. ³¹

Figure 7.

Simulation results for critical speed: (a) vehicle model in the VI-Rail, (b) equivalent conicities of applied wheel profiles with 60E1 rail, (c) M_O, (d) M_N and (e) M_C.

The ride comfort was estimated for the M_O, M_N and M_D. When the vehicle models travelled on a tangent track of 10 km with ERRI low rail irregularities, the accelerations were measured at the front, centre and rear positions of the carbody floor. The front and rear positions were on the centre of two bogies. In this simulation, the wheel profile for the W1 condition was applied to the vehicle models. Figure 8(a) shows the measured lateral accelerations at the front and rear positions in the M_N and M_D with a vehicle speed of 350 km/h, and Figure 8(b) to (d) show the estimated ride comforts for each model in accordance with EN 12299. ³² N_MV indicates the mean comfort of the test track section; the smaller the value, the better the ride quality. As shown in Figure 8(b), in the M_O, the differences between the indices of the front and rear positions are quite large because of the coupled mode of the bogies and the carbody. The ride comfort improved as the vehicle speed increased, which is consistent with the curve shape of ${\zeta }_{\mathrm{ld}}$ for the corresponding speed range in Figure 4(a). For instance, ${\zeta }_{\mathrm{ld}}$ takes values of 5.00%, 4.77%, 5.13% and 5.75% at vehicle speeds of 200, 250, 300 and 350 km/h, and N_MV of the centre position is evaluated as 2.79, 2.86, 2.72 and 2.51, respectively. The N_MV of rear position is 4.27–4.61, which can be assessed to be ‘uncomfortable’ or ‘very uncomfortable’ according to EN 12299. Other positions are 2.51–3.35, which are ‘medium’ scale. The ride comfort for M_N and M_D improved significantly compared with that of the M_O as shown in Figure 8(c) and (d). The N_MV of both models has the maximum on the rear position at the vehicle speed of 350 km/h. Their values are 1.18 and 1.60, indicating ‘very comfortable’ and ‘comfortable’ scales, respectively, despite the worst case. Therefore, considering uncertainties, the lowest damping ratio of 10% in the low-frequency mode is reasonable as a minimum target performance when designing the suspension parameters of the vehicle model used in this study. Among the 200, 250, 300 and 350 km/h speeds, the difference in ride comfort between M_N and M_D models was smallest at 200 km/h because the difference in ${\zeta }_{\mathrm{ld}}$ was the smallest. Additionally, as the vehicle speed increased, the difference increased because the ${\zeta }_{\mathrm{ld}}$ value of the M_D decreased slightly, whereas that of the M_N remained constant. Ride comfort at the centre position is better than at other positions in all models because the lateral accelerations caused by the carbody’s yaw or two bogies’ out-of-phase lateral motion cancel each other out at the centre position.

Figure 8.

Simulation results for ride comfort: (a) lateral accelerations of the carbody floor for M_N (top) and M_D (bottom), and ride comforts for (b) M_O, (c) M_N and (d) M_D, respectively.

Conclusions

Uncertainties in the parameters of the suspension elements widely used in railway vehicles, such as rubber bushes and hydraulic dampers, are often caused by the manufacturing process, service period and the environment. Therefore, designing suspension parameters such that the railway vehicle satisfies the minimum required performance despite their variations owing to uncertainty is necessary. Satisfying this requirement is guaranteed when the performance minimised by the uncertainty is equal to or higher than the minimum required performance. Hence, a design method for suspension parameters securing the minimum required performance of a railway vehicle despite the uncertainty is proposed. This method comprises two loops: the sub-loop finds the performances of the railway vehicle minimised by the parameter uncertainties, and the main loop seeks solutions where those minimum performances reach Pareto optimality. The proposed method also requires no information for parameter uncertainties apart from their range.

The low-frequency damping ratio in the low equivalent conicity condition and critical speed in the high equivalent conicity condition, which are dynamic characteristics representing the ride comfort and running stability of a railway vehicle, were selected as design targets. Since the damping ratio changes according to the vehicle speed and suspension parameters, the lowest low-frequency damping ratio in the service speed range was used instead of the damping ratios at certain specified vehicle speeds. Consequently, the overall low-frequency damping ratio in the service speed range could be improved.

The distribution of the parameter uncertainties acquired from the Pareto set suggests that the longitudinal suspension elements whose parameters significantly deviate from their nominal values in both directions and the lateral suspension elements whose parameters excessively change in the decreasing direction be devoted particular attention within the scope of maintenance.

A performance analysis using 100,000 random samples further verified the viability of the proposed method. The designed vehicle performed better than the minimum required performance in most parameter uncertainties, and the probability of being less than the minimum required was zero, indicating that the proposed method worked well as expected.

Because of the otherwise excessive computational cost of the proposed design method, a linearised vehicle model and its dynamic behaviour in the frequency domain were used for optimisation of suspension design, although they could only represent some characteristics of a railway vehicle. However, the simulation results for the target performance of this study using the railway vehicle dynamics software VI-Rail demonstrate the validity of the proposed method and the linearised vehicle model. The result for $v_c$ agreed well with the linearised model, and the ride comfort index (N_MV), according to EN 12299, significantly improved in the designed vehicle model, as expected from the low-frequency damping ratio of the linearised model.

Notably, the quantitative design results are limited to the vehicle model used in this study. However, the proposed design procedure can be applied to designing suspension parameters for almost all railway vehicles with appropriate objective functions calculated from various methods. It can also be extended to robust design of suspension parameters for other uncertainties, such as changes in wheel-rail contact characteristics due to wear or changes in vehicle mass and moment of inertia due to passenger/cargo weight, by simply adding uncertain parameters and their normalised ranges to the sub-loop.

In future studies, the actual range of suspension parameter uncertainties will be investigated from a number of high-speed railway vehicles in operation. Also, shapes of their wheels will be measured to determine the actual uncertainty range of wheel-rail contact parameters. Then, the suspension parameters will be re-designed with their actual uncertainty ranges and actual uncertainty range of wheel-rail contact parameters. Finally, the vehicle performance will be compared with this study.