Abstract
A hybrid friction model was recently developed by Azizian and Mureithi (2013) to simulate the friction behavior of tube-support interaction. However, identification and validation of the model parameters remains unresolved. In previous work, the friction model parameters were identified using the reverse harmonic method, where the following quantities were indirectly obtained by measuring the vibration response of a beam: friction force, sliding speed of the force of impact, and local displacement at the contact point. In the present work, the numerical simulation by the finite element method (FEM) of a beam clamped at one end and simply supported with the consideration of friction effect at the other is conducted. This beam is used to validate the inverse harmonic balance method and the parameters of the friction models identified previously. Two static friction models (the Coulomb model and Stribeck model) are tested. The two models produce friction forces of the correct order of magnitude compared to the friction force calculated using the inverse harmonic balance method. However, the models cannot accurately reproduce the beam response; the Stribeck friction model is shown to give the response closest to experiments. The results demonstrate some of the challenges associated with accurate friction model parameter identification using the inverse harmonic balance method. The present work is an intermediate step toward identification of the hybrid friction model parameters and, longer-term, improved analysis of tube-support dynamic behavior under the influence of friction.
1. Introduction
The friction model is an essential element in the detailed analysis of the dynamics of steam generator tubes in the nuclear industry [1]. Most of the friction models currently used to simulate tube-support interaction are cited in [2]. These are special cases of static friction models. The velocity-limited friction model (VLFM) [3] is a continuous Coulomb model without Stribeck effects [4]. The force balance friction model (FBFM) [5] and the spring damper friction model (SDFM) [6] are two models based on springs and dampers. The principle of the Karnopp friction model [7] is used in the FBFM friction model. Recently, a hybrid friction model [1, 8] has been developed to model tube-support interaction. In the hybrid model, all the properties and benefits of other friction models are included, namely, the dynamics of the Dahl model [9] (i.e., hysteresis effect), the dynamics of bristles [10], and the Stribeck effect [4] (transition from the static friction limitation to the kinetic friction limitation). This model also indirectly takes into account the distribution of the stresses in the contact area according to the principle of Cattaneo–Mindlin [11] to model the presliding phenomenon.
In finite element modeling (FEM) codes, the Coulomb friction model and the decay friction model are widely used; these are combined with the principle of the velocity limit [3] to model the presliding (sticking) regime [12]. Most of the numerical algorithms for analyzing dynamic friction in FEM codes are presented by Oden and Martins [13]. Diehl [14] has also numerically investigated the friction effect of a circular rigid body in sliding contact with a flexible beam.
In the present work, we validate the friction coefficient identified experimentally in previous work [15] using the finite element method (FEM) using the Abaqus software.
Four standard friction model parameters have been identified by direct friction force measurement. In the present ongoing work an indirect approach, based on acceleration measurement, has been proposed [15–17]. To identify the friction model parameters experimentally using this indirect method, accurate nonlinear normal modes are needed. Nonlinear normal modes (NNMs) and the principle of inverse harmonic balance (IHB) based on the harmonic balance (HB) method [18] were developed and presented in previous work [15–17].
In the previous work [15–17], the parameters of the Dahl [9] and LuGre [19] friction models, which are, respectively, based on Coulomb and Stribeck friction models were reported. In this paper, therefore, we analyze these two friction models to validate their friction coefficients identified also in the work [15], where the basic friction models identified parameters have been used to identify the parameters of the more complex models.
2. Principle of the Inverse Harmonic Balance (IHB) Method
Figure 1 shows a schematic of the test rig used to extract the contact force (friction and impact) and displacement. The specimen is a beam simply supported at one end (but allowing sliding hence friction effects) and clamped at the other. In Figure 2 all the forces acting at the contact point are represented; for more details see [15–17, 20].
The beam equation of motion iswith the following boundary conditions:
The characteristics of the beam are presented in Table 1 below. The beam response is measured by six accelerometers, where their type and position are cited in the previous works [15–17].
Based on modal superposition, the response of the system may be written in the following form:where are the generalized coordinates, are the normal modes of the system.
Using the modal superposition principle, equation (3), the boundary conditions equations, equation (2), and the Galerkin method, equation (1) can be rewritten as follows:where
Modal orthogonality may be expressed as follows:and is the half cylinder mass (contact element).
In our previous work, different methods were proposed to identify the nonlinear normal modes (NNMs). Here, we present only a summary of the most accurate method, where the NNMs are computed based on the coupled subharmonic approach. The NNMs are calculated in the following manner:
Firstly, the Fourier series of the response signal is written as follows:
Next this equation is written aswhere is the first subharmonic form (or cosines subharmonic form) and is the second subharmonic form (or sine subharmonic form). The modal equation of motion (equation (4)) becomeswhere and are the subharmonic form indices; represents the subharmonic form index used to normalize the modal equation of motion (equation (10)) (as the test function of Galerkin method [21]), and l is the second subharmonic form. The values of k and l can be chosen as follows: (first subharmonic form), (second subharmonic form) or , . In the present work, and were chosen. The modal orthogonality equation (5) becomeswhere
Equations (8) and (9) are used to calculate and by the normalization of from and . Equation (13) is used in the harmonic forms normalization.
Finally, the nonlinear normal modes can be represented by a natural smoothing spline or a series of trigonometric functions (equation (14)):where is chosen equal to 3 to minimize the fitting error.
3. Friction Models
Coulomb and decay friction models are widely used in finite element (FEM) software. In previous work [22], we analyzed these two friction models and validated the friction coefficient identified experimentally in the works [16, 17] using a 1D element. The results showed that the 1D element is not capable of modeling the contact problem with friction. Following the improvement of the inverse harmonic balance method [15], five new parameters of the hybrid friction model were reidentified using the newly calculated experimental results. In the present work, we analyze the different implementations of these models to validate the friction coefficient identified experimentally in [15] using a plane stress element. Firstly, we present the principles and formulation of FEM friction models.
3.1. Coulomb Friction Model
In the Coulomb friction model, the friction force is computed in an analytical formulation:where is the kinetic friction coefficient, is the resultant of the normal forces at the contact point, the velocity at the contact point, and is the resultant of the tangential force at the contact point in the sliding direction.
Generally, in the sticking regime, the resultant of the tangential force is less than the kinetic friction force . The sliding regime begins when the resultant of the tangential force reaches and surpasses this limit (). This formulation is computed using the Lagrange method in the FEM software Abaqus [12, 23]. The formulation of the Coulomb friction model using the penalty method in Abaqus is similar in form to the Karnopp friction model [7] or the velocity-limited friction model [3]. The latter is used to model the sticking (or presliding) regime (called the “elastic slip” regime in Abaqus [12, 24]). The Coulomb friction model with a FEM penalty formulation is given bywhere is the critical friction force and the critical elastic slip.
In general, or is defined by the user in most FEM software. However, in Abaqus, there are two options to define the limiting values: in the first, or are defined by the user (as above), while in the second, automatic option, is calculated during the simulations using the formula:where is the characteristic contact surface length and is the slip tolerance; the default value, [12].
The slip tolerance is also defined as the maximum allowable “elastic slip” (presliding displacement), expressed as a fraction of a characteristic length. The characteristic length is the average length of the contact surface elements (half cylinder). A value of 0.5% is considered typical and used by default in the code. Alternatively, an absolute value (e.g., based on known experimental values) of the characteristic length may be specified. Note that the “elastic slip” is, strictly, not slip at all but rather refers to a relative presliding displacement between the contact surfaces. The relative surface displacement is possible due the contact area asperities which deform elastically without breaking. Gross slip occurs once all the asperities are broken.
3.2. Decay Friction Model (Stribeck Friction Model)
The decay friction model is a special case of the Stribeck friction model [4] (equation (20)), with the exponent is equal to one. Equation (20) represents the general Stribeck friction model [4] in the slip regime. However, in the presliding (or elastic slip) regime, the friction force takes the same value in equations (16)–(18) with the replacement of by in these equations (16) and (18):where is the static friction coefficient, is the Stribeck velocity, and is the Stribeck exponent.
To validate the friction model parameters using a FEM code, we must validate the other parameters used as well as the contact algorithms (surface to surface or nodes to surface) and the slip tolerance . In the present work, the surface-to-surface contact algorithm with the hard contact model is considered. One of the important parameters investigated is the slip tolerance .
3.3. Hybrid Friction Model
This paper presents intermediate results in a longer-term project on hybrid friction model parameter identification and validation using the simple friction models used to create the general hybrid friction model.
In the previous work [15], five parameters of the hybrid friction model were identified using Dahl [9] and LuGre friction models [19]. Some of these parameters are the same as those of the Coulomb and decay (Stribeck) friction model. Starting then with these simpler models, some of the challenges associated with FEM-based parameter identification have been highlighted. Importantly, it is clear that friction models, as implemented in commercial FEM codes, should be used with caution. When quantities intimately related to the details of the friction model are considered, the resulting physical outputs may be significantly affected and far from true values. For instance, the slip tolerance is an adjustable parameter that is found in the numerical friction models. This parameter has, however, been found to have little effect on the computed friction force and displacement when physically realistic values are used.
The hybrid friction model [1, 8] (Figure 3) is based on the principle of the LuGre model [19]. All properties and benefits of other friction models are included, namely, the dynamics of the Dahl model [9] (i.e., hysteresis effect), the dynamics of bristles [10], and the Stribeck effect [4]. This model also takes into account the distribution of the stresses in the contact area according to the principle of Cattaneo–Mindlin [11] to model the presliding phenomenon.
The formulation of the hybrid friction model [1] is represented in the following equations:where , , and are, respectively, the elastic, plastic, and partial slip relative displacements. The presliding or sticking regime is modeled by stiffness , plastic damping and presliding damping , and two transition stiffness coefficients, elastic-plastic and plastic presliding , and the slip regime is modeled by the coefficients of the Stribeck damping which is given bywhere is the Stribeck function (or Stribeck friction model [4]), is the kinetic friction coefficient, is the static friction coefficient, and is the Stribeck velocity.
4. Numerical Simulation
4.1. Finite Element Modeling
In previous work [22], a 1D FEM model (Figure 4(a)) was used to validate the friction coefficient identified experimentally in the works [16, 17]. This model was also compared and validated using the 2D and 3D models. However, in this work, the beam of Figure 1 is modeled by 2D elements FEM model (Figure 4(b)). Then, in conjunction with the friction models presented above, it is used to validate the friction coefficient identified experimentally in the previous work [15]. The beam is discretized using 2D plane strain elements. The validation of this FEM model is done using a multistep mesh analysis. This includes (i) cantilever beam end static deflection calculation, (ii) modal analysis, and (iii) model size optimization. However, only the modal analysis and the model size (calculation time) were used in the choice of the meshing method and the element size of the half-cylinder (contact element). The results of this analysis are presented in Figures 5 and 6. The description of the elements presented in legend of those figures is defined in Table 2.
(a)
(b)
(a)
(b)
(a)
(b)
In general, the triangular elements are the most commonly used to mesh the cylindrical form because of geometrical compatibility. On the other hand, we have solved the problem of the quadrilateral elements’ incompatibility for this geometry type by dividing the circular section of the cylinder into four parts. This means that our contact element (half cylinder) should be divided in two parts (two quarter cylinders) (Figure 4(b)). The other contacting surface is a rigid surface for both models.
4.2. Finite Element Model Validation
These analyses are necessary in this work to optimize our model mesh to reduce the simulation time. For example, using our optimal mesh (Figure 4(b)), the simulation for each iteration (one frequency excitation of 1.2 seconds using direct implicit dynamic integration) takes 90 minutes using a PC having an Intel i7 3 GHz CPU and 24 GB RAM. To obtain one case of parameters analysis results this translates to 30 hours.
Figures 5 and 6 show that the quadrilateral elements give better results than triangular elements, in both objectives cases analysis (deflection the beam free end, and a modal analysis) using the cantilever beam. The CPS4I (4-node bilinear quadrilateral with incompatible modes) element is the best among all the elements for both cases. For this reason, the element type was chosen to model the beam and the half cylinder (the contact element). The triangular element form is the most commonly used to model the cylindrical structure. In the half cylinder element (contact element) meshing analysis, the solutions adopted to obtain the best results are to divide the half cylinder in two parts making it possible to use the CPS4I element for this geometry, and maintain a space ratio less than two to optimize the number of elements with bias mesh controller concentrated in the contact zone.
The characteristics of the beam modeled by FEM are presented in Table 1. The friction model parameters are presented in Table 3. The parameters , , and in equation (23) were identified experimentally using the Dahl friction model [9] and LuGre friction model [19]. These models are detailed in the previous works [15–17]. The objective function previously used for parameter identification is that given in equation (24). Optimization based on an objective function is used to validate the simulation results here as well. The objective function is more general than the more commonly used simpler form employed previously in [17]. This objective function is more useful to compare signals with different behaviors such as the appearance of the beating phenomenon in one response signal or of signals having different offsets (difference between the min/max values).where is the objective function, the experimental friction force, and the simulated friction force.
All the results are obtained at the excitation force amplitude of 6.5 N. A proportional (Rayleigh) damping is used to model the beam structural damping. The mass proportional damping coefficient is given by the following equation:where is the angular natural frequency of mode and the corresponding damping ratio.
The value of the damping ratio in Table 1 is identified experimentally using an impact test and the logarithmic decrement method.
5. Parameter Validation
All the results of the simulation are obtained using the software Abaqus version 6.12-3. Simulations are carried out for 1.2 seconds, with a 10−4 second fixed time step, using the implicit integration method. Data storage was done in 8 × 10−4 second time intervals in order to optimize the total simulation time. Two models were tested, the Coulomb friction model with two different formulations and Lagrange and Penalty formulations, and the decay friction model. We begin by presenting the modal analysis results.
Table 4 presents the first resonance frequency. The simulations using the two Coulomb (Lagrange and penalty) models give approximately equal first resonance frequencies (with a difference of 0.2 Hz); the predicted frequencies are approximately 3 Hz higher than the experimental frequency of 50 Hz. However, the decay friction model is closest to the experimental results when compared to the Coulomb models. For the resonance frequency, the Lagrange formulation is better than penalty formulation, while the decay friction model frequencies are the closest to the experimental values.
In the following figures, from Figures 7–13, the simulations of two cases for each Coulomb friction model formulation are presented. For the Lagrange formulation, simulations using the two friction coefficients identified via Coulomb and LuGre friction models (Table 3) are presented. However, for the penalty formulation, we present the simulation results obtained using the friction coefficient identified via the Coulomb model with two values of the slip tolerance : the default value (0.005) and the optimal value (0.0015) identified during this analysis using the beam response envelope validation criterion.
(a)
(b)
(a)
(b)
(a)
(b)
(c)
(a)
(b)
(a)
(b)
(a)
(b)
(a)
(b)
Figure 7 presents the first nonlinear normal mode (NNM). Simulation results are compared with the experimental normal modes identified using the inverse harmonic balance method (IHB) and the linear normal mode (LNM) of the clamped-simply supported beam. The simulation results using both friction models have approximately the same NNMs; the NNM form is also very close to the first NNM identified experimentally using the IHB. However, the NNM of the simulation results using the penalty formulation are even closer to the experiments NNM compared to the Lagrange formulation results.
Note that while the simulated 1st nonlinear modes (NNMs) appear similar to the experimental NNM, the small difference is important. This is because the modal spatial derivatives are particularly critical for accurate determination of the friction dynamics at the contact point; this is evident in equation (10). Hence, even though the NNMs look similar, the slight differences lead to significantly different responses (Figure 8). This is largely due to the importance of the first and the second derivatives of the mode form at the contact point which directly affects the resulting computed friction force in equations (4) and (10).
To better compare the results, we can use data shown in Figure 8. In this figure, the beam response envelopes (maximal deflections) using both formulations of the Coulomb friction model are presented. From Figure 8(a), Coulomb friction model with the Lagrange formulation, we find that the friction coefficient identified using the LuGre model yields the closer results to experiments than the value identified using the Coulomb friction model. Contrarywise for the penalty formulation (Figure 8(b)), the value identified using the Coulomb friction model () gave the best results but with the optimal value of the slip tolerance. The default value of the slip tolerance is unable to produce the correct results.
To provide more details on the beam response, the frequency response (FRF) of the system at the driving point (position 555 mm from the clamped end) is presented in Figure 9 for the Coulomb and decay friction models. Similarly to the observations made from Figure 8, the case with for the Coulomb friction model with Lagrange formulation (Figure 9(a)) and the case with for the penalty formulation (Figure 9(b)) are the closest results to the experimental response. These two cases can be considered as the optimal cases (having optimal parameters) for each model formulation.
From Figure 9, we can also observe that the penalty formulation gives smooth system response (without multiple-resonance peaks phenomenon), in the interval 45–56 Hz, compared to the Lagrange formulation, where this phenomenon is not observed experimentally in the interval (from 45 Hz up to the second natural frequency (160 Hz)). However, for both formulations, there are other resonances above 56 Hz, while not presented in the results, we can see the beginning of a second possible resonance after the 1st natural frequency in Figure 9(b). Furthermore, the experimental response occurs over a wider frequency bandwidth. The large response predicted by the Coulomb friction model with both formulations shows that the models do not correctly capture the energy dissipation due to friction. This is partly expected due to the known inadequacies of the Coulomb friction model.
Note, however, that there is gradual improvement in the predictive behavior between the more basic Coulomb and the improved decay friction model. The preliminary result of decay friction model peak frequency, in general, is closer to the experimental resonance frequency value than Coulomb model (Table 3); however, the beam response amplitude (Figure 9(c)) is higher than the experimental result. The results of Figure 9(c) are obtained using the optimal parameters of the penalty formulation and the optimal value of presented in Table 3 from the previous work [15] with the Abaqus decay friction model. In the next step, a complete analysis of the decay and Stribeck models must be done to improve our decay friction model.
Figures 10 and 11 present the slip displacement and slip velocity at the contact point. For both formulations, the simulated resonance slip displacement (Figure 10) of the optimal parameters defined above: for Lagrange formulation and for penalty formulation, leads to the closest values to experimental results. However, all the simulated resonance slip velocities at the contact point (Figure 11), for both formulations, are higher than the experimental values. This may come from the integration and the derivative methods used in Abaqus and/or due the inability of the Coulomb friction model to represent the real behavior of friction especially in the sticking condition. Thus, for both formulations at low frequencies, below the resonance, the system is in the sticking condition (zero slip displacement for the Lagrange formulation and slip displacement equal to the critical displacement for the penalty formulation). For the higher frequencies above the 1st resonance frequency, the slip displacement (Figure 10) values are closer to the experimental values. The slip velocity (as shown in Figure 11) for the Lagrange formulation is higher than the penalty formulation case. Both formulations thus give higher results with the penalty formulation being 1.5 times the experimental values.
We can also make the same observations in Figure 8, where for the beam frequency response at the driving point, for the low frequencies, the difference between the simulations and the experimental results is quite large. However, for the higher frequencies (above the resonance frequency), there is a good match between the simulations and the experimental results. We conclude that the Coulomb friction model is valid and applicable only for the cases with slipping conditions. Despite the penalty formulation improvement, the Coulomb friction model accuracy remains insufficient to represent the real behavior of the sticking or presliding condition.
Figures 12 and 13 present the contact force, tangential (friction) force, and normal (impact) force results. The friction force results, Figure 12, support the conclusion above for both formulations. The simulated friction force levels at the resonance frequency have closest values to the experimental results when the optimal parameters defined above are used. However, the simulated normal force levels obtained using the optimal parameters (Figure 13) at the resonance frequency are slightly higher than the experimental values, with an error of 4.39 N (4.93%) for the Lagrange formulation and 3.56 N (4%) for the penalty formulation. Thus, even when the results of the friction force, the system response, and the slip displacement are close to the experimental values, both formulations of the Coulomb model lead to slightly different values for the normal force. The continuous nature of the penalty formulation slightly improves the normal force result compared to the discontinuous Lagrange formulation.
Finally, we can use these results (Figures 8, 9, and 11) to obtain the optimal Coulomb friction model parameters for the FEM model. The optimal parameters are ( between 0.0015 and 0.005 with ) for the penalty formulation. For the Lagrange formulation, is less than 0.49 for the response system validation criterion and between 0.52 and 0.49 for the others validation criteria. In conclusion, we note that the accurate optimal parameters of the different Coulomb friction model formulations affect the system response, the friction force, and the slip displacement in the same manner and differently to the normal force.
6. Parameter Sensitivity Analysis
In this section, we investigate the parameter sensitivity of the two Coulomb friction model formulations. In Figures 14–19, we present the simulation results using the parameters of the two cases presented in Figure 7 to Figure 13 plus other cases near the optimal parameters of both formulations. However, we present only the figures (or the validation criteria) where the parameter variation effect is clear. Furthermore, the sensitivity of the two parameters in the penalty formulation is presented separately: subpart “b” in Figures 14–19 for the slip tolerance analysis and subpart “c” in Figures 14–19 for the friction coefficient analysis.
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
In the Lagrange formulation (Figure 14(a)), when the optimal parameter () is increased by 2%, the beam response envelope decreases by 17% relative to the optimal parameter response. However, the beam response envelope increases by 28% of the optimal parameters response when is decreased by 2% below the optimal value .
For the penalty formulation (Figures 14(b) and 14(c)), firstly, the slip tolerance sensitivity was tested by fixing the friction coefficient at 0.52 (Figure 14(b)). It was found that the beam response envelope increased with increasing slip tolerance (which means with increasing critical displacement ). Additionally, it was noted that the optimal slip tolerance could be improved within the interval (0.0011 to 0.0015) by using the experimentally identified . The friction coefficient sensitivity was also tested (Figure 14(c)) by fixing the slip tolerance to two values: the default value and the optimal value . These two cases are analyzed to understand the behavior of the variation. It was found that the beam response envelope increased with decreasing for both cases.
The same approach was used to analyze the other criteria presented in Figures 13–17.
For the following two criteria, FRF at driving point (Figure 15) and slip displacement (Figure 16), the observations of the parameters sensitivity are similar to the beam response envelope (Figure 14) where, the simulation result values increase with decreasing for both criteria using both formulations (Figures 15(a), 15(c), 16(a), and 16(c)), and the simulated values increase with increasing for the penalty formulation. Furthermore, for the penalty formulation, using the slip tolerance optimal value , the simulation result approaches the experimental value, unlike the case using the default value , with the appearance of a multi resonance, confirming the stability response of the slip tolerance optimal value.
For the slip velocity criterion (Figure 17), the observations of the parameters sensitivity are similar as the FRF at driving point (Figure 15) and slip displacement (Figure 16) for the penalty formulation (Figures 17(b) and 17(c)). However, for the Lagrange formulation, the simulated slip tolerance (Figure 17(a)) behaves differently around the optimal friction coefficient, . The simulated value increases away from the experiment at value for increasing or decreasing friction coefficient . This is another confirmation that is the optimal value for the Lagrange formulation.
For the following two criteria, the friction force (Figure 18) and the normal force (Figure 19), the simulation values decrease for decreasing in the Lagrange formulation (Figures 18(a) and 19(a)) and also by decreasing using the penalty formulation (Figures 18(b) and 19(b)). Furthermore, for the penalty formulation, using the default slip tolerance leads to values far from the experiments, and there is a small difference for the values around the optimal value . For the penalty formulation, the friction force (Figure 18(c)) increases with increasing friction coefficient . However, the normal force (Figure 19(c)) decreases with increasing using the optimal value , which leads to the closest values in experiments for both criteria.
Finally, we conclude from this analysis that using the optimal parameter values ( for the Lagrange formulation and with for the penalty formulation) yields the smoothest (the most stable) simulated results and the closest to experimental values for the six validation criteria presented above. This validates and confirms the accuracy of the friction model parameters identified using the inverse harmonic balance method IHB with LuGre model for the Lagrange formulation and with Coulomb friction model for the penalty formulation. We can conclude also that the accurate optimal parameters of the different Coulomb friction model formulations affect the system response, the slip displacement, and the slip velocity in the same manner. The friction force and normal forces are affected differently.
We can use the results of Figures 7–18 to obtain the best (optimal) Coulomb friction model parameters for the FEM model. The optimal parameters are for between 0.0015 and 0.0020 with for the penalty formulation and between 0.48 and 0.49 for the Lagrange formulation.
The friction coefficient identified experimentally using the Coulomb model is the optimal value for the penalty formulation of the Coulomb model. However, the friction coefficient identified experimentally, using the LuGre friction model, is the optimal value for the Lagrange formulation of the Coulomb model. This can be generalized only after testing different cases in future work.
Modeling the Stribeck effect also significantly improves the resonance frequency prediction. The next step involves moving to more sophisticated models—the Dahl, LuGre, and, eventually, the Hybrid friction models. Prior to this, the role played by the numerical implementation of these models in the FEM codes needs to be verified. This is important because the implementation of the friction models in FEM codes is not mathematically trivial and involves adjustments and additional parameters.
Finally, a word concerning the slip displacement: this displacement is an important quantity used to calculate the wear work rate (26). The slip displacement value at the contact point is plotted versus frequency in Figures 10 and 16 for the Coulomb model. The values of the slip displacement generated by the Coulomb friction model are close to the experimental ones at resonance and too small (even equal to zero for the Lagrange formulation) for frequencies lower than the resonance frequency. This leads to smaller computed work rates compared to the experimental values in this frequency interval.where is the work rate, is the normal force, is the sliding velocity, and is the duration of contact.
7. Conclusion
This paper presents the first step to validate, using 2D plane strain FEM elements, the parameters identified in previous work [15] using the inverse harmonic balance method results and the basic models of the hybrid friction model.
The accuracy of the friction model parameters identified using the inverse harmonic balance method (IHB) with the LuGre model is validated using the Lagrange formulation. The parameters identified with Coulomb friction model are validated using the penalty formulation. The accuracy of the optimal parameters of the different Coulomb friction model formulations affects the system response, the slip displacement, and the slip velocity in the same manner but differently from the friction and normal forces.
Both models produce friction forces, normal forces, and slip displacements of the correct order of magnitude compared to the friction force calculated using the inverse harmonic balance method at the resonance frequency. However, their FRF bandwidths are significantly far from the experimental results. Furthermore, Coulomb and decay friction models yield modal parameters (resonance frequency and NNMs) that are close to those of the experiments. The decay friction model yields the results closest to experiments for the NNMs and resonance frequency.
The Coulomb friction model is valid and applicable only for cases with slipping conditions. Despite the improvement gained using the penalty formulation in the Coulomb friction model, it remains insufficiently accurate to represent the real physical behavior for sticking or presliding condition.
The present work demonstrates that both static friction models produce results closest to the experiments at resonance in the slipping regime, but are incapable of accurately representing all the behaviors of friction, especially in the sticking regime, ultimately leading to incorrect work rate estimates in the complete frequency range. The results presented here represent a part of ongoing work aimed at modeling the detailed physics underlying friction.
Nomenclature
| : | Area of the cross section of the beam |
| : | Young’s modulus |
| : | Excitation force by a shaker |
| : | Kinetic friction force limitation |
| : | Excitation force amplitude |
| : | Slip tolerance (fraction of characteristic contact surface length) |
| : | Static friction force limitation |
| : | Quadratic bending moment |
| : | Length of the beam |
| : | Resultant of the normal force |
| : | Optimization objective function |
| : | Static load to ensure permanent contact between the beam and support |
| : | Support reaction (“impact” force) |
| : | Resultant of the force at the contact point in the sliding direction |
| : | Frictional force |
| : | Critical friction force |
| : | Experimental friction force |
| : | Simulated friction force |
| : | Position of the excitation point |
| : | Work rate |
| : | Plastic damping |
| : | Presliding damping |
| : | Stribeck damping |
| : | Thickness of the beam |
| : | Stribeck function |
| : | Elastic stiffness |
| : | Elastic-plastic stiffness |
| : | Plastic presliding stiffness |
| : | Characteristic contact surface length |
| : | Half cylinder mass (contact element) |
| : | Mode index |
| : | Linear generalized coordinates |
| : | Nonlinear generalized coordinates |
| : | Radius of the contact element (half cylinder) |
| : | Time dependence |
| : | Displacement at the contact zone |
| : | Elastic slip |
| : | Position of the accelerometer j |
| : | Position dependence |
| : | System response (beam deflection) |
| : | System response of the accelerometer j |
| : | Elastic slip relative displacements |
| : | Plastic slip relative displacements |
| : | Partial slip relative displacements |
| : | Proportional mass damping |
| : | Proportional mass damping of mode |
| : | Stribeck exponent |
| : | Damping ratio |
| : | Damping ratio of mode |
| : | Kinetic friction coefficient |
| : | Static friction coefficient |
| : | Density of the beam material |
| : | Velocity at the contact area |
| : | Linear normal modes |
| : | Nonlinear normal modes |
| : | First subharmonic form or cosines subharmonic form |
| : | Second subharmonic form or sine subharmonic form |
| : | Natural frequency of mode . |
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This manuscript is based on Dr. Hadji’s thesis [20].