Effects of the Seal Wire on the Nonlinear Dynamics of the Aircraft Engine Turbine Blades

Aircraft engines are complicated engineering systems containing thousands of components made of various materials in complex geometries [1–3]. As one of these components, metallic wires are currently widely used in aircraft engines to ensure aerodynamic sealing between turbine blades, which also contribute to the reduction of vibration amplitudes through frictional damping introduced by the contact between the wires and the blades [4,5]. Thus, the contacts between those components are not only functional for the aerodynamics of the system but also crucial from a purely structural point of view as they can influence the overall system behavior, just like the contacts between other components of these systems. This brings the necessity of predicting/modeling the nonlinear behavior due to these contacts accurately to make the design process easier and more efficient in terms of time and cost [6].

As one of the early works of modeling turbine blade vibration, taking the isotropic and anisotropic friction forces into account, Goodman and Klumpp dealt with frictional damping analytically [7]. Zmitrowicz developed a mathematical formulation to capture the physics of the bladed disk systems observed in experiments [8]. Armand et al. formulated a semi-analytical approach to model the contact capable of representing the bladed disks [9]. Charleux et al. performed numerical and experimental investigations on friction damping in bladed disk assemblies [10]. Pesek et al. utilized experiments and finite element method to carry out three-dimensional modeling of bladed disks [11]. Szwedowicz et al. conducted numerical modeling and experiments to explore the energy dissipation of a thin-walled frictional damper designed for turbine blades [12]. Afzal et al. inspected the plausible damping capabilities of a strip damper on a turbine blade [13]. As a common problem, the prediction of vibrational energy dissipation by under-platform dampers has been in the spotlight of several research items available in the literature [12,14–16].

Usually, the computational costs of these models are quite significant. Therefore, model order reduction methods are developed for generating more feasible models in terms of computational load and accuracy [17]. Thus, researchers have generated techniques to perform reduced-order modeling. Taking advantage of cyclic symmetry is a commonly used modeling method for aircraft engine applications. Grolet and Thouverez presented a modeling technique for geometrically nonlinear systems with cyclic symmetry [18,19]. Petrov developed an efficient method to model cyclic symmetric structures containing nonlinearities under multiharmonic excitations [20]. Again Petrov conducted reduced-order modeling on a bladed disk application case with frictional contact nonlinearities [21]. Dastani et al. investigated the behavior of a dovetail joint representing the root of a blade using the Craig–Bampton approach and compared the results to experiments [22]. Renson et al. also employed Craig–Bampton to build and solve the reduced-order models of complex aerospace structures [23]. Krack et al. presented a detailed literature review about vibration modeling of bladed disks taking friction into account [24].

There are several ways to model contacts with friction. Coulomb and Jenkins models are among the most common approaches to model friction. Each has its advantages and disadvantages. To model the contacts more cost-efficiently and yet accurately, Şanl itürk and Ewins developed a contact formulation for two-dimensional point contacts that can employ various contact models, as well as experimental data [25]. Petrov and Ewins presented an analytical formulation for frictional contact interface elements [26]. Firrone and Zucca presented some of the most common numerical contact modeling approaches [27]. These studies employed the harmonic balance approach to derive a solution [28].

Even though there are studies present in the literature focusing on the problem of modeling the dynamics of aircraft engine turbine structures, including contact nonlinearities, current nonlinear prediction methods struggle to accurately replicate the observed effect of the seal wire on the blade amplitude and frequency response. Inspection of contact and wear patterns on engine hardware indicate that contact assumptions in the prediction methods may not be fully representative [5,29,30]. Besides, there is a very limited number of published studies available in the literature concerning sealing applications directly. One of them is Nakane et al., which presented a methodology to develop seals for high-performance applications in gas turbines investigating various parameters, including leakage and rigidity of the seal [31]. Based on the results acquired in that study, leaf seals were found to be efficient and could be recommended for potential future use. Moreover, Demiroğlu et al. took advantage of a beam formulation to model the dynamics of the sealing to explore the effects of several factors related to the operational conditions and the shape of the seal on the tip force of seals [32]. The linear and nonlinear cantilever beam formulations were solved and compared. As one of the unique studies regarding sealing, Petrov dealt with the modeling of nonlinear vibrations of bladed disks in contact with seal wires and under-platform dampers employing explicit large-scale finite element models with various types of frictional dampers calculating the damping for several contact conditions [33].

The objective of this work is to investigate the modeling and dynamic analysis of this complex contact problem to determine the source of existing prediction errors of the current conventional industrial models. The shrouded blades of the intermediate-pressure turbine of a high bypass ratio turbofan engine are modeled using the finite element method taking advantage of the cyclic symmetry in the static and frequency domain models. The seal wires and the so-called lock plates used to maintain the blades in place are also modeled using the same method. The finite element model includes the contacts at the root of the blade, the lock plate-blade-disk contact as well as the contact between the seal wire and blade hub. The geometries used in this study are representative of the different components of such an engine structure. The simulation of the response consists of time domain and frequency domain analysis. The time domain calculations are mainly based on finite element analyses conducted in abaqus, while the frequency domain studies are performed using an in-house harmonic balance algorithm. Time domain simulations include transient implicit simulations to help determine the precise position of the seal wire between the blades and improve the modeling of the contact region. Different contact parameters and initial conditions at the interface between components are considered, and the energy dissipation associated with each individual contact is investigated.

This study is organized into four sections following this introductory section. In Sec. 2, the numerical models built are described. The numerical models are categorized as frequency and time domain simulations, and the different and common assumptions are specified in Secs. 3 and 4. The results are presented and discussed in these sections as well. Finally, the conclusions drawn from this research are expressed in Sec. 5.

All the numerical simulations performed within the scope of this study are based on the large-scale finite element model of the bladed disk system of an aircraft engine's turbine. The geometry consists of one sector of the bladed disk system (Fig. 1). Since it is a one-sector model, the lock plate only interacts with the components on that sector and not with any other externally predefined item like springs and dampers. The model takes into account the effects of the steady-state temperature distribution on the material properties, which is presented in Fig. 2. Besides, the rotational body forces at 7500 rpm and representative static aero-loads (pressure) acting on the system are also included in the models as shown in Fig. 3.

The finite element model, which is the basis of the entire modeling procedure, consists of 117,964 quadratic tetrahedral elements of type C3D10, 200,112 nodes, and 600,336 degrees-of-freedom. There are three components in the assembly that are in contact with the bladed disk, as displayed in Fig. 4. The number of nodes that are involved in the contact definition (or tie constraint definition) is presented in Table 1.

Here, K(u) denotes the stiffness matrix, u stands for the displacement vector, and f(u) indicates the force/load vector. The stiffness matrix is a function of the deformation state at that instant leading to nonlinearity. The displacement vector u can be computed via Eq. (1) using a numerical method like the modified Newton–Raphson method.

The nonlinear general static analysis in abaqus accounts for the cyclic symmetry of the system by solving the equations of equilibrium for a single sector and then applying the cyclic symmetry boundary conditions to obtain the solution for the complete system. In these simulations, cyclic symmetry is imposed on the disk to represent the boundary conditions on the lateral surfaces of the disk shown in Fig. 1 [20].

Apart from the mechanical loading, abaqus can also consider thermal effects in the analysis, which are modeled using the coefficient of thermal expansion and temperature-dependent material properties. Even though the static analysis takes these nonlinearities into account, the contact nonlinearities are excluded, and the contacting surfaces are tied rigidly.

In the analysis of the cyclically symmetric bladed disk structure with the seal wire, the geometry and loading of the system display perfect symmetry around the axis of rotation of the turbine. This allows one to reduce the size of the finite element model by only considering a single sector of the system, which is then expanded for the other sectors to get the results for the whole system. The boundary conditions applied to the sector must also exhibit cyclic symmetry, which means that they are periodic with respect to the angle of rotation, which is achieved by imposing the cyclic symmetry boundary conditions.

To implement cyclic symmetry boundary conditions in abaqus, a cyclic symmetry interaction must be defined, which involves specifying the sector, rotational symmetry axis, the cyclic symmetry boundary surfaces, and the type of analysis. This option is only available in static and frequency domain simulations. abaqus then solves the equations of equilibrium for the sector and applies the cyclic symmetry boundary conditions to obtain the solution for the complete system.

where $ωnat$ is the natural frequency, $ϕ$ is the mode shape vector, and $M$ is the mass matrix.

The results of the static analysis are then used as input for the modal analysis step. The results can provide insight into the dynamic characteristics of the structure under the combined effects of mechanical loading and temperature. Furthermore, in the modal analysis, no material damping is taken into account while constructing these models. None of these simulations considers the nonlinearities due to the contacts of the components at this stage. The contacting component interfaces are tied rigidly. This provides the baseline for the next step, which is the computation of the forced response of the system, including the nonlinearities.

Here, C stands for the damping matrix, u denotes the vector of displacement, $u˙$ denotes the vector of velocity, $u˙˙$ denotes the vector of acceleration, t indicates time, and $f(t,u)$ indicates the nonlinear forcing/excitation function.

The harmonic balance method is used to approximate the periodic response/solution of Eq. (3) as a sum of harmonic functions

where N is the number of harmonics included in the approximation, ω is the frequency of the periodic excitation, and a_i and b_i are the amplitude coefficients.

Substituting Eq. (4) into Eq. (3) and equating the coefficients of each harmonic term, a set of nonlinear algebraic equations for the unknown amplitude coefficients a_n and b_n are obtained.

Then numerical continuation method is used to solve these algebraic equations for a range of excitation frequencies, starting from a known solution at a reference frequency. By gradually increasing or decreasing (predictor–corrector) the excitation frequency, the continuation method traces out a path of solutions in the frequency domain, allowing one to explore the system's response over a range of frequencies.

Using this approach, the forced response of the bladed disk structure contacting the seal wire is calculated. The nonlinear contact elements are inserted between the rigidly tied nodes to achieve an accurate representation of the contact nonlinearities that are present in the current system [26]. The normal and tangential stiffness values for the contacts are taken as $30,000 N/mm$⁠, and the friction coefficient as 0.3 [34–36]. Only the contact between the seal wire and the blade hub is considered in the frequency domain model within the scope of this study. To finalize the frequency domain model, the system is excited with a point force acting on the pressure side of the blade in the surface's normal direction. Figure 5 displays the excitation and the measurement location for the frequency domain forced response analysis. Then, the forced response of this frequency domain model is calculated using the forced response suite, an in-house software developed at the Rolls-Royce Vibration University Technology Center, Imperial College London [6,21,26].

The results of the frequency domain calculations are presented in nondimensional form. The frequencies are divided by the first natural frequency of the system, whereas the displacements are divided by the static displacements calculated in the static step of the original finite element model.

Figure 6 shows the first mode shape of the bladed disk, which corresponds to the first flap mode of the blade for the nodal diameter 0, and Fig. 7 displays the first five modes for all the nodal diameters calculated with the described finite element model taking only the effects of the temperature and the prestress into account and no contact nonlinearity.

The results of the forced response analysis, including only the contact nonlinearity between the seal wire and the bladed disk, are presented in Figs. 8 and 9. The results in Fig. 8 show the results for a constant excitation force amplitude (F_exc) calculated with the first and the first ten harmonics, whereas Fig. 9 is computed using only the first harmonic with different levels of excitation forces. Figure 8 shows that both sets of calculation yield almost overlapping results while the calculations taking the first ten harmonics into account are significantly more expensive. The only noticeable difference between these sets of results can be noted as the overshoot of the displacement amplitude of the simpler model that only includes the first harmonic in the calculations. However, considering the differences between the computational costs of the models and the similarity between these sets of results, further investigations are carried out with the simpler model with only the first harmonic. Figure 9 reveals that the system displays a softening type of nonlinearity with increasing amplitudes of the excitation forces. Under smaller amplitudes of excitation forces, the displacements also remain relatively small, which causes the seal wire to be stuck with the bladed disk assembly, and no significant slip is observed. Increasing the excitation force amplitude leads to larger displacements/deformations in the system. Eventually, the friction forces between the seal wire and the bladed disk are overcome by the inertia of the seal wire and cause a slip. Thus, the stiffness of the system slightly drops, and the softening nonlinearity mentioned earlier is observed. Another aspect to note about these results would be that the effects of the nonlinearity observed here are quite weak. Considering the mass and stiffness properties of the seal wire, it needs to be stated that the contribution of the seal wire to the entire system's stiffness and mass is quite limited. Therefore, the nonlinearity caused by the seal wire does not change the behavior of the system dramatically.

where $fc$ is the contact force vector due to contact with other bodies. The contact force vector has two components: $fnorm$⁠, the normal force vector, and $ftang$⁠, the tangential force vector, respectively.

Here, F_f is the frictional force (concentrated), μ is the coefficient of friction, F_norm is the normal force (concentrated), and $u˙rel$ is the relative velocity between the two surfaces. The $sign(u˙rel)$ function is used to determine the direction of the frictional force, which is opposite to the direction of motion.

The frictional force is limited to a maximum value of $μFnorm$⁠. If $|Ff|$ is less than or equal to $μFnorm$⁠, the frictional force is given by the first case in Eq. (6). If $|Ff|$ is greater than $μFnorm$⁠, the frictional force is limited to $μFnorm$ and the direction is determined by the second case in Eq. (6).

The solution of this equation of motion is then obtained using the abaqus solver, which is an iterative process utilizing the Newmark time integration with the Newton–Raphson method.

where u is the nodal displacement vector, r and z are the radial and axial positions/coordinates, respectively, and n is the number of cyclically symmetric sectors.

To implement periodic boundary conditions in a finite element model, the boundary conditions are imposed on the nodes of the model on the opposite surfaces in the circumferential (θ) direction. Thus, the nodal displacement values of the finite element model on one surface are related to the corresponding nodal values on the opposite surface in the θ direction.

The bladed disk model displays cyclic symmetry, which makes it convenient to formulate the boundary conditions (Eq. (7)) in cylindrical coordinates [39,40]. In this coordinate system, if the periodic boundary conditions only apply on the opposite surfaces in the θ circumferential direction, the system is assumed to be invariant under rotations of $2π$ around the axis. This is equivalent to a unit cell with periodic boundary conditions in the circumferential direction.

These models are setup to investigate the effect of the initial position of the seal wire on the equilibrium and the dynamics of the entire system. Considering only one sector of the bladed disk assembly is modeled, it is assumed that the seal wire's initial conditions remain in the same plane as the original model assembly. Therefore, the initial position of the seal wire is parameterized based on the two translational degrees-of-freedom and one rotational degree-of-freedom of its center of gravity, as shown in Fig. 11. The sets of initial conditions for the seal wire are given in Table 2.

The contact pressure distributions on the seal wire are presented with a simplified two-dimensional visualization approach which is inspired by the idea of pressure films. The nodes on the top surface of the seal wire are taken and projected onto a two-dimensional surface to color-code the contact pressures as shown in Fig. 12.

Like the frequency domain model, the results of the time domain models are also nondimensionalized using reference values. The reference value for the power dissipation is selected as the time average of the power dissipation of the model with the original sets of initial conditions of the seal wire. Similarly, the maximum contact pressure value of the model with the original sets of initial conditions of the seal wire is chosen as the reference.

The nondimensionalized contact pressure distributions on the top half surface of the seal wire at static equilibrium calculated for the different sets of the initial conditions of the seal wire are given in Fig. 13. These results are extracted from the simulations where all the contact sets are taken into consideration. The data presented suggest that the maximum contact pressure can vary up to 20% with even the slightest changes in the initial position of the seal wire. The different static equilibrium locations and varying contact pressure distributions prove the nonuniqueness of the static equilibrium for the system.

Figures 14 and 15 display the steady-state displacement responses of the measurement points at the blade tip and on the seal wire for the three sets of initial conditions of the seal wire. Both of these figures show that there are different equilibrium positions for each set of initial positions of the seal wire, which is already supported by the results of the quasi-static simulations. Considering that the system is nonlinear, it is expected that changes in the initial conditions would lead to different solutions. Together with the findings from the contact pressure distribution, the results clearly show that there are multiple equilibrium positions of the system for minor changes in the initial conditions. Similarly, including different sets of contacts in the simulations leads to various equilibrium positions for the components. The involvement of different contact sets changes the nonlinearity of the system. Thus, the components settle to slightly different steady-state solutions. However, the difference in the steady-state static positions can be identified as minor. This outcome is aligned with the results of the frequency domain simulations that show the weak nonlinear characteristic caused by contact with the seal wire. In Fig. 14, the oscillatory displacements are more significant as the measurement point is at the tip of the blade. Also, due to the form of the motion, the effects of the frictional damping of each contact set can be observed at the tip of the blade. Thus, a slight phase difference can be observed between the simulations with different sets of contacts. Whereas, in Fig. 15, the displacement amplitudes and the phase lags between responses are much smaller. Since the displacements are relatively small, the frictional damping levels are not significant. This agrees with the frequency domain simulations, which point out to the weakness of the effects of the nonlinearity due to the contact between the seal wire and the bladed disk.

Furthermore, Fig. 16 shows the steady-state power dissipation by frictional contacts for the three different sets of initial conditions of the seal wire. Like the contact pressure distributions and the steady-state displacements presented earlier, the power dissipation also shows slight variations for the sets of initial conditions of the seal wire. These variations occur due to the changes in the contact nonlinearities in the system. A more significant source of variation in power dissipation is the sets of contacts that are included in the simulations. Individual contacts and all contacts are investigated as separate test cases. Furthermore, the sums of power dissipation of the individual contacts are also presented in the plots. The first thing to emphasize is that summing the power dissipation of the dissipation of individual contacts does not yield similar curves as the dissipation curve of the cases where all contacts are involved. This points out a significantly nonlinear characteristic for the entire system. Besides supporting the claim made when presenting the frequency domain results, the damping contribution of the contact between the seal wire and the blade hub is quite limited. For all the sets of initial conditions, the root is found as the major source of frictional damping and nonlinearity.

Finally, Fig. 17 visualizes the ratio of energy dissipated per period for each set of contacts handled within the scope of this study. Each set of initial conditions for the seal wire is nondimensionalized individually. The energy dissipation of all contacts being considered is selected as the reference value for the nondimensionalization. The data reinforces the idea of the seal wire and lock plate not contributing significantly to frictional energy dissipation in the aircraft engine turbine's bladed disk assembly, and it shows the significance of the initial conditions of the seal wire as the data for the three sets of initial conditions are not notably distinct. The seal wire dissipates slightly less than $10%$ of the energy dissipated when all contacts are taken into account for the analysis. Meanwhile, the lock plate damps roughly 40% more energy compared to the seal wire on average. Furthermore, the primary contribution to damping appears to come from the contacts at the blade root. Over 80% of the total energy dissipation is achieved by the contacts at the root region. Finally, it is important to note that the sums of the energy dissipation fractions for the seal wire, the lock plate, and the root are not equal to 1. The reason for that is already explained earlier.

This study deals with the numerical modeling of a bladed disk assembly of an aircraft engine's turbine, including the nonlinearities due to frictional contacts. The contact between the bladed disk and the seal wire is brought to the spotlight, and the effects of that particular contact are on the dynamics of the whole system.

The results show that the damping contribution of the seal wire is limited compared to the energy dissipated at the interface of other components (blade-disk, blade-lock plate). However, the results also suggest that the seal wire has multiple equilibrium positions depending on its initial conditions in the simulations, which significantly affects contact conditions and hence the overall system response. First, the forced response in the frequency domain is analyzed and shown that the contact with the seal wire alone does not introduce nonlinearity of significant strength. The time domain simulations prove that the static and dynamic equilibrium of the system is dependent on the initial positioning of the seal wire, and there is no single equilibrium state for the system. In connection with that, the contact pressure distributions are not the same in different equilibrium states. The damping contributions of the contacts present in the system are also evaluated, and the contribution of the seal wire is found to be limited compared to the contributions of the root and the lock plate.

This study presents the first steps in the evaluation of the nonlinear behavior of complex mechanical interfaces and highlights the need for further research and method development. For future works, it is necessary to consider two sectors of the system to explore the true contact behavior of the wire within a single vibration cycle. Furthermore, the effect of segmented components (e.g., the lock plate) having a different sector size on the number of blades can be further investigated. In a real engine, it is known that the components wear over time, and thus the contact characteristics can change. The presented methodology can be used for modeling to explore the sensitivity of the system to such evolution of the contacts.

• Rolls-Royce (Innovate UK and the ATI for supporting this research as part of the MALIT programme (Grant No. 113180; Funder ID: 10.13039/501100000767).

• Imperial College Research Computing Service.

Letters and Symbols

A =

excitation amplitude

a_i =

excitation amplitude

b_i =

excitation amplitude

C =

damping matrix

f =

force/load vector

F =

force (N)

K =

stiffness matrix

M =

mass matrix

n =

number of sectors

N =

number of harmonics

P =

pressure (N m⁻²)

r =

spatial coordinate (mm)

t =

time (s)

u =

displacement (vector) (mm)

$u˙$ =

velocity (vector) $(mm s−1)$

$u¨$ =

acceleration (vector) $(mm s−2)$

x =

spatial coordinate (mm)

y =

spatial coordinate (mm)

z =

spatial coordinate (mm)

θ =

spatial coordinate (deg)

μ =

friction coefficient

$ϕ$ =

mode shape

ω =

(nondimensional) frequency

Superscripts and Subscripts

exc =

excitation

f =

friction

n =

numbering index

nat =

natural

norm =

normal

p =

pressure side

rel =

relative

s =

suction side

tang =

tangential

0 =

static