Simulation of Transient Nonlinear Friction-Induced Vibrations Using Complex Interface Modes : Application to the Prediction of Squeal Events

During the past decades, the problem of friction-induced vibration and noise has been the subject of a huge amount of works. Various numerical simulations with finite elements models have been largely investigated to predict squeal events. Although a nonlinear analysis ismore predictive thanComplex Eigenvalues Analysis, one of themain drawbacks of the time analysis is the need of large computational efforts. In view of the complexity of the subject, this approach appears still computationally too expensive to be used in industry for finite element models. In this study, the potential of a new reduced model based on a double modal synthesis (i.e., a classical modal reduction via Craig and Bampton plus a condensation at the frictional interface based on complex modes) for the prediction of self-excited vibrations of brake squeal is discussed.The effectiveness of the proposed modal reduction is tested on a finite element model of a simplified brake system. It will be shown that numerical results of times analysis by applying the proposed reduction correlate well with those of the nonlinear analysis based on a reference model, hence demonstrating the potential of using adapted modal reductions to predict the squeal propensity and to estimate self-excited vibrations and noise.


Introduction
Brake squeal is the result of friction-generated vibration and noise.It corresponds to a very common problem in industrial application.Thereby friction-induced vibration has become one of the most costly problems and important concerns in automotive engineering.Indeed this instability which results in a harsh noise, especially in the frequency range between 0.1 kHz and 15 kHz, is annoying to passengers' hearing which leads to consistent customer complaints.Comprehensive overviews on friction-induced vibration and much classic understanding of the mechanism of the disc brake squeal phenomenon can be found in [1][2][3][4][5][6][7].
Finite element models and numerical simulations have been extensively investigated to detect brake squeal instabilities and to predict amplitudes during squeal events.
Nowadays three kinds of analysis can be performed in the design and optimization process of brake systems in order to avoid or at least to reduce the squeal occurrence: eigenvalue analysis to detect frequency instabilities, time analysis, and nonlinear methods [8] to determine self-excited vibrations during the squeal event and acoustic analysis to estimate acoustic radiation.While the Complex Eigenvalue Analysis (CEA) represents a fundamental tool that is widely used in industry to predict unstable vibration modes, this approach may lead to under-or overestimation of unstable modes [9][10][11][12][13][14], and more realistic results can be achieved via the time nonlinear analysis [15].While a brake squeal is characterized primarily as an acoustic phenomenon, the prediction of the acoustic response during squeal events has been usually ignored.There have been only few attempts to calculate and analyze the acoustic radiation of brake squeal [16][17][18][19].
Indeed the calculation of acoustic events necessarily involves an estimate of the nonlinear amplitudes and velocities of brake squeal via a time analysis.One of the major drawbacks of the nonlinear time simulation, nevertheless, is its high computational cost which explains that this approach based on time analysis has mostly been neglected in the past.As a consequence, it is difficult to calculate and track the acoustic radiation during squeal events.So the development of numerical methods to facilitate the prediction of amplitude vibrations during squeal events is one of the crucial issues and challenges.
As previously explained working in the time domain and performing simulations on a full finite element model require prohibitive computational costs.In order to overcome these issues, one solution is to perform model reduction.The use of appropriate reduction bases for the prediction of squeal frequencies and self-excited vibrations during the squeal event and its associated numerical implementation are complex tasks.Recently some efforts have been spent to overcome this problem by considering the potential of using reduction techniques of finite element models.Vermot des Roches [20] proposed a reduction method adapted to large models with local nonlinearities.The reduction method allows keeping the exact system real modes.Even if the computations remain costly, the numerical results were very encouraging.These results contributed greatly to new developments on the use of different modal reductions for the community of brake squeal.Then Loyer et al. [21] proposed using spatial model reductions for friction destabilized systems in order to study the nonlinear behaviors of an elastic layer.They concluded that the inclusion of static modes in the reduction basis is mandatory to approximate the self-excited vibrations of the mechanical system.They also indicated that the model reduction may lead to an alteration of the transient evolution, often leading to an advance phenomenon.
Recently, Besset and Sinou [22] proposed a new efficient reduction method based on a generalized double modal synthesis method.This numerical technique can significantly reduce the size of the original finite element model of the brake system by combining a classical modal reduction and a condensation at the frictional interface via a reduced complex mode basis.One of the main limitations of this previous study is the restriction of the problem on the prediction of frequency instabilities.One objective of this paper is, therefore, to address the use of the generalized double modal synthesis method in conducting rigorous numerical simulations of brake squeal for both the detection of squeal frequencies and the calculation of self-excited vibrations during squeal events.This paper consists of four major sections: introduction, finite element brake model and reduction method adapted to numerical brake system, numerical simulations and results based on the Complex Eigenvalue Analysis and the transient analysis, and conclusions.

Numerical Brake Model
In this section, the finite element model of the brake system is first presented.Then the double modal reduction based on the Craig and Bampton reduction and the interface reduction via complex eigenmodes is briefly introduced.

Finite Element Model.
A simplified model of a brake system (i.e., "pad-on-disc": a steel pad in sliding contact with an annular disc) is depicted in Figure 1.The finite element model used in this paper was explained previously in detail in [22].At the upper surface, the pad is constrained in the tangential and radial directions (with no rotational degrees of freedom).The inner surface of the disc is clamped.The finite frictional surface-to-surface contact is composed of 220 contact nodes.It can be noted that a convergence study versus the interface refinement (i.e., the number of contact points at the frictional interface) has been carried on in order to ensure the independence of the simulation results versus the interface meshing.Details on the geometric and material properties of the finite element model are given in Table 1.
The nonlinear equation of the numerical model can be written in the following form: where M, K, and C are the classical mass, stiffness, and damping matrices, respectively.A modal damping (with the damping coefficients   fixed at 0.1% for each mode) is considered.X is the generalized displacement vector while the dot denotes derivative with respect to time.F ext represents the hydraulic pressure that is directly applied to the backplate of the pad.The pressure force applied on the upper surface of the pad is equal to 2.10 7 Pa in our case.A velocity field corresponding to the disc rotation is imposed on the disc surface.F NL contains contributions from both the contact nonlinear forces and frictional forces at the pad/disc interface with a classical Coulombs law.Considering previous experimental investigations [23], we assume that the nonlinearities at the friction interface between the pad and the disc are due to the cubic nonlinear terms and the possible loss of contact.The components of the nonlinear contact force vectors denoted by   contact,pad and   contact,disc (  friction,pad and   friction,disc , resp.) along the normal direction (along the tangential direction, resp.) are given at th contact node by where   corresponds to the relative displacement between the pad and the disc for the th contact node.  and  NL define the linear and nonlinear stiffness at the friction interface.Values of these linear and nonlinear stiffness are fixed at 9.10 2 N⋅m −1 and 4.10 4 N⋅m −3 , respectively.V  is the relative velocity between the pad and the disc and  defines the friction coefficient that is assumed to be constant.

Reduction Method Adapted to Numerical Brake Model.
In order to reduce the finite element model and more specifically  to condense the pad and disc, the Craig and Bampton technique is first used.This reduction technique is one of the most classical processes to reduce the finite element by retaining all the contact nodes at the pad/disc interface and generalized degrees of freedom.The Craig and Bampton technique [24] consists in building a reduced model using fixed interface modes and constraint modes.This general process with the application on the brake model under study was explained previously in detail in [22].Convergence of Craig and Bampton method for a similar finite element brake system has also been discussed in [25].Secondly, a reduction method using complex interface modes is applied in order to condense the pad/disc section (i.e., the remaining physical interface degrees of freedom that are kept unreduced via the Craig and Bampton method).This reduction method is now briefly explained.
In the process, the mass and stiffness matrices (denoted by M and K, resp.) are firstly implicitly reduced.They are obtained by using the transfer matrix T 2 : where (u   , u   ) and (  ,   ) define the interface degrees of freedom and the generalized degrees of freedom, respectively (see [22] for more details).
By considering the reduced mass and stiffness matrices via the Craig and Bampton method (denoted by M and K, resp.) and the transfer matrix T 2 , we obtain For the reader comprehension, the reduced stiffness matrix K via the Craig and Bampton method contains not only the structural stiffness matrix K but also the Jacobian matrix of the linearized expression of the nonlinear contact and friction forces F NL at the disc/pad interfaces around the nonlinear equilibrium (see [22] for more details).
Then, right and left eigenvectors (denoted by X   and X   ) are calculated by resolving the following equations: By introducing the matrices X  and X  of the first modes X   and X   , respectively, and the two matrices T  and T  , Finally we obtain the reduced brake model: where  defines the generalized degrees of freedom of the reduced brake system.We have Algorithm in universal formal pseudocode is provided in [22].

Numerical Simulations
The numerical simulations and results are decomposed into two main parts: the stability analysis via CEA and the transient nonlinear analysis.
The stability analysis that is the first step in order to predict the squeal propensity of a numerical brake model is only briefly presented in this study.This is explained by the fact that a comprehensive study was conducted previously by the authors on this subject in [22].
In the second part of this section, numerical simulations of the transient nonlinear self-excited vibrations will be presented.More specifically, efficiency of the reduction technique based on complex modes will be illustrated.All the results will be analyzed and compared to a reference calculation based on the finite element model without any interface condensation.

Stability Analysis.
The stability analysis of the numerical brake model can be performed by using the CEA on the reduced finite element model (see (7)).As a reminder, the nonlinear steady-state operating point and the Jacobian matrix of the linearized expression of the nonlinear contact and friction forces F NL at the disc/pad interfaces are taken into account in this calculation.
The stability of a brake system is determined by the sign of real part of eigenvalues: an equilibrium is asymptotically stable if all eigenvalues have negative real parts; it is unstable if at least one eigenvalue has positive real part.In this last case instability is generated from the equilibrium point and selfexcited vibration occurs.
Figures 2(a) and 2(b) display the evolutions of real parts and frequencies as a function of the friction coefficient, respectively.As previously demonstrated in [22], the coalescence pattern is well predicted by the reduction process based on complex interfaces modes (in comparison with the reference numerical model).The curves highlight the phenomenon referred to as mode coupling: as the friction coefficient increases, the real parts increase or decrease slowly and the Hopf bifurcation point occurs for a friction value of 0.48.The associated frequencies of the two modes involved in the coalescence (i.e., the stable and unstable modes) tend to get closer.It is observed that the two modes reach not exactly the Hopf bifurcation point (for  = 0.48) due to the non-equally damped configuration between the two modes (see [26,27] for more details on the effect of damping for the coalescence patterns).After the Hopf bifurcation point, the frequencies of the two modes remain close and their associated real parts become opposite: one mode is stable (i.e., the associated real part is negative and decreases versus the increase of the friction coefficient) whereas the other one is unstable (i.e., the associated real part is positive and increases versus the increase of the friction coefficient).

Transient Simulation and Self-Excited Vibrations.
Based on the previous CEA results (see Figure 2 and Section 3.1), the following friction coefficients were chosen for the transient nonlinear simulation and analysis: (i)  = 0.4 with no instability.
Firstly the time series (transient dynamics and selfexcited vibrations), then the spectral content by using short time Fourier transform (to visualize faster time-varying signals), and finally the limit cycles and the Fast Fourier Transform (to represent the nonlinear stationary signal in the frequency domain) are studied in detail.The displacements selected to illustrate the results correspond to (i) one degree-of-freedom of the disc in the direction normal to the disc surface (i.e., the direction normal to the frictional interface), (ii) one degree-of-freedom of the disc in the direction tangential to the disc surface.
For each case, results from the reduction process based on complex interfaces modes are compared with those from the reference model.Convergence versus the number of interface modes has been conducted (results not presented in the following in order not to overload the paper content).In the following, it has been chosen to show results with 15 interface modes.
Figures 3(a), 3(b), 4(a), and 4(b) depict the time series and the spectrograms for  = 0.4 for one degree-of-freedom of the disc in the direction normal and tangential to the contact surface, respectively.It appears that the initial vibration amplitude decreases until no more movement occurs.The same behavior is obtained for the reduction process based on complex interfaces modes and the reference numerical model.It is also noted that the time simulation agrees with the CEA results that predict no instability for  = 0.4.
As the friction coefficient is increased to  = 0.5, the vibration amplitudes grow (see Figures 3(c), 3(d), 4(c), and 4(d)).It can be seen that a fast increase of the transient oscillations (i.e., the initial growing transient oscillations with the highest amplitudes) is followed by a slight decrease until the nonlinear stationary periodic motions.Evolution of the transient vibrations and the associated spectrograms during the transient oscillations and stationary self-excited vibrations for the reference model and the reduced system are very similar (compare Figures 3(c) and 3(d) and Figures 4(c) and 4(d), resp.).Moreover, Figures 5(a), 5(b), 6(a), and 6(b) give the limit cycles and Fast Fourier Transform on the stationary self-excited vibrations for both the normal and tangential directions.It is clearly shown that the results based on the reduction technique are in perfect agreement with the reference performed via the finite element model without any interface condensation.As indicated in Figures 5(b) and 6(b), the FFTs exhibit two significant peaks: the dominant peak is related to the contribution of the fundamental frequency  = 394 Hz.It is noted that this frequency for the stationary self-excited oscillations is very close (but not equal) to the frequency of the unstable mode predicted by the CEA at  = 0.5, as indicated in Table 2. Also the frequency of the second predominant peak (at 788 Hz) is equal to a multiple (i.e., a superharmonic) of the fundamental frequency  and corresponds to the second harmonic component (i.e., 2).The contribution of 2x harmonic component is also visible on the limit cycles with the appearance of a "small-inside loop."Moreover it is also observed that the FFTs exhibit Increasing the friction coefficient further to  = 0.6 does not show major differences, as illustrated in Figures 3(e), 3(f), 4(e), and 4(f).The squeal event corresponds to nonlinear harmonic self-excited friction-induced vibrations of the mode-coupling mechanism with the participation of one single instability: we have only one fundamental frequency (at  = 367 Hz) and its superharmonics participate in the nonlinear oscillations.However, it appears that the nonlinear behavior is more complex with the major contributions of not only the 1x and 2x harmonic components (at about 367 Hz and 734 Hz) but also a strong participation of the 3x harmonic component (at about 1101 Hz) during the transient time simulation of squeal events (see Figures 5(d) and 6(d)).Considering more specifically the nonlinear steady-state friction-induced motions, the FFTs exhibit other superharmonics .However, all the combinations with  > 3 are generally less important than the fundamental frequencies  as well as the second and third harmonic components (i.e., 2 and 3).Considering more specifically the FFT for  = 0.6 in the tangential direction (see Figure 6(d)), two additional comments can be made: firstly, the major contribution comes from the 3x harmonic component.Secondly, the contributions of the 5x and 2x harmonic components are greater than that of the fundamental frequency.It can also be noted that the fundamental frequency of the stationary selfexcited oscillations is slightly different from the frequency of the unstable mode predicted by the CEA at  = 0.6, as indicated in Table 2.These results highlight the need to consider time analysis not only to estimate the contribution of each harmonic components of the self-excited vibrations but also to predict the value of the fundamental frequency and the main frequency contribution during squeal events.It can be noted that the values of the fundamental frequency calculated from the reference model and the reduced system for the transient time analysis are very close (as shown in Table 2 for  = 0.5 and  = 0.6).Finally, all the results via the reduction technique based on complex modes (for the time series, the spectral content, and the limit cycles with the Fast Fourier Transform) are in perfect agreement with the reference performed via the finite element model without any interface condensation.This illustrates the efficiency and robustness of the proposed reduction method.The computation time and the number of degrees of freedom of the numerical model are extensively reduced by using the reduction process based on complex interfaces modes.This demonstrates the potential of the proposed method (see Table 3 that highlights the less computational cost of the proposed method).

Conclusion
Extension of the double modal synthesis approach [22] has been proposed and used to study the self-excited vibrations and time analysis during squeal events.By conducting a direct comparison with a reference model, the proposed modal reduction strategy has been found to be highly effective not only to predict the transient and stationary amplitudes but also to estimate the nonlinear vibration signature during squeal events.

Table 1 :
Finite element details and material and geometric properties of the simplified brake system.

Figure 2 :
Figure 2: Stability analysis based on the reduced model (red marks) and the reference model (black lines).(a) Real parts.(b) Imaginary parts.

Figure 5 :
Figure 5: Limit cycles and Fast Fourier Transform on the stationary self-excited vibrations for one degree-of-freedom of the disc at the normal component: (a, b)  = 0.5, (c, d)  = 0.6, (a, c) limit cycles, and (b, d) FFT, black color for the reference and red color for the reduced model with 15 interface modes.

Table 2 :
Frequency of the unstable mode via CEA and time analysis.

Table 3 :
Computational costs of the proposed reduced method and the reference model.