Numerical example: results of nonviscous eigenvalues (rad/s) and eigenvectors.
Nonviscously damped vibrating systems are characterized by dissipative mechanisms depending on the time history of the response velocity, introduced in the physical models using convolution integrals involving hereditary kernel functions. One of the most used damping viscoelastic models is Biot’s model, whose hereditary functions are assumed to be exponential kernels. The free-motion equations of these types of nonviscous systems lead to a nonlinear eigenvalue problem enclosing certain number of the so-called nonviscous modes with nonoscillatory nature. Traditionally, the nonviscous modes (eigenvalues and eigenvectors) for nonproportional systems have been computed using the state-space approach, computationally expensive. This number of real eigenvalues is directly related to the rank of the damping matrices associated with the exponential kernels. The state-space approach has traditionally been used up to now as the only method to compute the nonviscous modes for nonproportionally damped systems. Motivated by this open problem, we propose in this chapter to describe the available numerical methods for classically damped systems and present the recent methods for nonclassically damped systems. It is shown that the problem of finding the nonviscous modes can be reduced to solve as a set of linear eigenvalue problems. The presented methods are compared through a numerical example.
- vibrating systems
- nonviscous damping
- eigenvalues and eigenvectors
- nonproportional systems
- numerical methods
1. Introduction and background
It has been always very difficult to model the physical fundamentals of damping in structural dynamics. In general, the proposed models depend on several parameters, which must be fitted according to experimental results. The viscous model, proposed by Rayleigh , is the most used representation of dissipative forces for vibrating systems as it predicts an exponential decay rate of displacements, something that can be observed experimentally in a great variety of structural materials such as metals, concrete, wood, glass, or masonry. However, damping models need to be updated for the mathematical modeling of the real behavior of the so-called viscoelastic damping materials, widely used for vibration control and energy dissipation devices. Although the term viscoelastic damping has traditionally been used, in the last years the concept nonviscous damping is also found in the bibliography, since this behavior can be considered as a generalization of the classic viscous damping. These materials, used in different areas of engineering as mechanical, civil, industrial, or aeronautics, are formed by polymer derivatives, rubbers, and rubber-like materials, and are characterized by a time-dependent constitutive model and by frequency-dependent Young and shear moduli.
Viscoelastic models of energy dissipation are introduced in the structure assuming that the damping forces are proportional to the history of the degrees-of-freedom (dof) velocities via kernel hereditary functions. These functions, also named damping functions, are the terms of the viscoelastic damping matrix in time domain, denoted by
where the dofs’ time-domain response is represented by
In this chapter, we will analyze Biot’s damping model with exponential kernels, a restriction commonly assumed in engineering applications. The expressions of the normalized damping functions in time and in frequency domain are, respectively
where with are the relaxation or nonviscous parameters and are the damping matrices (in general asymmetric) of the limit viscous model, obtained if the relaxation parameters tend to infinite, that is
The coefficients control the time and frequency dependence of the damping model while the spatial location and the level of damping are modeled via the matrices . From this property, it is not strange to find them also as the damping coefficients of Biot’s model. The level of damping is closely related to the magnitude of these limit-damping matrices, while the relaxation parameters give information on how far is our nonviscous model from a viscous behavior . It is also easy to demonstrate that the limit viscous damping and the time-domain kernel function are related by
Eqs. (2) and (4) clearly show the frequency dependence of the damping matrix, characteristic in this type of systems. This fact leads to a nonlinear eigenvalue problem whose eigenvalues are the roots of the equation
In general, the damping matrix admits a rational representation, so that the polynomial of the denominator is at least of one order less than that of the numerator . Thus, the determinant can be written as a polynomial, whose order is greater than and, therefore, the total number of roots of Eq. (7) can be expressed as and arranged as
where are complex conjugate pairs and are negative real numbers named nonviscous eigenvalues. The name is chosen precisely because they are characteristic of nonviscous or viscoelastic models. The number of these nonviscous eigenvalues will depend on the nature of the damping function, particularly on the number of hereditary exponential kernels. The complex conjugate pair forces the solution to be oscillatory, whereas the other eigenvalues are associated with overdamped, nonoscillatory modes. The latter modes decay rapidly and in general are not important for the system response.
The representation of the hereditary behavior was originally introduced by Boltzman  at the end of the nineteenth century. Its application to viscoelastic materials and to damping of vibrating systems was studied by different authors in the middle of the twentieth century. Among them, it is worth mentioning specially Biot [4, 5] whose multi-exponential hereditary model has widely been used for modeling viscoelastic damping materials. The fundamentals of viscoelasticity, a thorough study on the time-dependence constitutive models, and its application for modeling damping materials can be found in books such as Fluegge , Nashif , and Jones . Although this chapter is closely related to Biot’s damping model, we must not forget the other viscoelastic models based on the fractional derivatives and widely used for representing the frequency-dependent behavior of damping materials. This model allows to use less parameters than exponential-based models , although the mathematical treatment is more difficult to implement, especially in the time domain, which is computationally more expensive .
This chapter is focused on the study of the nonviscous modes of a nonviscously damped vibrating system. It is known that the effect of not considering these modes in the time-domain response is not important [2, 11]. Additionally, the exact calculation requires the use of the state-space approach, significantly increasing the computational effort and losing the physical insight of the involved internal variables [12, 13]. Maybe, these two reasons put together explain why they have not been analyzed in detail in the literature. Recently, Lázaro  derived a numerical approximated method to extract the nonviscous modes avoiding the state-space approach. In this chapter, we present a review of the nonviscous modes with nonoscillatory nature, giving their characteristics, mathematical properties, and the current available numerical methods for their computation.
2. Single degree-of-freedom systems
A single dof nonviscously damped vibrating system is dynamically characterized by a mass , a linear stiffness , and a nonviscous hereditary damping function
Checking solutions of the form
Multiplying this expression by , it results in a -order polynomial. If the system is lightly or moderately damped, the set of eigenvalues presents the form , where are a pair of conjugate-complex eigenvalues representing the modes of oscillatory nature. The rest roots are negative real numbers representing the nonviscous eigenvalues of nonoscillatory nature (overcritically damped). In this point, we focus on giving a mathematical characterization of these eigenvalues and to provide efficient methods to approximate the nonviscous eigenvalues avoiding solving the polynomial equation.
2.1. Mathematical characterization of eigenvalues
Let us see that the damping function evaluated at a nonviscous eigenvalue must always verify certain inequality related to the dynamic properties of the system, say mass and stiffness . Eq. (12) is rewritten in terms of the undamped natural frequency and of a new dimensionless damping function denoted by
The characteristic Eq. (11) now becomes
Reordering this equation, we can express it as
equivalent to . This inequality always holds for any real eigenvalue of Eq. (12). This result is a generalization of the well-known relationship between the dynamic parameters , , and of a single dof viscously damped oscillator for having real eigenvalues: (condition for critical damping).
As a direct consequence, we can define the following set:
assuring that every real eigenvalue of Eq. (12) lies inside
2.2. Numerical computation
It is known that the influence in the response of the nonviscous modes is much less important than that of the oscillatory complex modes [2, 16, 17]. For this reason, it is reasonable to look for closed-form approaches, avoiding the computational effort needed for solving the characteristic polynomial. Two methods based on the hypothesis of light damping can be found in the literature. They allow to approximate the nonviscous eigenvalues using closed-form formulas as function of the dynamic and damping parameters. The first one due to Adhikari and Pascual  approximates the nonviscous eigenvalues with the first iteration of Newton’s method applied to the characteristic polynomial. The second one, developed by Lázaro in his PhD Thesis  and published in the paper , is a perturbation-based approach. Both methods will be described in detail below and can be applied for both single dof systems and multiple dof systems with proportional (or classical) damping.
2.2.1. Adhikari and Pascual’s method
Let us denote by to the damping ratio associated to the th exponential kernel. Introducing these parameters, the characteristic equation can be written as
As mentioned before, the characteristic polynomial can be obtained multiplying the above equation by , resulting the -order polynomial
The method of Adhikari and Pascual  is based on the application of the first iteration of Newton’s method with as the initial point. Indeed, assuming that is close to the solution, can be explicitly calculated from the first-order expansion of around the initial point
After some simplifications, the expressions of Adhikari and Pascual published in Ref.  can be rewritten in terms of the current notation as
Under the hypothesis of light damping , the nonviscous eigenvalue lies close to ; therefore, it is expected that the solution from Eq. (21) accurately estimates the exact solution.
2.2.2. Lázaro’s method
Lázaro’s method [19, 20] is based on considering the th nonviscous eigenvalue , as a function of the th associated damping ratio . The damping ratio can be interpreted as a perturbation parameter of Eq. (18). Thus, we can write and Eq. (18) can be written for this eigenvalue as
Now, multiplying this equation by , we obtain
With this operation, the singularity associated to the th nonviscous eigenvalue can be avoided. The function
Eq. (24) explicitly defines as a function of . Assuming light damping, we can expand in terms of the damping parameter , considering the latter as a perturbation parameter within the equation. Thus,
The value can be obtained evaluating Eq. (24) at
The eigenvalue associated to the th hereditary kernel lies closely to the th relaxation parameter [15, 21]. Therefore, we are interested in the real solution . The first-order derivative can be calculated solving for after taking derivatives with respect to in Eq. (24). The rest of higher-order derivatives are derived following the same procedure and using the previously calculated results. In general, it is sufficient to take up the second-order term since this approximation accurately estimates the nonviscous eigenvalues within a wide range of the damping ratios, including lightly and moderately damped structures [15, 20]. After obtaining the coefficients and , the closed-form expression for remains as follows:
Both Lázaro’s and Adhikari and Pascual’s methods are presented as closed-form expressions. On one hand, numerical computation of polynomial roots is avoided, and on the other hand the analytical expressions allow to explicitly observe the dependence of the nonviscous eigenvalues as functions of the rest of the parameters of the problem.
3. Multiple degrees-of-freedom systems
This section deals with the properties of the nonviscous modes in asymmetric nonproportional viscoelastically damped vibrating systems. A generalization of the mathematical characterization proved for single dof systems in the previous point will be derived. Regarding numerical analysis, the available methods for computing nonviscous modes will also be presented. As mentioned in the introduction, we consider an -dof vibrating structure with mass and stiffness matrices denoted
The eigenvalues can be separated in conjugate-complex pairs with oscillatory nature and the real nonviscous eigenvalues . The number of nonviscous eigenvalues and the range of the damping matrices are directly related [12, 13]. In fact, Adhikari and Wagner  proved that
3.1. Mathematical characterization of eigenmodes
It is assumed that the damping matrix is not proportional, that is, does not verify the necessary conditions to be diagonal in the modal space of the undamped problem . As known, proportional damping matrices allow to reduce an -dof system to single dof systems due to the simultaneous decoupling capability. Thus, for these kinds of structures, the results of the previous point would apply. Assuming the nonproportionality, each nonviscous mode is characterized by a real eigenvalue
We define the following expressions for each nonviscous eigenmode:
These values can be interpreted as modal mass and stiffness, respectively, associated to the nonviscous mode. Using these new modal parameters, we can write that
We introduce functions
which can be interpreted as the dimensionless modal representation of the damping matrix at the th real mode. Introducing this relationship in Eq. (33), we obtain
We can identify in this equality the same form as that of Eq. (14), derived for single dof oscillators. Therefore, and using identical mathematical manipulations, we can deduce that
expression of which represents the generalization for multiple dof systems of the necessary condition derived for single dof systems in the previous point, Eq. (16). Additionally, Eq. (36) can also be considered as a generalization of the result published by Lázaro and Pérez-Aparicio  for symmetric systems.
3.2. The state-space approach
In this section, the general state-space representation of the dynamic problem will be described. This methodology allows to transform the general -dof system of integro-differential equations into a system of first-order differential equations through the introduction of internal variables. It was developed by Wagner and Adhikari  for symmetric systems and by Adhikari and Wagner  for asymmetric system.
It turns out that the final size of the extended state-space formulation is directly related to the rank of the damping matrices , . Because of that, it is appropriate to introduce the algebra associated to the matrix . Let us assume that
where and for are the right and left eigenvectors of the nonzero eigenvalues of . These two bases have special relevance in the developments of the state-space method and it is convenient to group them in the two rectangular matrices
so that the following relations are straightforward:
Let us return now to the system of integro-differential equations presented in Eq. (1) written in terms of the dof
For our purposes, we need the time derivative , which can be calculated using Leibniz’s rule for differentiation of an integral, yielding
With these new variables, Eq. (1) can be expressed as
In the above expression, the vector
Introducing this transformation into Eq. (43) and premultiplying by
Premultiplying by matrix and denoting by , we can write after some operations
In these expressions,
showing that the extra order of the state-space formulation of a nonviscously damped vibrating system is governed by the rank of the damping matrices. Hence, the total number of nonviscous eigenvalues is given by . As known, checking solutions of the form in the free-motion equations
The complete solution of this problem allows to construct the spectral set of nonviscously damped systems. On one hand, we have complex modes with oscillatory nature and, on the other hand, the nonviscous eigenmodes with their respective eigenvectors. A detailed study of the eigenvalue problem of Eq. (52) has been described in the work of Adhikari and Wagner . From a mathematical point of view, the problem of calculating the eigenmodes is totally solved. However, we can expose two reasons why it is worth to deepen in the numerical problem of the nonviscous modes: (a) to solve a linear eigenvalue problem as that shown in Eq. (52) requires in general
3.3. Approximate numerical method
As described above, we derive here the numerical method proposed by Lázaro  for the computation of nonviscous modes. We work under the generally accepted assumption of light damping, something that allows to predict that the nonviscous eigenvalues are close to the relaxation parameters . Let us consider the following decoupling of the damping matrix in the Laplace domain associated to the th relaxation parameter,
Something similar can be made for the dynamic stiffness matrix, yielding
where denotes the dynamic stiffness matrix without the th hereditary damping function. Note that under this manipulation, the function is now continuous and with continous derivatives at . Let us denote by
Let us define the matrix
and Eq. (57) can be written as
Since the damping is assumed to be light, is close to and consequently there exists certain , such that , with . Expanding the matrix around and neglecting second-order terms, we obtain
where . From the definition of , we have that
Substituting this result together with in Eq. (59) and rearranging
Following the same steps for the left eigenvectors from Eq. (56), we obtain the following relation between and
From Eqs. (62) and (64), , , and represent an eigensolution of the generalized linear asymmetric eigenvalue problem of matrix with respect to . Denoting by , then is eigenvalue of Eq. (26) with multiplicity . Consequently, there exist other no-null eigenvalues, which will be named . Hence, the complete spectral set of problem (64) can be listed as
and the nonviscous modes associated to the th relaxation parameter can be denoted by
We highlight two interesting results from this method: (i) the computation of the nonviscous modes has been reduced to solve linear eigenvalue problems of order and (ii) there is no need to previously calculate neither the modal space of the undamped model nor the eigenproblem of the matrices . We find a limitation because a hypothesis of light damping has been used in the linearization (60). For vibrating problems under a higher level of damping, the method can be adapted just taking the second-order term in the expansion of , that is
Introducing this expression in Eq. (59) and after some manipulations the resulting right and left eigenvalue problems are
In general, the second-order approximation will lead to better approximations, although in this case a larger problem must be solved; this will be confirmed in the numerical example. The reader who wants to deepen in detail in higher-order approximations and their associated computational cost can refer to the work of Lázaro . In this paper, it is proved that, from a computational point of view, it is profitable to increase the order of approximation up to certain limit order after which it is better to use the state-space approach. That limit value of the approximation order is
3.4. Numerical example
In this numerical example, the presented computational methods to calculate the nonviscous modes will be compared. For that, we use a five-degree-of-freedom discrete system with viscoelastic dampers, shown in Figure 1. Each dof represents the displacement of a mass kg. The linear stiffness between the masses is N/m. Two nonviscous dashpots are located between ground and dof 1 and between dofs 3, 4, and 5, whose constitutive relationships are expressed as the sum of exponential kernels.
The damping coefficients are Nms and Nms and the relaxation parameters rad/s. Since we have three relaxation parameters, the damping matrix in time domain yields
and according to the dashpots and rigidities distribution, the damping matrix coefficients and the stiffness matrix are
|1st order approx.||−9,767255360||−24,552342078||−44,490259517||−44,818527825|
|2nd order approx.||−9,762478350||−24,539554707||−44,480042888||−44,817178561|
|1st order approx.||0,924263325||0,992396067||0,000786787||0,001405038|
|2nd order approx.||0,920823393||0,991936757||0,000819965||0,001425786|
The rank of these matrices can easily be calculated obtaining
The number of nonviscous eigenvalues of this system is . The results of the four nonviscous eigenvalues and eigenvectors are shown in Table 1. Exact solutions based on the state-space approach are shown in the first rows. Below, we find the approximated solutions calculated with Lázaro’s method using both the first- and the second-order approximation (see Eqs. (62) and (68), respectively). The relative error is also shown below each result (in brackets) for both eigenvalues and eigenvectors. For the latter, the relative error is calculated in terms of the vector norms. Note that in general, the eigenvalues are calculated more accurately than eigenvectors. Indeed, the relative error of the former is one order of magnitude lower than that of the latter. As expected, the second-order approximation improves notably the solution, decreasing the relative errors two or three orders of magnitude respect to those computed from the first-order approximation. In general, since the effect of the nonviscous modes in the response is not relevant, it is justified to use the first-or second-order approximations presented in this text, even for moderately or highly damped vibrating structures .
In this chapter, the mathematical modeling of damping materials has been presented. These materials are characterized by presenting dissipative forces depending on the history of degrees-of-freedom velocities via exponential kernel functions (or Biot’s model). The free-motion vibration of these structural systems leads to a nonlinear eigenvalue problem. There exist two types of eigensolutions: on one hand, the complex eigenmodes, with oscillatory nature and considered as perturbations of the undamped natural modes, on the other hand, the so-called nonviscous modes, overcritically damped modes (without oscillatory nature), characteristic of the type of damping model. These latter modes are the main objective of the research of the present chapter.
The nonviscous modes behind a viscoelastic exponential-damping-based system are closely related to the relaxation parameter of the exponential functions. In general, their influence in the response of the system is several orders of magnitude less important than that of the complex modes. In this paper, we try to summarize some of the most relevant properties of these modes, both from a theoretical and from a numerical point of view. Nonviscous modes for both single and multiple dof systems are studied. For both cases, a necessary condition of nonviscous modes relating to eigenvector, eigenvalue, and dynamic matrices is provided. Additionally, numerical methods to extract nonviscous eigenvalues and eigenvectors, assuming asymmetric and nonproportional dynamic matrices, are reviewed. The results have been compared with a numerical example.