Free and forced vibrations with damping

This example verifies the frequency-dependent spring and dashpot elements available in Abaqus.

This page discusses:

Problem description
Results and discussion
Input files
References
Tables
Figures

ProductsAbaqus/Standard

There are several different mechanisms that can cause damping in a system. In linear viscous damping the damping force is directly proportional to the velocity. In many cases such simple expressions for the damping forces are not available directly. However, it is possible to obtain an equivalent viscous damping coefficient by equating the loss of kinetic and strain energy to the energy dissipation. Hysteretic and viscoelastic damping are two important damping mechanisms that are more complex than linear viscous damping. In the frequency domain these mechanisms can be simulated by using dashpots with viscous damping coefficients that depend on the forcing frequency. Frequency-dependent springs will also be needed for modeling viscoelastic damping.

To illustrate how to model viscous, hysteretic, and viscoelastic damping mechanisms, springs and dashpots with constant and frequency-dependent properties will be used in frequency domain dynamic analyses of one- and two-degree-of-freedom discrete mass-spring-dashpot systems. In addition, viscous damping is modeled in the time domain by using a constant dashpot coefficient.

Abaqus also allows for spring and dashpot properties that depend on temperature and user-defined field variables. This dependence provides an easy means to vary material properties of springs and dashpots during time-domain analysis. In doing perturbation analysis (such as frequency-domain steady-state dynamic analysis) with Abaqus, temperature and field variable variations are not permitted within an analysis step. However, since the base state temperature and field variable values for each perturbation analysis step can be changed, it is possible to perform a multiple-step perturbation analysis that uses different temperature- and field-variable-dependent material properties that correspond to the base state temperature and field variable values. This dependence feature will be illustrated in analyses 2 and 3 described below. These two analyses employ both the direct-solution and the subspace-based steady-state dynamic procedure in Abaqus.

The one- and two-degree-of-freedom mass-spring-dashpot systems are shown in Figure 1. The following dynamic analyses are performed: (1) free vibration of the one-degree-of-freedom system after it is given an initial displacement and then released; (2) steady-state response to applied harmonic loading of the one-degree-of-freedom model with viscous damping; (3) steady-state response to applied harmonic loading of the one-degree-of-freedom model with hysteretic damping; and (4) steady-state response to applied harmonic loading of the two-degree-of-freedom model with viscoelastic damping. In all cases the forcing function is applied to the point mass closest to the anchor point, and numerical results are compared to the exact solutions for the system.

Problem description

The basic constant parameters of the analysis models are as follows:


Spring constant, k	5253.8 N/m (30 lb/in)
Damping coefficient, c	21.02 N/m-s (0.12 lb/in-s)
Mass, m	4.536 kg (0.02588 lb-s²/in)

SPRING1 and DASHPOT1 elements are used in analyses 1–3. SPRING2 and DASHPOT2 elements are used in analysis 4.

In analysis 1 the model is the one-degree-of-freedom system shown in Figure 1. The initial displacement is 25.4 mm (1 in), so the force in the spring is initially 133.4 N (30 lb). The problem is run in two steps: a static step, wherein the initial displacement is imposed, and a dynamic step, during which the structure is allowed to oscillate. The dynamic step is run with automatic time stepping, using two different values for the half-increment tolerance: 44.48 N (10 lb) and 4.448 N (1 lb). The higher value of the tolerance should give moderately accurate results, while the lower value should result in a more accurate solution. Implicit Dynamic Analysis Using Direct Integration gives guidelines for choosing a value for the half-increment tolerance for realistic, multiple-degree-of-freedom systems.

In analysis 2 a harmonic loading of the form $P = \bar{P} \cos (Ω t)$ is applied to the single-degree-of-freedom system, where $Ω$ is the circular frequency. The equation of motion for this system is

m \ddot{u} + c \dot{u} + k u = P .

The direct-solution and the subspace-based steady-state dynamic procedures are used to calculate the steady-state vibrations in this system with low and high viscous damping coefficients, $c = $ 0.12 and 0.24. The dashpot coefficient in this model is defined as a function of the first field variable, and the change of the field variable value is carried out in a dummy general static step placed between two direct-solution steady-state dynamic steps.

Analysis 3 is identical to analysis 2 in all aspects except that hysteretic damping is modeled instead of linear viscous damping. Hysteretic damping, also known as structural or solid damping, is observed in the vibration of many solid materials and can be attributed to internal friction. This form of damping produces a hysteresis loop in the force-displacement plot for each loading cycle that is proportional to the amplitude and tends to stay constant with rising forcing frequency. The energy loss is proportional to the displacement amplitude squared for both viscous damping and for hysteretic damping. This fact suggests that structurally damped systems subjected to harmonic excitation can be modeled as viscously damped systems with an equivalent coefficient of viscous damping that is inversely proportional to the frequency: see Denhartog (1985). The equation of motion for this one-degree-of-freedom system is, thus, written readily as

m \ddot{u} + \frac{γ k}{Ω} \dot{u} + k u = P,

where $γ$ is a damping coefficient and $Ω$ is the forcing frequency. The equivalent viscous damping coefficient is $c_{e q} = \frac{γ k}{Ω} .$ For harmonic motion we have the relationship $\dot{u} = i Ω u$ ; therefore, $\frac{γ k}{Ω} \dot{u} = i γ k u$ , where $i = \sqrt{- 1}$ is the imaginary number. Hence, the equation of motion can also be written as

m \ddot{u} + k (1 + i γ) u = P .

Abaqus also allows direct specification of structural damping; however, this direct specification can be used only in modal-based analysis and is accurate only for small damping values. See Material Damping for further discussion. In this analysis the effects of low damping ( $γ$ =0.125) and high damping ( $γ$ =0.25) are compared, following the same procedure as used in analysis 2. The data set containing the frequency-dependent dashpot coefficients at intervals of 0.05 Hz over the frequency range of 0 to 10 Hz is included (file vibration_dampdata1.inp).

Analysis 4 involves a two-degree-of-freedom system with viscoelastic damping. Viscoelastic materials are often used in a structure to improve the damping characteristics of the structure or its components. In a one-dimensional test specimen made of linear viscoelastic material, an applied cyclic stress $σ^{*} = σ_{0} \exp (i Ω t)$ will result in a steady-state cyclic strain response, $ϵ^{*} = ϵ_{0} \exp (i (Ω t - δ))$ with the same frequency but out of phase by the phase angle $δ$ . The phase angle is also known as the loss angle and is a function of frequency. The damping ability of the material is dependent on it and not on the stress and strain amplitude. The ratio of the stress and strain defines the complex modulus, $E^{*} (Ω) = σ^{*} / ϵ^{*}$ , where the real part is termed the storage modulus and the imaginary part the loss modulus. The equation of motion for the steady-state forced vibration of a single-degree-of-freedom viscoelastic system of mass m is simply

m \ddot{u} + K^{*} (Ω) u = P,

where $K^{*} (Ω)$ is the complex stiffness proportional to the complex modulus $E^{*} (Ω)$ . Making use of the substitution $\dot{u} = i Ω u$ for harmonic motion, we can rewrite the equation of motion as

m \ddot{u} + c (Ω) \dot{u} + k (Ω) u = P,

where $k (Ω) = ℜ (K^{*} (Ω))$ and $c (Ω) = \frac{1}{Ω} ℑ (K^{*} (Ω))$ . Referring to Frequency Domain Viscoelasticity, we can identify that the equivalent viscous damping coefficient is $c (Ω) = ℜ (κ^{*}) k_{\infty}$ and the spring stiffness is $k (Ω) = (1 - Ω ℑ (κ^{*})) k_{\infty}$ , where $κ^{*} = κ^{*} (Ω)$ is the Fourier transform of the nondimensional relaxation function $κ (t)$ and $k_{\infty}$ is the long-term spring stiffness. The equation of motion for viscoelastic damping resembles the one for hysteretic damping to the extent that viscoelastic damping can also be simulated in discrete mass-spring-dashpot systems using frequency-dependent springs and dashpots. This form of damping is simulated in the two-degree-of-freedom discrete mass-spring-dashpot system shown in Figure 1 with the following parameters: $m_{1} = m_{2} = m$ and $K_{1}^{*} = K_{2}^{*}$ , such that the real and imaginary moduli are $k_{1} (Ω) = k_{2} (Ω) = k (Ω) = (1 - Ω ℑ (κ^{*})) k_{\infty}$ , and $c_{1} (Ω) = c_{2} (Ω) = c (Ω) = ℜ (κ^{*}) k_{\infty} .$ The frequency dependence of $κ^{*} (Ω)$ assumes the power law formula $κ^{*} (Ω) = κ_{1}^{*} f^{- b}$ , where b is a real constant, $κ_{1}^{*}$ is a complex constant, and $f = Ω / 2 π$ is the frequency in cycles/time. The equation of motion for the two-degree-of-freedom system is now readily developed and can be written as

[\begin{matrix} m & 0 \\ 0 & m \end{matrix}] {\begin{matrix} {\ddot{u}}_{1} \\ {\ddot{u}}_{2} \end{matrix}} + [\begin{matrix} K_{1}^{*} + K_{2}^{*} & - K_{2}^{*} \\ - K_{2}^{*} & K_{2}^{*} \end{matrix}] {\begin{matrix} u_{1} \\ u_{2} \end{matrix}} = {\begin{matrix} P_{1} \\ P_{2} \end{matrix}} .

Since harmonic loading of the form $P_{j} = {\overset{ˉ}{P}}_{j} \exp (i Ω t)$ produces the harmonic oscillation $u_{j} = \bar{u_{j}^{*}} \exp (i Ω t)$ with the complex oscillation amplitude $\overset{ˉ}{u_{j}^{*}} = {\bar{u}}_{j R} + i {\bar{u}}_{j I}$ , the equation of motion for the loading parameters used in this analysis, $P_{1} = \bar{P} \cos (i Ω t)$ and $P_{2} = $ 0, can be rewritten in terms of the real and imaginary parts of the oscillation amplitudes as follows:

[\begin{matrix} 2 k - m Ω^{2} & - k & 2 Ω c & - Ω c \\ - k & k - m Ω^{2} & - Ω c & Ω c \\ 2 Ω c & - Ω c & - 2 k + m Ω^{2} & k \\ - Ω c & Ω c & k & - k + m Ω^{2} \end{matrix}] {\begin{matrix} {\bar{u}}_{1 R} \\ {\bar{u}}_{2 R} \\ - {\bar{u}}_{1 I} \\ - {\bar{u}}_{2 I} \end{matrix}} = {\begin{matrix} \bar{P} \\ 0 \\ 0 \\ 0 \end{matrix}} .

The frequency-dependent spring and dashpot properties are generated by a Fortran program using the basic model constants for the mass, m, and for the spring, $k_{\infty} .$ In addition, the parameters b=1.38366, $ℜ (κ_{1}^{*})$ =2.3508 × 10⁻², and $ℑ (κ_{1}^{*})$ =6.5001 × 10⁻² are used. This form of the power law dependence of frequency of $κ^{*} (Ω)$ does not describe the viscoelastic properties for all frequencies accurately. In particular, this formula is incorrect for low frequencies since the stiffness becomes negative. Therefore, values computed using this formula for frequencies below 0.77 Hz are discarded in this analysis. The frequency-dependent data for the dashpot coefficients and for the spring stiffness are written at intervals of 0.035 Hz over the frequency range of 0.77 to 14 Hz.

Results and discussion

For analysis 1 the exact solution for damped free oscillation is

u = u^{0} \exp (- η ω t) [\cos ω \sqrt{1 - η^{2}} t + \frac{η}{\sqrt{1 - η^{2}}} \sin ω \sqrt{1 - η^{2}} t],

where $ω = \sqrt{k / m}$ is the natural frequency of the undamped system (34.05 rad/s in this analysis), $η = c / 2 \sqrt{m k}$ is the ratio of damping to critical damping (0.068 in this analysis), and $u^{0}$ is the initial displacement of 25.4 mm (1 in).

The exact solution and the Abaqus solutions obtained using the different values of the half-increment tolerance are plotted in Figure 2. The tighter tolerance provides the more accurate solution, showing a slight phase shift later in the response. The looser tolerance shows considerably more phase shift, as expected.

At any time during the analysis Abaqus can provide a summary of the energy present in the structure, as well as quantities such as viscous and plastic dissipation. Summation of the various energy quantities yields an energy balance. Comparison of this balance with the initial strain energy of the system yields the energy lost due to numerical damping in the time integration operator. Table 1 is a summary of all the energy terms at the end of the problem ( $t = $ 0.7 seconds). Since the initial strain energy is 1.695 N-m (15 lb-in), the numerical damping loss is 1% for the smaller value of the tolerance and 9.1% for the larger value.

For analysis 2 the steady-state response of the viscously damped single-degree-of-freedom system subjected to a cosine forcing function is given by

u = \bar{u} \cos (Ω t - ϕ),

where the amplitude of oscillation is

\bar{u} = \frac{\bar{P}}{\sqrt{{(k - m Ω^{2})}^{2} + {(c Ω)}^{2}}} = \frac{\bar{P}}{k} \frac{1}{\sqrt{{(1 - {(Ω / ω)}^{2})}^{2} + {(c Ω / k)}^{2}}},

with $ω = \sqrt{k / m}$ , and the phase angle of the response is

ϕ = \arctan \frac{c Ω / k}{1 - {(Ω / ω)}^{2}} .

The response calculated by the direct-solution and subspace-based steady-state dynamic procedures are in exact agreement with these solutions. The amplitude and the angle of phase lag for the frequency range of 0 to 10 cycles/time are shown in Figure 3 and Figure 4, respectively.

For analysis 3 the steady-state response of the structurally damped single-degree-of-freedom system subjected to a cosine forcing function can be obtained from the solutions for the viscously damped case by replacing the constant dashpot coefficient with the equivalent frequency-dependent dashpot coefficient,

c = k γ / Ω,

so the amplitude of oscillation is

\bar{u} = \frac{\bar{P}}{k} \frac{1}{\sqrt{{(1 - {(Ω / ω)}^{2})}^{2} + γ^{2}}}

and the phase angle of the response is

ϕ = \arctan \frac{γ}{1 - {(Ω / ω)}^{2}} .

The Abaqus solutions obtained by the direct-solution and subspace-based steady-state dynamic procedures are again in exact agreement with these analytical results, as shown in Figure 5 and Figure 6. Comparing these results with those in Figure 3 and Figure 4 for viscous damping, two differences are apparent. First, resonance (maximum amplitude) occurs at $Ω = ω$ and not at $Ω < ω$ as in analysis 2. Second, the phase angle for $Ω = $ 0 is $ϕ = \arctan γ$ instead of zero as in analysis 2; therefore, motion with structural damping, where the energy dissipation is rate independent, will never be in phase with the forcing function.

For analysis 4 the steady-state response of the two-degree-of-freedom system with viscoelastic damping is obtained numerically by solving the system of four equations for the real and imaginary parts of the response at the two nodes, from which the response amplitudes

| {\bar{u}}_{j} | = \sqrt{u_{j R}^{2} (Ω) + u_{j I}^{2} (Ω)}

and the phase angles

ϕ_{j} = \arctan \frac{- u_{j I}}{u_{j R}}

are obtained. The solutions obtained by the direct-solution steady-state dynamic procedure in Abaqus are also in exact agreement with the numerical solutions. The results are presented in Figure 7 and Figure 8, which show the amplitudes and angles of phase lag of the response, respectively, at the free nodes for the frequency range of 0.77 to 14 cycles/time. It is important to realize that small intervals must be used in the frequency sweep to obtain results with high accuracy (in particular for the peak response). Furthermore, when the frequency dependence is nonlinear, such as exhibited in this system, the quality of the solution also depends upon the accuracy of the frequency-dependent spring and dashpot data used in the calculations. Abaqus assumes that the properties vary linearly over each frequency interval; consequently, a small interval size should be used in the discretization of the data to minimize interpolation errors.

Input files

vibration_1dof_dyn_haft1.inp: One-degree-of-freedom time-integration dynamic analysis with HAFTOL set to 4.448 N (1 lb).
vibration_1dof_ssdyn_viscous.inp: One-degree-of-freedom direct-solution and subspace-based steady-state dynamic analysis with viscous damping.
vibration_1dof_ssdyn_hyster.inp: One-degree-of-freedom direct-solution and subspace-based steady-state dynamic analysis with hysteretic damping.
vibration_2dof_ssdyn_visco.inp: Two-degree-of-freedom direct-solution steady-state dynamic analysis with viscoelastic damping.
vibration_1dof_dyn_haft2.inp: Problem with HAFTOL set to 44.48 N (10 lb).
vibration_dampdata1.inp: Frequency-dependent damping coefficients used in analysis 3.
vibration_dampdata2.inp: Frequency-dependent damping coefficients used in analysis 4.
vibration_springdata.inp: Frequency-dependent spring stiffness used in analysis 4.

References

Denhartog, J. P., Mechanical Vibrations, Dover, 1985.

Tables

Table 1. Energy balance at 0.7 seconds.
		Solution with half-increment tolerance = 4.448 N (1 lb)	Solution with half-increment tolerance = 44.48 N (10 lb)
Kinetic energy	N-m	0.0472	0.0033
Kinetic energy	lb-in	0.418	0.029
Strain energy	N-m	0.0490	0.1943
Strain energy	lb-in	0.434	1.720
Dissipated energy	N-m	1.5817	1.3445
Dissipated energy	lb-in	14.000	11.900
Total energy	N-m	1.6780	1.5421
Total energy	lb-in	14.852	13.649
Energy loss through numerical damping	N-m	0.0167	0.1526
Energy loss through numerical damping	lb-in	0.148	1.351

Figures