## Abstract

The effect of damping on the behaviour of oscillations in the vicinity of bifurcations of nonlinear dynamical systems is investigated. Here, our primary focus is single degree-of-freedom conservative systems to which a small linear viscous energy dissipation has been added. Oscillators with saddle–node, pitchfork and transcritical bifurcations are shown analytically to exhibit several interesting characteristics in the free decay response near a bifurcation. A simple mechanical oscillator with a transcritical bifurcation is used to experimentally verify the analytical results. A transcritical bifurcation was selected because it may be used to represent generic bifurcation behaviour. It is shown that the damping ratio can be used to predict a change in the stability with respect to changing system parameters.

## 1. Introduction

The development of methods for predicting impending bifurcations in dynamical systems has recently gained considerable attention in the literature. One well-known indicator of an approaching bifurcation is *critical slowing down*, where a system responds sluggishly to external perturbations [1]. This phenomenon has been observed in biological systems as an early warning sign of extinction in species [2]. It has also been used to forecast climate shifts [3,4] and to predict phase transitions in statistical mechanics [5,6].

In the study of mechanical oscillators, stability reversals are a common concern. However, because these are second-order systems, the prediction of these bifurcations is made more difficult as inertial forces may mask (or exacerbate) an underlying static instability. One common method of predicting an oncoming instability is to track the natural frequency of the free response [7]. The undamped natural frequency of small (linear) oscillations about an equilibrium will tend to zero as the equilibrium approaches an unstable (buckling) point. This is because the undamped linearized natural frequency is proportional to the square root of the local tangential stiffness of the system. For experimental (necessarily damped) systems, however, the undamped natural frequency is not directly observable, and one instead needs to track the system eigenvalues to predict the onset of an instability. The point at which the real part of the largest eigenvalue of the system becomes positive indicates the onset of divergent/unstable behaviour. However, in viscously damped nonlinear systems, this instability (with respect to a change in parameters) is frequently preceded by a transition from an underdamped to an overdamped response when the imaginary parts of the eigenvalues vanish *prior to* the stability reversal. This is made apparent by considering the definition of the standard damping ratio for a spring–mass–dashpot system (), which is , where *β* is the viscous damping coefficient, *k* is the local spring stiffness, *ω*_{n} is the undamped natural frequency and *m* is the system inertia for linearized oscillations about a static equilibrium. Clearly, a loss of stability (*k*→0) with the change of a system parameter will be preceded by a transition from underdamped (*ζ*<1) to overdamped (*ζ*>1) behaviour. This phenomenon has been shown to occur for saddle–node and supercritical pitchfork bifurcations [7,8,9].

This paper presents a survey of the effects of different generic bifurcations on the damping ratio of nonlinear oscillators, along with experimental results for the transcritical bifurcation. The transcritical is of particular interest because it may be used to generate saddle–node bifurcations as well as equilibrium path intersections similar to pitchfork bifurcations. This is made evident by consideration of the generic transcritical bifurcation within the context of a lightly damped single degree-of-freedom (SDOF) mechanical oscillator,
1.1The equilibrium paths for this system are given by
1.2and hence for a given control parameter *λ* we obtain corresponding equilibria,
1.3These values are stable for 2*x*_{e}>*λ*. Figure 1*a* shows equilibrium paths for a number of initial imperfections, *γ*. For the perfect case (*γ*=0), the trivial path (*x*_{e}=0) has a simple exchange of stability with the non-trivial path (*x*=*λ*) at the origin. For small initial imperfections, we obtain the characteristic dependence of behaviour on the sign of *γ*. For negative values of *γ*, the equilibrium path reaches a local maximum (a saddle–node bifurcation), and the system loses stability completely. For positive values of *γ*, the equilibrium path smoothly veers to larger *x*-values under increasing *λ*. The shaded region of figure 1*a* represents unstable behaviour. Thus, in the context of damped linear vibrations, we expect small-amplitude oscillations to decay exponentially within the unshaded region.

The natural frequency and damping ratio associated with these small oscillations will depend, however, crucially on the value of *λ*. Returning to the governing equation of motion, equation (1.1), we can consider the behaviour of linear vibrations by perturbing the equilibrium: *x*=*x*_{e}+*δ*, where *δ* is a *small* value. Placing this back in the equation of motion and linearizing, we obtain
1.4Thus, the linear natural frequency (squared) is given by
1.5This relation is also shown in figure 1*b* for the same set of values. We note the linear relationship for the geometrically perfect case where the natural frequency, as expected, drops to zero (with stiffness) at buckling (*λ*=0). This has been used as a non-destructive predictor of buckling, although for more realistic (multi-degree-of-freedom and continuous) systems the relation is not necessarily so simple. However, for *γ*≠0, the natural frequency merely drops close to zero in the vicinity of buckling, but then increases beyond buckling *for positive values of γ*. For negative values of

*γ*, the natural frequency does indeed drop to zero, and at values below the critical load for the geometrically perfect case. The shaded region indicates imaginary values of the natural frequency, with perturbations close to the equilibrium position growing exponentially with time. Finally, we observe the role of damping under the variation of the control parameter. Introducing the non-dimensional damping ratio discussed earlier for linear vibration theory

*ζ*=

*β*/2

*mω*

_{n}to equation (1.1), one may get 1.6which is plotted in figure 1

*c*. Suppose the physical system has a linear damping coefficient of

*β*=0.67 Ns m

^{−1}, and is geometrically perfect (

*γ*=0). Then, according to equation (1.6), the damping ratio increases with

*λ*, starting from

*λ*=−1.0, say, until the motion becomes critically damped at 2

*ζ*/

*β*=3, and this occurs when

*λ*reaches , as indicated by the vertical dashed line in figure 1

*c*. The system will then experience an increase in damping with a maximum value when

*λ*=0 (in this case infinitely damped). Upon subsequent increase in

*λ*, the damping ratio reduces, passing back through critical damping at , and then onto light damping as the system enters its highly buckled regime. For

*γ*≠0, the range and extent of critical damping is diminished, and, for this specific case, it can be shown that the damping never becomes overdamped provided

*γ*>0.003086

## 2. Damping near generic bifurcations

The generic bifurcations, i.e. saddle–node, sub- and supercritical pitchfork and transcritical bifurcations [10], may be used to completely categorize the behaviour of SDOF oscillators with a nonlinear restoring force, at least for small oscillations about an equilibrium. Consider the SDOF oscillator given by
2.1where *β* is a constant viscous damping coefficient, and the restoring force *f* depends on the parameters *λ* and *γ*, which will be considered forcing and imperfection parameters, respectively, herein. The static equilibria are described by *f*(*x*_{e},*λ*,*γ*)=0, whereas their stability may be obtained by investigation of the linearized dynamics about these equilibria [7]. However, one may determine much more from the linearized equations of motion about an equilibrium than just stability. By equating the constant damping coefficient *β* to the standard coefficient of velocity of a SDOF oscillator 2*ζω*_{n}, with undamped natural frequency *ω*_{n} and damping ratio *ζ*, one may get the damping ratio for small oscillations about an equilibrium,
2.2

Table 1 shows the equilibrium relationship, undamped natural frequency and damping ratio for saddle–node, sub- and supercritical pitchfork and transcritical bifurcations.

To better investigate the equations in table 1, the equilibrium and damping relationships have been plotted along with the eigenvalues of the linearized equation of motion that are given by . Figure 2 shows the saddle–node static equlibria (*λ* versus *x*_{e} in the shaded plane where a solid line indicates a stable equilibrium, whereas a dashed line indicates an unstable equilibrium) plotted against the damping ratio *ζ*, along with the eigenvalues of small oscillations about the equilibria. What is interesting about this figure is that the damping ratio crosses from below to above *ζ*=1 (square) before the equilibrium becomes unstable (circle) along the path indicated by the arrow. It also shows that the damping ratio tends to infinity as the critical point is approached. This will cause the system to behave more like a first-order quasi-static system near the critical point [11]. When looking at the eigenvalues, the transition from under- to overdamped occurs where the imaginary components vanish leaving two real eigenvalues. Eventually the larger eigenvalue (dots) becomes positive, and the equilibrium becomes unstable. Note that the saddle–node bifurcation is a co-dimension one bifurcation, and hence the imperfection parameter has no effect other than shifting the location of the bifurcation; therefore, the case *γ*=0 is sufficient to completely characterize the behaviour.

The subcritical bifurcation is shown in figures 3 and 4 for the perfect and imperfect cases, respectively. For the perfect case (*γ*=0), the primary path, A, again transitions from under- to overdamped *prior to* becoming unstable, whereas the secondary paths, B, never stabilize. In the presence of an imperfection, the primary path in fact goes through a saddle–node bifurcation. Note that the distribution of the dots in figures 3*b* and 4*b* is different because they were plotted using *λ* and *x* parametrization, respectively, given their different descriptions.

Figures 5 and 6 show the equilibria and damping ratio in the vicinity of a supercritical bifurcation. For the perfect case, the primary path, A, is unchanged from the subcritical bifurcation; however, in this case, the secondary paths, B, are stable. All of these paths also exhibit a transition from under- to overdamped response. The imperfect case again produces another saddle–node bifurcation (paths B and D for two different *β*-values); however, the primary path is perhaps more interesting. Along this path, the damping ratio at first grows to a peak (inverted triangle) as the equilibrium nears the (now inaccessible) critical point, and then begins to reduce thereafter. The peak value may or may not be above *ζ*=1 depending on the value of *β* as seen by comparing paths A and C in figure 6. By setting the expression for *ζ* in table 1 equal to unity, one may produce the depressed cubic
2.3which has only one solution for , which implies that *β*_{crit} is the amount of damping needed to cause the system to become critically damped for a given imperfection *γ*.

Finally, figures 7 and 8 show the transcritical bifurcation for the perfect and imperfect cases. In the perfect case, the primary path (trivial solution), A, is again identical to the primary paths in the perfect sub- and supercritical pitchfork bifurcations. The secondary path, B, also yields the same type of behaviour, although rotated and reversed with respect to the chosen path direction. Figure 7 also includes a plot of the potential energy surface to highlight the geometric underpinning of the equilibria, and their stability as turning points on the potential energy surface. Unlike the sub- and supercritical pitchfork bifurcations, the transcritical bifurcation is not symmetric with respect to the imperfection parameter *γ*. Paths C and D in figure 8 are the equilibrium paths for the case *γ*<0, which again exhibit saddle–node bifurcations with transitions between under- and overdamped behaviour. Paths A and B are the equilibrium paths for *γ*>0, in which case no bifurcations occur. Path B behaves very much like the secondary path for the imperfect subcritical bifurcation, as it never stabilizes. Path A remains stable; however, it exhibits the same increasing and decreasing damping ratio as the imperfect supercritical bifurcation. Once again by equating the damping ratio to unity, and forcing a single solution, one is able to obtain an expression for the critical value of the damping coefficient below which the system remains underdamped, which in this case is *β*_{crit}=(64*γ*)^{1/4}.

The preceding discussion shows that two primary types of behaviour may occur in SDOF systems. If the equilibrium path crosses through a critical point, the response will transition from under- to overdamped and then finally to a quasi-static response (or critical slowing down). If an equilibrium path does not cross a critical point, it will still exhibit a changing damping ratio that will grow larger as it passes nearby to a critical point. This indicates that measurement of damping may prove useful in finding traces of unobservable critical points in experimental systems.

## 3. The transcritical bifurcation in detail

In order to verify these findings, a SDOF oscillator with a restoring force that was able to mimic a transcritical bifurcation was constructed. Figure 9 shows a photograph and schematic of the experimental apparatus, where *θ* is the single state of the system. The system parameters are given by:

— mass of roller assembly

*m*_{s}=1.26 kg, mass of each of the two (inextensible) beams*m*_{b}=0.041 kg and mass of centre hinge assembly*m*_{c}=0.121 kg;— spring stiffnesses of the three linear springs

*k*_{1}=18.6 N m^{−1},*k*_{2}=27.3 N m^{−1}and*k*_{3}=28.4 N m^{−1};— length of each of the two beams

*L*=32.4 cm;— spring lengths (corresponding to

*θ*=0)*L*_{2}=24.8 cm and*L*_{3}=24.2 cm;— unstretched spring lengths

*L*_{2us}=*L*_{3us}=10.0 cm;— lengths of the (inextensible) strings attaching springs

*k*_{2}and*k*_{3}to the structure*L*_{s2}=15.2 cm and*L*_{s3}=17.8 cm; and— spring tilt angle (corresponding to

*θ*=0)*ϕ*=45^{°}.

The two system parameters are the displacement-induced loading parameter *a* (the control parameter) and imperfection parameter *γ*. These were adjusted by moving the magnetic blocks connected to the springs. The imperfection parameter was ‘zeroed’ such that springs *k*_{2} and *k*_{3} were in mutual equilibrium, when *θ*≈0 (it is a practical impossibility to remove all bias in an experimental system). The forcing spring *k*_{1} was attached last such that the structure felt no force from this spring when *a*=0. The initial displacement required to achieve a linear spring stiffness for spring *k*_{1} was negligible when compared with the displacement range of the testing; hence, it was assumed to be linear, starting from zero displacement. The angle *ϕ* introduced the asymmetry that allowed the system to mimic a transcritical bifurcation for small displacements.

### (a) Theoretical

The equation of motion describing the dynamics of the system in figure 9 is best derived by using a Lagrangian approach. The kinetic and potential energy of the system are given by
3.1where the kinetic energy in the springs has been assumed negligible as the springs were much lighter than the structure. Note that the third term in the kinetic energy is from the contribution of the beam that is attached to the roller; this beam rotates and translates causing this term to be more complex. Applying these equations to the Euler–Lagrange equation, and noting that the potential energy is a function of *θ* only and that the kinetic energy may be written as , yields
3.2Finally, simplifying the expression for , and adding a damping term, yields the final equation of motion
3.3where the damping coefficient will be determined by fitting to the experimental data. This form of damping does, however, come with the assumption that the system's energy dissipation is primarily viscous, which will be elaborated on later.

The static equilibria of the system are therefore given by the expression *V* ′=*V* ′_{a}+*V* ′_{R}=0. Applying this to the potential energy in equation (3.1) yields
3.4where, in the interest of brevity, the derivative *V* ′_{R} is not shown explicitly. The stability of these equilibria may be investigated by determining the undamped natural frequency of small oscillations about an equilibrium position. Linearizing equation (3.3) for small displacements, *δ*, about an equilibrium position, *θ*_{e}, yields
3.5Therefore, the undamped natural frequency is given by
3.6The equilibrium and undamped natural frequency for the experimental system in the case *γ*=0 is shown in figure 10. The equilibrium relationship (figure 10*a*) clearly exhibits a transcritical bifurcation in the vicinity of the critical point (*a*=*a*_{crit}=32.8 cm, *θ*=0^{°}). The stability of an equilibrium path is indicated by the square of the undamped natural frequency (*ω*_{0}=15.9 rad s^{−1} is undamped natural frequency corresponding to *a*=0), where a positive value indicates a decaying oscillatory (i.e. stable) response, whereas negative values indicate a divergent unstable response. Therefore, given that is positive definite, equation (3.6) clearly agrees with the theorem of minimum potential energy. For larger displacements, the unstable path in figure 10 restabilizes; however, this domain was outside of the region of experimental focus as only the local transcritical behaviour was of interest. An interesting effect shown in figure 10*b* is that, just as in the simplified normal form shown in figure 1, the square of the undamped natural frequency for the trivial solution (*θ*=0) decays linearly even though this system is more complex.

The damping ratio may again be determined by comparing with the standard damping coefficient 2*ζω*_{n}, which yields
3.7

### (b) Experimental

Experimental data were taken for two different values of the imperfection parameter *γ*. As discussed earlier, the perfect case *γ*=0 was very difficult to calibrate experimentally; therefore, tests were instead performed for a small positive and negative imperfection. The test data were taken by fixing the control parameter *a* then perturbing the structure to small oscillations using a short pulse, and recording the free decay with a high-speed camera (Prosilica GC640) that captured the location of a light-emitting diode that was attached to one of the beams. The data processing was done with LabVIEW software.

To ensure that the free response was behaving linearly, data from oscillations were not used until the magnitude was below 5^{°}. This is a somewhat arbitrary choice; however, it is likely to be smaller than needed to ensure that the linear terms dominate the behaviour. For each value of the control parameter *a*, the structure was perturbed a number of times (at least nine times). A typical experimental time series is shown in figure 11 for *a*=5 cm (note that *a*_{crit}=32.8 cm is the buckling point in the case *γ*=0). This time series was used to determine the experimental equilibrium position, undamped natural frequency and damping ratio. The equilibrium position was taken as the average of the final resting positions of the structure after the perturbations.

The undamped natural frequency and damping ratio were determined by fitting the analytical solution of the standard SDOF spring–mass–dashpot model (for ),
3.8where *θ*_{0} was taken as the distance between the first peak of an experimental free decay and the final resting position in order to account for the linearization about a non-zero resting position. The experimental data *prior to* the first turning point were omitted from the analysis in order to ensure that the perturbing force was no longer present. This also simplified the initial conditions giving . Figure 12 illustrates the system identification using a single free decay from figure 11. In figure 12*a*, the experimental data are compared against equation (3.8) for *ζ*=0.047 and *ω*_{n}=15.64 rad s^{−1}. These values were obtained by minimizing the average error between equation (3.8) and the multiple experimental free decay samples for each load parameter value. This was done by selecting *ζ* and *ω*_{n} values, then finding the average squared difference between the experimental and theoretical curves, and then finally averaging the error for the selected *ζ* and *ω*_{n} values. This process was repeated over a fine grid of *ζ* and *ω*_{n} to produce an error surface. The region around the minimum error for figure 11 is shown in figure 12*b*. There is a clear minimum at *ζ*=0.047 and *ω*_{n}=15.64 rad s^{−1}; hence, these values were selected as the experimental damping ratio and undamped natural frequency. Figure 12 shows that the optimal *ω*_{n} is especially easy to find as the error grows rapidly up to 20 times the minimum error as one moves away from the minimum. In the *ζ* direction, however, the error grows more slowly, which made selecting *ζ* more difficult.

The procedure discussed above was repeated for multiple *a* values for two different imperfections. Figures 13 and 14 show a comparison between the experimental and analytical (equations (3.4), (3.6) and (3.7)) results for *γ*=−0.2 cm and *γ*=0.4 cm, respectively. The experimental results are shown as black data points, whereas the analytical expression is a solid line. The control parameter has again been normalized by *a*_{crit}=32.8 cm, whereas the undamped natural frequencies have been normalized by the unforced natural frequency *ω*_{0}. For the analytical curves, *ω*_{0}=16.3 rad s^{−1} and *ω*_{0}=15.8 rad s^{−1} for *γ*=−0.2 cm and *γ*=0.4 cm, respectively, were obtained by inverting the expression *a*_{e}(*θ*)=0, and substituting this angle into equation 3.6. For the experimental data, this was more difficult because the system could not be tested in the unloaded configuration, as the spring *k*_{1} would bind without any prestress. Therefore, the normalization was obtained by extrapolating the first four data points above *a*=0 (which appear approximately linear by visual inspection) to yield *ω*_{n}=16.7 rad s^{−1} and *ω*_{n}=16.4 rad s^{−1} for *γ*=−0.2 cm and *γ*=0.4 cm, respectively, which is in close agreement with the analytical results. These figures both show an excellent agreement between the experimental and analytical equilibrium relationship along with the natural frequency. As shown in the generic transcritical bifurcation discussed earlier, the squared undamped natural frequency exhibits an approximately linear trend for the primary path. For *γ*=−0.2 cm, the primary path remains stable, and hence the squared natural frequency remains positive, whereas for *γ*=−0.4 cm the system exhibits saddle–node bifurcations where the squared natural frequency becomes negative.

The analytical damping coefficient in was determined empirically by visually fitting *ζ* from equation 3.7 to the data points in figures 13*c* and 14*c*. Interestingly, the curve that seemed to best fit both of these datasets is given by . This function does not properly capture the damping coefficient for *θ*=0, because . Therefore, a function of the form where *β* is some constant, may better catch the overall dissipation of the system. This is, however, an overdetermined problem because many forms of damping could be used to fit the experimental data for a finite set of tests. However, the chosen function accurately matches the experimental data over the domain studied. The damping ratio, as expected, clearly exhibits a peak for load values near *a*_{crit}, although values above unity were not observed experimentally.

It is not immediately clear what effect Coulomb damping had on the selection of the function . However, it is known that the Coulomb damping effect will grow with the equilibrium angle. This is because the primary source of dissipation, the roller mechanism, is increasingly displaced by perturbations as the equilibrium angle grows. Therefore, the behaviour of the damping ratio, particularly the decreasing trend at load levels above *a*_{crit}, is entirely owing to the geometric effect exhibited in this paper.

## 4. Conclusions

It is shown that bifurcation in the equilibria of a system can cause drastic changes in the damping ratio of linearized oscillations about a nonlinear equilibrium. A SDOF experimental link model is used to show that the phenomena obtained theoretically may also be observed in real systems. Two interesting phenomena were observed. First, a stable equilibrium path will exhibit an increasing damping ratio, as a critical point is approached. In the case that the equilibrium path passes through a critical point, the damping ratio will tend to infinity as the critical point is approached, which is the critical slowing down effect. The second, and perhaps more interesting, observation is that, along stable equilibrium paths that do not exhibit any bifurcations (owing to imperfections), the damping ratio will still grow, perhaps even above unity, as the equilibrium path passes near the critical point.

## Acknowledgements

The authors acknowledge the support of Air Force Office of Scientific Research (AFOSR) grant no. FA9550-09-1-0204.

## Footnotes

One contribution of 17 to a Theme Issue ‘A celebration of mechanics: from nano to macro’.

- © 2013 The Author(s) Published by the Royal Society. All rights reserved.