Numerical and Experimental Study of Friction Damping in Blade Attachments of Rotating Bladed Disks

In order to mitigate high cycle fatigue risks in bladed disks, the prediction of the vibration levels early in the design process is important. Therefore, the different sources of damping need to be modeled accurately. In this paper the impact of friction in blade attachments on forced response is investigated both numerically and experimentally. An efficient multiharmonic balance method is proposed in order to compute the forced response of bladed disks with contact and friction nonlinearities in blade roots. For experimental validation purposes, a rotating bladed disk was tested in a vacuum chamber, with excitation being provided by piezoelectric actuators. A model of the rig was built and numerical results were obtained with a normal load dependent coefficient of friction and a constant material damping ratio. Nonlinear behavior observed experimentally at resonances was well reproduced and an acceptable correlation was found with experimental resonant frequencies, amplitudes, and amount of damping throughout the spinning speed and excitation level range. The proposed numerical method can therefore serve to enhance the prediction of the alternating stresses in bladed disk assemblies.


INTRODUCTION
In gas turbine engines, the nonhomogeneous air flow causes high-frequency vibrations of rotating bladed disks.Excessive amplitudes can be responsible for high cycle fatigue failures (Srinivasan [1]).Thus there is a need for tools able to estimate with sufficient accuracy the dynamic stresses, especially near resonances.The paper of Seinturier et al. [2] presents an industrial application of many recent advances in the computation of the forced response and points out some of the remaining difficulties in this domain.
The amount of damping in the system has a direct influence on the vibration levels, but it is still very difficult to estimate it precisely.Three main phenomena are responsible for energy dissipation in bladed disks: material hysteresis, aerodynamic damping, and dry friction in joints.In his review article, Srinivasan [1] provides a comparison of typical Q factors due to the main sources of damping.Material damping is very small for the alloys used and is sometimes neglected in the computations.The level of aerodynamic damping is generally higher but its prediction is difficult and requires computational fluid dynamics (CFD) analyses that take into account blade motion (see Seinturier et al. [2], Berthillier et al. [3]).Finally, energy dissipations by friction can occur in many locations, for instance between shrouded blades or in the joints between the blades and the disk.Moreover dry friction dampers are widely used to decrease the vibration level of dangerous resonances.The efficiency of such devices was studied experimentally on a rotating bladed disk in Tokar' et al. [4] and theoretically, for instance, in the early work of Griffin [5].The harmonic balance method (HBM) is the most widely used technique for the computation of the steady-state response with friction damping.Many efforts have been done in the recent years to propose numerically efficient implementations of this method (Berthillier et al. [6], Petrov and Ewins [7], Poudou and Pierre [8], Nacivet et al. [9]).Among them, the dynamic lagrangian frequency time (DLFT) method, proposed by Nacivet et al. [9], is characterized by its ability to handle directly the nonsmooth contact laws of unilateral contact and friction.
Experimental studies of blade root damping include the paper of Rao et al. [10].Their test rig features a beam with "T" roots at both ends.A special device based on the thermal expansion phenomenon is used to pull apart the beam, thus inducing traction stresses that mimic centrifugal stresses.The results confirmed the expected decrease of damping as the traction force was increased.Experiments on real rotating components have been reported by Kielb and Abhari [11] and Tokar' et al. [4].In the former study, the same turbine bladed disk was tested in vacuum with piezoelectric excitation as well as in a shock tunnel.The comparison between the results revealed that for this particular case the aerodynamic damping was largely predominant over the structural damping.The structural damping, attributed essentially to friction in the firtree attachments, was found to be inversely proportional to the square of the rotational speed.Tokar' et al. [4] measured friction damping in two different blade attachments (dovetail and streamline joints) and provided solutions to increase energy dissipation.It was shown in particular that the damping ability of dovetail joints could be improved by the use of antifriction coatings.
The modeling of blade root friction presents some difficulties.For one blade, the contact interfaces include two or more flanks which are heavily loaded because of rotation.On one flank, the contact pressure is neither constant in space nor in time.High-stress gradients have to be expected near the edges of the contact zone (Sinclair et al. [12], Beisheim and Sinclair [13]).Additionally, the normal pressure varies with time according to the dynamic forcing of the blade and separation of part of the interface can occur.It is therefore necessary to include a normal degree of freedom in the contact model.The relative displacements in the tangential directions are of very small amplitude and a microslip behavior is expected.Microslip models have already been used in the analysis of blades with friction dampers (Menq et al. [14], Sanliturk et al. [15]), but such models are not applicable in the case of blade root friction because the normal loads are not known in advance.The first numerical studies devoted to the effects of root friction on the forced response of blades have been published recently (Charleux et al. [16], Petrov and Ewins [17]).In both cases three-dimensional finite element models are used and the laws of unilateral contact and friction are enforced within a certain number of contact elements distributed over the root flanks, thus allowing for partial slip and partial separation.
In this paper, a DLFT method is used to compute the forced response of a bladed disk with friction in blade roots.The main objective of the study is to verify if the method is able to correctly predict the nonlinear dynamic behavior of the system.To this end, an experimental validation is proposed.In the next section, the formulation employed to treat the contact and friction nonlinearties is described.The test rig is presented in Section 3 and the corresponding numerical model is detailed in Section 4. Finally in Section 5 the numerical results are compared to the experimental results and the accuracy of the proposed numerical method is discussed.

NUMERICAL METHOD
The DLFT method detailed below allows for the computation of the steady-state response of structures with contact and friction.It is a harmonic balance formulation and therefore displacements and forces are assumed to be periodic.Moreover, our presentation in this paper is restricted to small deformations and small displacements, which means in particular that the sliding distances are small.We further assume that meshes are matching on the contact zones, which authorizes the use of node-to-node contact elements.For each contact element, relative displacements and contact forces are computed during resolution according to the Coulomb law of friction and the unilateral contact law.The input data required comprise the mass, stiffness, and material damping matrices as well as the vector of external forces for each solid.

Equations of motion
Let us consider the case of two flexible solids in contact with friction, for instance a disk and a blade.The equations of motion obtained after discretization can be written for each body l as M l , C l , and K l are mass, material damping, and stiffness matrices for the body l. u l , ul , ül are the displacement, velocity, and acceleration vectors, respectively.f l ex stands for the external forces.The vector of unknown nonlinear contact forces, f l c , includes normal and tangential components.The contact zones are defined by the node-to-node contact elements chosen by the user.Rotation of the local axes may be necessary in order to have one degree of freedom (dof) in the contact normal direction and two dofs in the tangential plane for each contact node.
Assuming that the steady-state motion is periodic, a Galerkin procedure is performed in order to formulate the equations in the frequency domain.Displacements and forces are then represented by multiharmonic vectors which are formed by the Fourier coefficients of the harmonics retained for the resolution.At this point, the size of the problem can be reduced by performing two exact reductions in the frequency domain (see Nacivet et al. [9]).In the first one only the degrees of freedom involved in the contact elements are retained.Further factor two reduction may be obtained by presenting the problem in terms of relative displacements.The equations of motion finally take the following form: where u r , λ, and f r are the multiharmonic vectors of relative displacements, Lagrange multipliers, and reduced external forces, respectively.Z r represents the reduced dynamic stiffness matrix.The Lagrange multipliers are equal to the unknown contact forces.A nonlinear solver is used to determine the zero of the following function: (3)

Computation of contact forces
For the computation of function f it is necessary to determine the contact forces λ for a given vector of relative displacements u r .Here, λ is formulated as a penalization of the equations of motion in the frequency domain where ε is a penalty coefficient and x r is a new vector of relative displacements.It will be seen in the next section that the pair ( λ, x r ) is determined through an alternating frequency time (AFT) method so that the contact and friction conditions are fulfilled in the time domain.Equation (3) reduces to f ( u r ) = ε( u r − x r ).The convergence thus ensures that (2) is verified and that contact and friction are taken into account since at this point u r = x r .The value of the penalty coefficient should be chosen so as to balance the contributions of the equation of motion to be solved with the contributions of the contact constraints.Taking the spectral radius of the dynamic stiffness matrix ρ(Z r ) generally gives good results.

Prediction-correction in the time domain
The contact forces are computed in the time domain where transition criteria between the three possible states (stick, slip, and separation) are easily expressed.Equation( 4) can be reformulated as where The period is split into N time steps.Let us focus now on one given contact element.For this contact element, the counterparts of λ, λ u , λ x , x r in the time domain are {λ n } n=1,...,N , {λ n u } n=1,...,N , {λ n x } n=1,...,N , {x n r } n=1,...,N , respectively.At each iteration of the nonlinear solver, λ u is calculated and the set {λ n u } n=1,...,N is obtained with an inverse FFT algorithm.A prediction-correction strategy is used at each time increment t n to compute the contact forces λ n .A predicted contact force λ n pre is calculated assuming that the vector of relative tangential displacements x T r remains the same between t n−1 and t n , and that the normal relative motion x N r is zero.This leads to Exponents N and T refer to the normal component and to the vector of tangential components, respectively.Vectors (i.e., λ n,T u ) are boldfaced while scalars (i.e., λ n,N u ) are not.The next step is to correct the contact force so that the complete contact law is verified.This is done by computing λ n x and the corrected contact force is then given by In order to calculate λ n x , three cases must be distinguished.
(i) Separation: λ n,N pre ≥ 0. The predicted force is a traction force which means that the contact is lost.The contact force has to be set to zero which requires (ii) Stick: The normal relative displacement is zero and the tangential relative displacement is constant.This leads to Again, there is no normal relative displacement.The correction is made assuming that the normal force does not change and that the tangential contact force has the same direction as the tangential predicted force.The tangential relative speed, implicitly defined by must be in the direction of the predicted tangential force, while the contact force must be on the Coulomb cone.Therefore, λ n x is given by Once λ n x has been determined, the contact force λ n is obtained by (8).This operation is repeated along the period for each contact element.Then, an FFT algorithm is used to obtain the expression of the contact forces in the frequency domain and the residual f ( u r ) from ( 3) is returned to the nonlinear solver.
The numerical procedure used is based on a specific estimation of forces λ (( 5)-( 12)) which has been summarized by a flow chart shown in Figure 1.

Test rig
The experimental data were obtained with a test facility designed to gain a better understanding of friction phenomena in bladed disks.In the present study, the only friction zones considered are the dovetail joints between the disk and the blades.A schematic view of the test bench used can be seen in Figure 2. The rotating disk supports four blades and is placed in a vacuum chamber so as to minimize the effects of aerodynamic forces.The disk is rotating on a hollow shaft supported by two ball bearings and driven by an electric motor.The blades considered here are the compressor blades of an aircraft engine.Two of them are shown in Figure 3.

Instrumentation
The four blades are placed around the disk as two diametrically opposite pairs as seen in the mesh of Figure 6.Two adjacent blades are equipped with piezoelectric actuators and strain gauges, while the other two remain bare.The two equipped blades, herein called blade 1 and blade 2, are shown in Figure 3.The lead zirconate titanate (PZT) ceramics are 1 cm square flat layers with a thickness of 1 mm.An insulating epoxy adhesive is used to bond them to the airfoils and to adapt to the slightly curved surface.Two ceramics are used for each blade, one on the pressure side and one on the suction side.A parallel wiring together with a suitable orientation of polarization directions allow them to work out of phase and to generate a flexion motion with a common tension signal.The PZT layers are placed on high-strain regions of the first bending mode so as to maximize the effectiveness of the excitation.A slip ring, located at the extremity of the shaft, is used to transmit excitation and measurement signals.

Measurements
The first step of the experimental procedure is to pump the air out of the vacuum chamber.The results presented in this paper were obtained with a stabilized pressure of 20 mbar.Then the bladed disk is set into rotation.Five rotational speeds were studied from 1000 rpm to 5000 rpm.In order to study the frequency response around the modes of interest, a swept sine excitation was provided to the actuators with a sufficiently slow sweep rate so as to avoid artificial distortions of responses in case of low damping levels.For each spinning speed, voltage levels up to 100 V can be applied to the piezoelectric ceramics.In this study, the retained levels were 10 V, 20 V, 40 V, and 80 V. Frequency response functions were provided by measuring only the fundamental component and the half-power bandwidth method was used on an isolated resonance peak to evaluate the global amount of damping in the system.

Experimental results
In the present paper, attention is focused on a narrow frequency range which includes two resonance peaks.An example of frequency responses for blades 1 and 2, to an excitation on blade 2, is shown in Figure 4 for a given excitation level of 80 Volts and for a 3000 rpm rotating speed.For these two modes, the strain energy is mainly localized in a single blade which vibrates according to its first bending mode.Such a phenomenon is due to mistuning between the blades.In the first resonance (mode 1), blade 1 has the greatest vibration amplitude among the four blades.The amplitude of blade 2 is much lower and blade 1 and 2 vibrate in phase.In the second resonance (mode 2) at a higher normalized frequency (value around 4), strain energy is mainly localized in blade 2 with blade 1 and 2 vibrating out of phase.Note that strains are in phase (resp., out of phase) while blades vibrate out of phase (resp., in phase) as far as strain gauges are not located on the same side for each blade for this measurement.
Figure 5 shows a set of experimental frequency responses obtained at 2000 rpm.The excitation was applied on blade 2 only and the response was measured on the same blade.The resonance of mode 2 is clearly visible while mode 1 is not detected.It can be noted that the peak is shifted to the left as the excitation voltage is increased, which is evidence of a nonlinearity.This effect is more clearly seen in Figure 5(b), where the strain amplitude was divided by the excitation level.Had the system been linear, all the curves would have coincided on this graph.The behavior is linear for low-vibration amplitudes but nonlinear in the proximity of the resonance.It is presumed that as the amplitude increases, the amount of slip in blade attachments increases resulting in a higher damping level.Further results will be presented in Section 5 and compared to numerical results.

Finite element model
A numerical model of the bladed disk was created.The four blades and the disk were meshed using second-order tetrahedral elements (i.e., 10 nodes per element).The mesh can be seen in Figure 6.The shaft was modeled with beam elements and rigid body elements were used to link the shaft and the three-dimensional mesh of the disk.The two bearings were modeled with axial and radial linear springs.The commercial FE code Samcef was used to carry out the preliminary linear calculations.The model was divided into five substructures: one for each blade and one for the disk-shaft assembly.The number of degrees of freedom (dofs) was reduced in accordance to the Craig & Bampton component modes synthesis method [18].This method was applied using a static mode for each contact dof plus additional ones for piezoelectric actuators and stain gauges dofs and finally the reduction basis is completed by dynamic clamped modes.The reduced basis representing the diskshaft assembly encompasses the disk contact nodes and 18 modal dofs.In the present study, 8 × 3 = 24 contact elements are defined on each flank of the blade roots, as shown in Figure 7.The natural frequencies obtained with the reduced model were compared with those obtained with the full-size model and the accuracy of the reduced model was found to be acceptable since the relative discrepancy remained under 1% for the first 11 modes, including the two studied modes.At the present stage, rotation is taken into account but the frictional contact nonlinearity is not considered: the contact interfaces are assumed to be perfectly welded.Under this assumption the mass matrix M and the tangent stiffness matrix K t for each substructure are computed and retrieved.The tangent stiffness matrix can be decomposed as where K is the structural stiffness and K c is the centrifugal stiffness.The term −K c is a stiffness correction due to the expression of the equations of motion in the rotating frame and is responsible for a softening effect.K g is the geometric stiffness matrix which represents the stiffening of the substructure under the rotation-induced stresses.For each substructure, the equations to be solved are written in the reduced basis as where f ex stands for the external forces, f c represents the contact forces, and C is a Rayleigh damping matrix.The gyroscopic matrix is not taken into account.Indeed, the results shown in this paper are computed for particular modes of the structure where there is no shaft bending.In this case, it was found that the forced response computed with the gyroscopic matrix is only very slightly modified and that the natural frequencies were shifted by less than 0.01%.These results are in accordance with those of Sinclair et al. [12], who found with a Ritz method that the modes of a rotating radial beam were not affected by the Coriolis acceleration.

Piezoelectric excitation modeling
The piezoelectric actuators were modeled with the simple unidirectional model shown in Figure 8.In this model, the thin adhesive film is not taken into account and the load is assumed to be transferred to the structure via the two extreme points or nodes A and B. The force produced by the actuator when a voltage U a is applied to the electrodes can be written as (Jemai [19]) where d 31 is the piezoelectric charge constant of the PZT ceramic, S 1 is the area of the cross section, Y 11 is the Young modulus of the ceramic in direction 1. e, l 1 , and Δl 1 are, respectively, the thickness, the length, and the elongation in direction 1 of the ceramic layer.The stiffness of the piezoelectric actuator between points A and B is denoted by K p : If K s is the local stiffness of the structure in direction 1 at the location of the actuator, then F 1 can also be expressed as Combining ( 15), (16), and (17), a simple expression of the excitation force as a function of the applied voltage is obtained: This convenient proportional relationship was used to model blade excitation.For each piezoelectric actuator, two nodes located at the position of the ceramic were kept in the Craig-Bampton reduced basis of the blade.The stiffness of the structure K s was obtained by a simple FE computation: opposing static forces were applied in direction 1 at nodes A and B and the calculated relative displacement served to evaluate K s .The validity of this modeling was verified a posteriori by the comparison of the experimental and numerical results (see Section 5).

Frequency mistuning measurement and model updating
Mistuning between the four blades is unavoidable because of manufacturing tolerances and material discrepancies.Furthermore, in the presented experiment, two blades are equipped with excitation and measurement devices (blades 1 and 2), while the other two are not.It is important to include  mistuning in the numerical calculations, since it is known to substantially affect the modes and the forced response of bladed disks (Srinavasan [1], Seinturier et al. [2]).
In order to quantify the frequency mistuning, the four blades were tested on the small bench shown in Figure 9.The blades were excited by a hammer blow while a laser vibrometer was used to measure the impulse responses.Piezo actuators and strain gauges were not used in this test.As the frequency depends on boundary conditions, a calibrated tightening torque of 10 Nm was imposed in each case by means of a M12 (12 mm diameter) screw that squeezes the blade root, which resulted in a good reproducibility of the measurement.The mean frequencies gathered in Table 1 were normalized because of some confidentiality requirements.These results suggest that the piezoelectric ceramics are responsible for a slight stiffening effect.
This data was used to update the Young modulus of each blade finite element model.The densities were also adjusted after weighing the blades.Finally, impulse responses of the shaft-disk assembly without the blades were measured at rest and the found natural frequencies were used to update the corresponding model.In particular, the bearing spring constants had to be adjusted.

Nonlinear analysis
The steady-state solution of ( 14) was sought with the DLFT method presented in Section 2. Three-dimensional node-tonode contact elements were used.They take into account the variable normal load at the blade root and the 2D stick-slip motion.There are 192 of them which represent 384 contact nodes and 1152 nonlinear degrees of freedom.With the linear dofs added, the total size of the problem is 1290.Using the reductions in the frequency domain mentioned in Section 2, the size of the system to be solved by the nonlinear solver is actually only 576.
In order to facilitate convergence at the first frequency step of the forced response, the static equilibrium of the bladed disk subjected to centrifugal loads only is computed.The relative displacements obtained in the blade root are then used as a starting point for the calculation of the dynamic response.
The steady-state responses presented in this paper were determined with the first three harmonics and a penalty coefficient ε of 2.5 × 10 4 .The maximum amplitudes were compared with those obtained with one harmonic.The greatest variation found was 10% for a rotational speed of 1000 rpm and an excitation level of 80 V. Computations with five harmonics were also performed, but the three harmonics results were improved by less than 1%.Beyond, the calculations became very expensive and brought no further significant improvements.So retaining three harmonics offers a good compromise between speed and accuracy.For the studied case, converged results are obtained if the penalty coefficient ε is chosen from 5×10 3 to 10 5 .Variations of less than 1% in maximum amplitudes are encountered if ε is varied within this range.Moreover, the spectral radius of the dynamic stiffness matrix ρ(Z r ) is approximately 2.5 × 10 4 at the studied resonances.This value is actually located in the above-mentioned range and was chosen for the DLFT computations.

Nonfrictional damping
In the studied experiment, damping is not due solely to blade root friction.Energy is also dissipated in materials and in the instrumentation.Moreover, a certain level of aerodynamic damping might remain.All these other sources of damping were represented by an equivalent viscous damping ratio ζ.In practice, they were introduced in the reduced model in the form of a Rayleigh damping matrix: where α and β are coefficients that were adjusted so as to obtain the desired ζ for the studied resonance.The same value of ζ was used for all rotational speeds and all excitation levels.It was therefore assumed that any change in the global amount of damping was due to dry friction in dovetail joints.
In order to choose the value of ζ we considered a case where energy dissipation by dry friction in root joints is minimal, that is, at the highest possible spinning speed (5000 rpm) and at the lowest excitation level.In these conditions, a damping ratio of 5.5 × 10 −4 was identified in the response.This value was assigned to ζ.

Friction parameters
In the blades studied, an antifriction coating is used to protect the blade roots from wear and fretting fatigue.Under these conditions, the coefficient of friction is expected to be low (about 0.1) and to decrease as the normal load is increased (Rabinowicz [20], Ibrahim [21]).The friction coefficient was therefore assumed to be governed by the equation given in Ibrahim [21]: where γ is constant, N is the normal load, and ν has a positive value.
The practical implementation of such a relation for blade-disk contact requires some precision.It is recalled that the normal contact pressure is not uniform on the flanks and has a time-varying component due to vibration.Hence, according to (20), the coefficient of friction should vary too.This was not taken into account however and the coefficient of friction μ was assumed to be the same for all contact elements and constant in time.For each rotational speed of interest, the value of μ was chosen thanks to a simplified 2D static analysis.As can be seen in Figure 10, the centrifugal force causes normal and tangential static contact forces N and T because of the flank angle α.Assuming that rotation brings about macroslip in the blade root-which was confirmed later by simulations with the 3D FE model-then As the experiments and simulations were conducted at five rotational speeds from 1000 rpm to 5000 rpm, these are denoted by Ω j = 1000 j, for j from 1 to 5. N j and μ j are the corresponding normal force and the coefficient of friction, respectively.As the centrifugal force is proportional to the square of the rotational speed, the static equilibrium is given by where c is constant.If we define a normal force ratio by then ( 20) can be reformulated as With ( 23) and ( 24) and using the j = 1 case to express constant c, (22) becomes The parameters involved in this equation have the following values: (i) flank angle α = 50 • , (ii) exponent of normal load dependence ν = 0.34, (iii) coefficient of friction at 1000 rpm μ 1 = 0.15 (adjusted numerical value as explained below).Exponent ν was chosen using Seitzman et al. [22], where molybdenum disulfide coatings for titanium alloys are studied.The solid lubricant used on the blade roots being of a similar type, a similar value of ν was assumed.Equation ( 25) is nonlinear in ψ j and is solved with a Newton method for j from 2 to 5 for a given μ 1 .Therefore, the simplified 2D static analysis permits choosing a single unknown parameter (here μ 1 ) to be adjusted.The value of this latter is updated so that frequency responses obtained by means of the DLFT analysis at several rotating speed match with experimental ones.The corresponding coefficients of friction obtained by (24) are listed in Table 2.
The elastoplastic deformations of asperities were not taken into account in the model because of lack of precise information about the local characteristics of the contact studied.This means that when a contact element is in the sticking state, no relative displacement is permitted.

RESULTS
In this last section of the paper the experimental and numerical frequency responses are compared.Attention is mainly focused on a particular mode of the bladed disk (mode 2 described in Subsection 3.4).A comparison between experimental and numerical frequency responses is presented in Figure 11.Figures 11(a), 11(c) are the results for blade 1, and Figures 11(b), 11(d) represent those for blade 2. A good agreement between numerical and experimental results is observed both in amplitude and phase curves.Resonance of mode 1 is visible in the response of blade 1.By dividing the amplitudes by the excitation level, as done in Figure 5, a linear behavior could be observed at this resonance.This was also the case in the simulations since all the contact elements remained in stick state.Hence, the observed discrepancy most likely means that material damping was underestimated for this particular rotational speed.For mode 2, the experimental peaks are slightly shifted to the left as U a is increased from 10 V to 80 V, which is not well reproduced by the simulation.Among the suspected causes are the level of discretization and the absence of tangential contact stiffnesses in the model.

Resonant frequencies
Figure 12 shows the evolution of the resonant frequency for mode 2 over the speed range.These results are obtained for the lowest excitation level, which minimizes the frequency shift due to the nonlinearity.The numerically predicted centrifugal stiffening is in good agreement with the experiment.But one can notice that the two curves do not have exactly the same shape.The experimental curve is almost straight, while the numerical curve bends a little.One possible explanation for this may be that the model does not include contact stiffnesses.Indeed, experiments reported in Ferrero et al. [23] and Crassous et al. [24] have shown that the tangential contact stiffness increases with normal load.Including such normal load dependent contact stiffnesses in the model would straighten the numerical curve.

Amplitudes and damping
The maximum amplitude at which the blades vibrate is a key information for the designer interested in preventing high cycle fatigue problems.Dynamic stresses together with static stresses are necessary to verify if the blades remain sufficiently under their endurance limit.
Figure 13 shows the maximum vibration levels of the bladed disk for the resonant response of mode 2. First of all, it is important to stress that for a speed of 5000 rpm and an excitation of 10 V, the amplitudes obtained numerically are very close to the experimental ones.In this case, energy dissipation in blade roots is minimal and behavior is almost linear.Furthermore, one can note that in Figure 14(a) the measured and numerical damping ratios are very close.So with the same amount of damping, both computations and experiments give the same amplitude, demonstrating the acceptable accuracy of the simple model used for piezoelectric actuation.Figure 13 also reveals that, to a good extent, the numerical model reproduces the main trends observed in the experiments.Nonetheless, one point of interest seen in (b) is that the slope of the strain versus speed curves obtained with the numerical model can be negative for high speeds and low amplitudes, whereas the experimental curves always exhibit a positive slope.Kielb and Abhari [11] tested a rotating bladed disk in a similar facility at speeds up to 20 000 rpm and also found ever increasing amplitudes.The fact that structural damping model does not take into account the speed and prestress influence due to rotating speed can explain the discrepancies between experimental and experimental slopes.
International Journal of Rotating Machinery  For the highest recorded amplitude (blade 2, 5000 rpm, 80 V), the maximum computed alternate Von Mises stress in the airfoil part of the blade is approximately 20 MPa.This is more than twenty times below the endurance limit.Nevertheless, the excitation level was sufficient to produce a significant nonlinearity for each speed.
The results presented in Figure 14 confirm that the total damping decreases as the spinning speed increases.The damping levels are of the same order of magnitude as those found by Tokar' et al. [4] and Kielb and Abhari [11].Experimental curves are rather well reproduced by simulation even if damping is quite underestimated for low-vibration levels (curves (a) and (b)).The benefit of taking into account friction in blade root is obvious in this case, as linear calculations would have produced horizontal lines with a damping ratio of 5.5 × 10 −4 .
Among the possible causes of the observed discrepancies is the assumption made in the model that friction in dovetail joints is the only source of nonlinearity.For instance as said before, the possible variation of material damping with vibration amplitude was not taken into account.Moreover, as the rotational speed increases, centrifugal stiffening could be 1000 2000 3000 4000 5000 Rotational speed (rpm) 3.9 responsible for a decrease in material damping.Smith and Wereley reported such a phenomenon in the case of rotating composite beams (Smith and Wereley [25]), but no similar experimental work on titanium alloys has been found in the literature.Another source of error relates to the number of contact elements used for simulation.Sinclair et al. [12] showed that very refined grids were necessary to obtain converged contact stresses in dovetail joints, but such discretizations are not compatible with dynamic computations.In this study, a compromise was found and 24 contact elements per flank were used.This is insufficient to accurately predict the contact stress distribution and therefore also brings about errors in the estimation of slip regions and slip amplitudes.

CONCLUSIONS
An attempt to predict the nonlinear forced response of bladed disks with friction in blade attachments has been presented.A DLFT method was used to compute the steadystate response taking into account nonsmooth Coulomb friction and unilateral contact laws.To validate this method, experimental testing of a bladed disk rotating in a vacuum chamber was conducted.With harmonic excitations, the measured frequency response exhibited a nonlinear behavior in the vicinity of resonances, attributed to blade root friction.A three-dimensional finite element model of the rig was constructed and updated taking into account mistuning of the blades.
The simulations were carried out with a normal load dependent coefficient of friction and a constant material damping.Under these assumptions, experimental damping levels have been reproduced with reasonable accuracy for different rotational speeds and different excitation levels.The method can therefore be used to better predict the alternating stresses in bladed disk assemblies.The observed discrepancies between numerical and experimental results can be due in particular to the number of contact elements retained which was insufficient to accurately predict the distribution of contact stress.Contact stiffnesses, that may vary with normal load, were not included in the model presented in this paper and could also play an important role.There is a need for more experimental data to enhance the level of modeling and the precision of numerical results.
In particular, material damping should be determined alone as well as its variation with amplitude and centrifugal load.

Figure 1 :
Figure 1: Flow chart, computation of the function for the nonlinear solver.

Figure 2 :Figure 3 :
Figure 2: Schematic view of the test bench.

Figure 4 :
Figure 4: Sample of experimental frequency response to an excitation on blade 2 at 3000 rpm.

Figure 11 :
Figure 11: Numerical and experimental frequency responses at 3000 rpm.

Figure 12 :
Figure 12: Numerical and experimental resonant frequencies for U a = 10 V.

Table 1 :
Measured frequencies of the first bending mode and deviation from average value.

Table 2 :
Coefficients of friction used in the model.