Sage Journals: Discover world-class research

Abstract

The complex damping model is only applied in frequency-domain. Based on the constitutive relation of viscoelastic materials, a vibration equation is equivalent to the constitutive relation of complex damping theory and the frequency response function satisfying the causality constraint. Compared with the existing equivalent complex damping model, the proposed equivalent damping model can consider more natural frequencies. The calculation formulas of Gauss precise integral method and improvement constant average acceleration method are derived. By the equivalent complex damping, the seismic time-history response of the typical example could carry on the numerical calculation. Compared with the Gauss precise integral calculation results of complex damping vibration equation, some conclusions can be analyzed. The proposed two methods could avoid divergent phenomenon in the time-domain numerical solution of complex damping vibration equation. The time-domain integral method of equivalent complex damping theory is stable and convergent. The Gauss precise integral method has strong equivalence and big computational complexity. The improvement constant average acceleration method has weak equivalence and small computational complexity.

Keywords

Equivalent complex damping time-domain integral Gauss precise integral method convergent

Introduction

Complex damping model is one of the common damping models,^1–3 which is widely used in the dynamic analysis of engineering structures. The complex damping model assumed that the damping force is proportional to the displacement in complex-domain.⁴ The advantage is that the energy dissipation is independent of frequency in harmonic steady-state vibration responses. The frequency-domain calculation results of the complex damping model are consistent with the experiment results of the steady-state due to harmonic waves.⁵ The damping matrix of the complex damping model depends only on the loss factor of the material and the stiffness matrix of the structures.⁶ Therefore, the advantage of the complex damping model is that the damping matrix is easier to construct. However, the complex damping model has some disadvantages. Firstly, the complex damping model is constructed in complex-domain, and the external excitation has to match the imaginary duality of the load.⁷ The calculated process is complex. Secondly, the free vibration responses are divergent.⁸ The complex damping model cannot be applied in time-domain, which is only in frequency-domain. Thirdly, the frequency response function is non-causal.^9–11 The divergent terms in the general solution is abandoned in time-domain method.^3,12,13 This method can obtain the stable results, however, which is essentially a subjective approximation algorithm. The calculating errors are difficult to estimate.

Compared with complex damping model, the viscoelastic damping model is the most common damping model and the frequency response function is causal.^14–16 The viscoelastic damping model assumed that the damping force can meet viscoelastic constitutive relation in real-domain. However, the damping matrix of the traditional viscoelastic damping model is difficult to construct, which needs to consider the modal damping ratios of the structures.¹⁷ The reason is that compared with the loss factor of complex damping model, the modal damping ratio is more difficult to obtain. It is important to study a new damping model, which combined with the advantages of complex damping model and viscoelastic damping model.

To ensure complex damping can be used in time-domain analysis, a series of viscoelastic damping models are obtained. Based on frequency-domain constitutive relation, some equivalent complex damping models of the single-degree-of-freedom systems are constructed, such as the Kelvin-Voigt damping model,¹⁸ Maxwell-Kelvin damping model,¹⁹ and Biot damping model.²⁰ The Kelvin-Voigt damping model can be further applied for the multi-degree-of-freedom (MDOF) systems, which only consider the natural frequency of the 1st mode. The Maxwell-Kelvin damping model and Biot damping model can consider multiple natural frequencies. However, the corresponding time-domain motion equation is difficult to solve and the calculated process is complex. It is necessary to introduce a simple damping model that considers more natural frequencies for the MDOF systems.

On this basis, the purpose of this paper is to propose an equivalent complex damping model, which can combine with the advantages of complex damping model and viscoelastic damping model. Based on the characteristic of equivalent complex damping model, the Gaussian precise integration method is obtained by aid of state place model. At the same time, combined with Newmark-β method, an improved constant average acceleration method is proposed based on equivalent complex damping model. Then the proposed time-domain integration methods are compared with the improved time-domain method and the Gaussian precise integration method based on traditional complex damping model in numerical examples. The calculated results show that the two time-domain integration methods can be used to calculate the convergent results.

An equivalent complex damping model based on viscoelastic material constitutive relation

Construction of equivalent complex damped motion equation

Some standard linear bodies are composed as a general linear body in series and parallel. For the MDOF systems, the corresponding differential equation of the constitutive relationship is

\sum_{i = 1}^{S} P_{i} \frac{d^{i} f (t)}{d t^{i}} = \sum_{j = 0}^{S} Q_{j} \frac{d^{j} x (t)}{d t^{j}}

(1)

where

f (t)

is the external excitation vector;

x (t)

is the displacement vector;

P_{1}, \dots, P_{S}

and

Q_{0}, Q_{1}, \dots, Q_{S}

are undetermined parameter matrix; and S is the total number of standard linear bodies

Based on Fourier transform method, equation (1) is expressed as

Ω + i ϖ P_{1} Ω + \dots + {(i ϖ)}^{S} P_{S} Ω = Q_{0} Γ + i ϖ Q_{1} Γ + \dots + {(i ϖ)}^{S} Q_{S} Γ

(2)

where

ϖ

is vibration frequency;

Ω

is Fourier transform of

f (t)

; and

Γ

is Fourier transform of

x (t)

Equation (2) is rewritten as

Ω = (U + i D) Γ

(3)

where

{\begin{cases} U = [(Q_{0} - ϖ^{2} Q_{2} + \dots) (1 - ϖ^{2} P_{2} + \dots) + (ϖ Q_{1} - ϖ^{3} Q_{3} + \dots) (ϖ P_{1} - ϖ^{3} P_{3} + \dots)] \\ {[{(1 - ϖ^{2} P_{2} + \dots)}^{2} + {(ϖ P_{1} - ϖ^{3} P_{3} + \dots)}^{2}]}^{- 1} \\ D = [- (Q_{0} - ϖ^{2} Q_{2} + \dots) (ϖ P_{1} - ϖ^{3} P_{3} + \dots) + (ϖ Q_{1} - ϖ^{3} Q_{3} + \dots) (1 - ϖ^{2} P_{2} + \dots)] \\ {[{(1 - ϖ^{2} P_{2} + \dots)}^{2} + {(ϖ P_{1} - ϖ^{3} P_{3} + \dots)}^{2}]}^{- 1} \end{cases}

(4)

The constitutive relation of complex damping model in the frequency-domain is

Ω = (1 + i v) K Γ

(5)

where

v

is loss factor;

K

is stiffness matrix; and

i = \sqrt{- 1}

Equation (3) is equivalent to equation (5). Then they are obtained as

{\begin{cases} K = U \\ v K = D \end{cases}

(6)

The elastic force does not change with vibration frequency, namely, $U$ does not change with vibration frequency. So they are obtained as

{\begin{cases} P_{1} = P_{2} = \dots = P_{S} = 0 \\ Q_{2} = Q_{4} = \dots = 0 \end{cases}

(7)

Then

K = Q_{0}

(8)

v K = ϖ Q_{1} - ϖ^{3} Q_{3} + ϖ^{5} Q_{5} - ϖ^{7} Q_{7} + \dots + ϖ^{2 S - 1} Q_{2 S - 1}

(9)

In order to approximate the complex damping model, the damping force does not change in a wide frequency range. However, it is difficult to guarantee equation (8) in a wide frequency range. The natural frequencies are adopted, and equation (9) is rewritten as

{\begin{cases} v K = ω_{1} Q_{1} - ω_{1}^{3} Q_{3} + ω_{1}^{5} Q_{5} - ω_{1}^{7} Q_{7} + \dots + ω_{1}^{2 S - 1} Q_{2 S - 1} \\ v K = ω_{2} Q_{1} - ω_{2}^{3} Q_{3} + ω_{2}^{5} Q_{5} - ω_{2}^{7} Q_{7} + \dots + ω_{2}^{2 S - 1} Q_{2 S - 1} \\ \dots \\ v K = ω_{S} Q_{1} - ω_{S}^{3} Q_{3} + ω_{S}^{5} Q_{5} - ω_{S}^{7} Q_{7} + \dots + ω_{S}^{2 S - 1} Q_{2 S - 1} \end{cases}

(10)

where

ω_{1}, \dots, ω_{S}

are natural frequencies.

It is assumed that

A = [\begin{array}{l} \begin{array}{l} ω_{1} & - ω_{1}^{3} & ω_{1}^{5} & - ω_{1}^{7} & \dots & ω_{1}^{2 S - 1} \end{array} \\ \begin{array}{l} ω_{2} & - ω_{2}^{3} & ω_{2}^{5} & - ω_{2}^{7} & \dots & ω_{2}^{2 S - 1} \end{array} \\ \dots \\ \begin{array}{l} ω_{S} & - ω_{S}^{3} & ω_{S}^{5} & - ω_{S}^{7} & \dots & ω_{S}^{2 S - 1} \end{array} \end{array}]

(11)

The natural frequencies are not equal to each other for most MDOF systems, and S is the number of degrees-of-freedoms. However, some natural frequencies are equal to each other for the systems with multiplex eigenvalues. The modal damping ratio of the structure is equal at the same natural frequency. So the equivalent complex damping model only considers different natural frequencies. S is defined as the number of non-multiplex eigenvalues. Here natural frequencies are not equal to each other. Then $A$ is invertible matrix.

Then equation (10) is expressed as

[\begin{array}{l} Q_{1} \\ Q_{3} \\ \dots \\ Q_{2 S - 1} \end{array}] = A^{- 1} [\begin{array}{l} v K \\ v K \\ \dots \\ v K \end{array}]

(12)

The constitutive relationship is expressed as

f (t) = Q_{0} x (t) + Q_{1} \frac{d x (t)}{d t} + Q_{3} \frac{d^{3} x (t)}{d t^{3}} + \dots + Q_{2 S - 1} \frac{d^{2 S - 1} x (t)}{d t^{2 S - 1}}

(13)

Based on equation (13), an equivalent damping model is obtained, and the corresponding time-domain motion equation is expressed as

M \frac{d^{2} x (t)}{d t^{2}} + Q_{0} x (t) + Q_{1} \frac{d x (t)}{d t} + Q_{3} \frac{d^{3} x (t)}{d t^{3}} + \dots + Q_{2 S - 1} \frac{d^{2 S - 1} x (t)}{d t^{2 S - 1}} = - M I u_{g} (t)

(14)

where

M

is the mass matrix;

u_{g} (t)

is the external excitation acceleration; and

I

is the distribution vector.

Initial values of motion equations

In the calculation process of structural dynamic response, only the initial displacement and initial velocity are adopted. So the initial values need to be further determined. When the natural frequency of the first mode is considered, the equivalent complex damping constitutive relationship is

{\begin{cases} Q_{0} = K \\ ω_{1} Q_{1, 1} = v K \end{cases}

(15)

Then the corresponding time-domain motion equation is expressed as

M \frac{d^{2} x (t)}{d t^{2}} + Q_{0} x (t) + Q_{1, 1} \frac{d x (t)}{d t} = - M I u_{g} (t)

(16)

When the natural frequencies of the first and second mode is considered, the equivalent complex damping constitutive relationship is

M \frac{d^{2} x (t)}{d t^{2}} + Q_{0} x (t) + Q_{1, 1} \frac{d x (t)}{d t} + Q_{2, 3} \frac{d^{3} x (t)}{d t^{3}} = - M I u_{g} (t)

(17)

By aid of trigonometric series, the external excitation acceleration can be expanded as

u_{g} (t) = A_{0} + \sum_{j = 1}^{N_{a} / 2} A_{j} \sin (θ_{j} t) + \sum_{j = 1}^{N_{a} / 2} B_{j} \cos (θ_{j} t)

(18)

Based on equation (18), it is further obtained as

\frac{d u_{g} (t)}{d t} = \sum_{j = 1}^{N_{a} / 2} θ_{j} A_{j} \cos (θ_{j} t) - \sum_{j = 1}^{N_{a} / 2} θ_{j} B_{j} \sin (θ_{j} t)

(19)

The derivations of both sides of equation (17) are taken, which is expressed as

M \frac{d^{3} x (t)}{d t^{3}} + Q_{0} \frac{d x (t)}{d t} + Q_{1, 1} \frac{d^{2} x (t)}{d t^{2}} + Q_{2, 3} \frac{d^{4} x (t)}{d t^{4}} = - M I \frac{d u_{g} (t)}{d t}

(20)

In the same way, the equivalences at natural frequencies of multiple vibration modes can be considered in turn. The corresponding equivalent complex damped motion equation and the corresponding derivation equation can be constructed, and the simultaneous equations can be obtained as

{\begin{cases} M \frac{d^{2} x (t)}{d t^{2}} + Q_{0} x (t) + Q_{1, 1} \frac{d x (t)}{d t} = - M I u_{g} (t) \\ M \frac{d^{2} x (t)}{d t^{2}} + Q_{0} x (t) + Q_{2, 1} \frac{d x (t)}{d t} + Q_{2, 3} \frac{d^{3} x (t)}{d t^{3}} = - M I u_{g} (t) \\ M \frac{d^{3} x (t)}{d t^{3}} + Q_{0} \frac{d x (t)}{d t} + Q_{2, 1} \frac{d^{2} x (t)}{d t^{2}} + Q_{2, 3} \frac{d^{4} x (t)}{d t^{4}} = - M I \frac{d u_{g} (t)}{d t} \\ M \frac{d^{2} x (t)}{d t^{2}} + Q_{0} x (t) + Q_{3, 1} \frac{d x (t)}{d t} + Q_{3, 3} \frac{d^{3} x (t)}{d t^{3}} + Q_{3, 5} \frac{d^{5} x (t)}{d t^{5}} = - M I u_{g} (t) \\ M \frac{d^{3} x (t)}{d t^{3}} + Q_{0} \frac{d x (t)}{d t} + Q_{3, 1} \frac{d^{2} x (t)}{d t^{2}} + Q_{3, 3} \frac{d^{4} x (t)}{d t^{4}} + Q_{3, 5} \frac{d^{6} x (t)}{d t^{6}} = - M I \frac{d u_{g} (t)}{d t} \\ \dots \\ M \frac{d^{2} x (t)}{d t^{2}} + Q_{0} x (t) + Q_{S - 1, 1} \frac{d x (t)}{d t} + Q_{S - 1, 3} \frac{d^{3} x (t)}{d t^{3}} + \dots + Q_{S - 1, 2 S - 3} \frac{d^{2 S - 3} x (t)}{d t^{2 S - 3}} = - M I u_{g} (t) \\ M \frac{d^{3} x (t)}{d t^{3}} + Q_{0} \frac{d x (t)}{d t} + Q_{S - 1, 1} \frac{d^{2} x (t)}{d t^{2}} + Q_{S - 1, 3} \frac{d^{4} x (t)}{d t^{4}} + \dots + Q_{S - 1, 2 S - 3} \frac{d^{2 S - 2} x (t)}{d t^{2 S - 2}} = - M I \frac{d u_{g} (t)}{d t} \end{cases}

(21)

x (t_{0})

and

\frac{d x (t_{0})}{d t}

are substituted into equation (21).

\frac{d^{2} x (t_{0})}{d t^{2}}, \frac{d^{3} x (t_{0})}{d t^{3}}, \frac{d^{4} x (t_{0})}{d t^{4}}, \dots, \frac{d^{2 S - 3} x (t_{0})}{d t^{2 S - 3}}

can be obtained.

An improved Gaussian precise integration method

State space equation

The Gaussian precise integration method has the advantages of high accuracy and short step length.²¹ An improved Gaussian precise integration method is used to solve the equivalent complex damping equation. Firstly, based on state space method,²² the state space equation of the equivalent complex damping model is constructed. Auxiliary equations are introduced. Based on equation (14), the equivalent equations are

{\begin{cases} \frac{d x (t)}{d t} = \frac{d x (t)}{d t} \\ \frac{d^{2} x (t)}{d t^{2}} = \frac{d^{2} x (t)}{d t^{2}} \\ \frac{d^{3} x (t)}{d t^{3}} = \frac{d^{3} x (t)}{d t^{3}} \\ \frac{d^{4} x (t)}{d t^{4}} = \frac{d^{4} x (t)}{d t^{4}} \\ \frac{d^{5} x (t)}{d t^{5}} = \frac{d^{5} x (t)}{d t^{5}} \\ \dots \\ \frac{d^{2 S - 1} x (t)}{d t^{2 S - 1}} = - Q_{2 S - 1}^{- 1} Q_{0} x (t) - Q_{2 S - 1}^{- 1} Q_{1} \frac{d x (t)}{d t} - Q_{2 S - 1}^{- 1} M \frac{d^{2} x (t)}{d t^{2}} - Q_{2 S - 1}^{- 1} Q_{3} \frac{d^{3} x (t)}{d t^{3}} \\ \dots - Q_{2 S - 1}^{- 1} Q_{2 S - 3} \frac{d^{2 S - 3} x (t)}{d t^{2 S - 3}} - Q_{2 S - 1}^{- 1} M I u_{g} (t) \end{cases}

(22)

Equation (22) is rewritten as

W^{(1)} = H W + R

(23)

where

W = {[\begin{array}{l} x (t) & \frac{d x (t)}{d t} & \frac{d^{2} x (t)}{d t^{2}} & \dots & \frac{d^{2 S - 2} x (t)}{d t^{2 S - 2}} \end{array}]}^{T}

(24)

R = {[\begin{array}{l} 0 & 0 & 0 & \dots & - Q_{2 S - 1}^{- 1} M I u_{g} (t) \end{array}]}^{T}

(25)

H = [\begin{array}{l} 0 & E & 0 & 0 & 0 & 0 & \dots & 0 \\ 0 & 0 & E & 0 & 0 & 0 & \dots & 0 \\ 0 & 0 & 0 & E & 0 & 0 & \dots & 0 \\ 0 & 0 & 0 & 0 & E & 0 & \dots & 0 \\ 0 & 0 & 0 & 0 & 0 & E & \dots & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & L & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & E \\ - Q_{2 S - 1}^{- 1} Q_{0} & - Q_{2 S - 1}^{- 1} Q_{1} & - Q_{2 S - 1}^{- 1} M & - Q_{2 S - 1}^{- 1} Q_{3} & 0 & - Q_{2 S - 1}^{- 1} Q_{5} & \dots & - Q_{2 S - 1}^{- 1} Q_{2 S - 3} \end{array}]

(26)

E

is the unit matrix.

The general solution of equation (23) is expressed as

W = e^{H t} W_{0} + \int_{0}^{t} e^{H (t - τ)} R (τ) d τ

(27)

where

W_{0} = {[\begin{array}{l} x (t_{0}) & \frac{d x (t_{0})}{d t} & \frac{d^{2} x (t_{0})}{d t^{2}} & \dots & \frac{d^{2 S - 2} x (t_{0})}{d t^{2 S - 2}} \end{array}]}^{T}

(28)

Gaussian precise integration method

The whole time-history is discrete according to the time step $Δ t$ , any time can be expressed as

t_{k} = k Δ t (k = 0, 1, 2, \dots)

(29)

At the time $t_{k}$ , the time-domain response vector is

W_{k} = e^{H t_{k}} W_{0} + \int_{0}^{t_{k}} e^{H (t_{k} - τ)} R (τ) d τ

(30)

At the time $t_{k + 1}$ , the time-domain response vector is

W_{k + 1} = e^{H Δ t} W_{k} + \int_{t_{k}}^{t_{k + 1}} e^{H (t_{k + 1} - τ)} R (τ) d τ

(31)

Based on three-point Gaussian integral format, the integral term at the right side of equation (31) is expressed as

\int_{t_{k}}^{t_{k + 1}} e^{H (t_{k + 1} - τ)} R (τ) d τ = \sum_{i = 1}^{3} A_{i} T_{i} R (t_{k} + 0.5 Δ t (1 + y_{i})) 0.5 Δ t

(32)

where

{\begin{cases} A_{1} = \frac{5}{9}, T_{1} = e^{0.5 Δ t (1 + \sqrt{0.6}) H}, y_{1} = - \sqrt{0.6} \\ A_{2} = \frac{8}{9}, T_{2} = e^{0.5 Δ t H}, y_{2} = 0 \\ A_{3} = \frac{5}{9}, T_{3} = e^{0.5 Δ t (1 - \sqrt{0.6}) H}, y_{3} = - \sqrt{0.6} \end{cases}

(33)

e^{H Δ t}

e^{0.5 Δ t (1 + \sqrt{0.6}) H}

e^{0.5 Δ t H}

, and

e^{0.5 Δ t (1 - \sqrt{0.6}) H}

are calculated based on the Gaussian precise integration method. Therefore, an improved Gaussian precise integration method can be realized based on the equivalent complex damping model.

An improved Newmark-β method

In the numerical method of dynamic responses, the Newmark-β method is widely adopted, which has the advantages of simple principle. Especially, the constant average acceleration method is based on the Newmark-β method, which is unconditionally convergent.

It is assumed that

\frac{d^{2} x (τ)}{d t^{2}} = \frac{1}{2} (\frac{d^{2} x (t_{k + 1})}{d t^{2}} + \frac{d^{2} x (t_{k})}{d t^{2}}) (t_{k} < τ < t_{k + 1})

(34)

Then the assumption is further expressed as

\frac{d^{2 S - 1} x (τ)}{d t^{2 S - 1}} = \frac{1}{2} (\frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}} + \frac{d^{2 S - 1} x (t_{k})}{d t^{2 S - 1}}) (t_{k} < τ < t_{k + 1})

(35)

Based on equation (35), they are obtained as

{\begin{cases} \frac{d^{2 S - 2} x (t_{k + 1})}{d t^{2 S - 2}} = \frac{d^{2 S - 2} x (t_{k})}{d t^{2 S - 2}} + \frac{1}{2} Δ t (\frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}} + \frac{d^{2 S - 1} x (t_{k})}{d t^{2 S - 1}}) \\ \frac{d^{2 S - 3} x (t_{k + 1})}{d t^{2 S - 3}} = \frac{d^{2 S - 3} x (t_{k})}{d t^{2 S - 3}} + Δ t \frac{d^{2 S - 2} x (t_{k})}{d t^{2 S - 2}} + \frac{1}{4} Δ t^{2} (\frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}} + \frac{d^{2 S - 1} x (t_{k})}{d t^{2 S - 1}}) \\ \dots \\ \frac{d^{2} x (t_{k + 1})}{d t^{2}} = \frac{d^{2} x (t_{k})}{d t^{2}} + Δ t \frac{d^{3} x (t_{k})}{d t^{3}} + \dots + \frac{1}{2} \frac{1}{(2 S - 3)!} Δ t^{2 S - 3} (\frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}} + \frac{d^{2 S - 1} x (t_{k})}{d t^{2 S - 1}}) \\ \frac{d x (t_{k + 1})}{d t} = \frac{d x (t_{k})}{d t} + Δ t \frac{d^{2} x (t_{k})}{d t^{2}} + \dots + \frac{1}{2} \frac{1}{(2 S - 2)!} Δ t^{2 S - 2} (\frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}} + \frac{d^{2 S - 1} x (t_{k})}{d t^{2 S - 1}}) \\ x (t_{k + 1}) = x (t_{k}) + Δ t \frac{d x (t_{k})}{d t} + \dots + \frac{1}{2} \frac{1}{(2 S - 1)!} Δ t^{2 S - 1} (\frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}} + \frac{d^{2 S - 1} x (t_{k})}{d t^{2 S - 1}}) \end{cases}

(36)

At the time $t_{k + 1}$ , the time-domain equation is expressed as

M \frac{d^{2} x (t_{k + 1})}{d t^{2}} + Q_{0} x (t_{k + 1}) + Q_{1} \frac{d x (t_{k + 1})}{d t} + Q_{3} \frac{d^{3} x (t_{k + 1})}{d t^{3}} + \dots + Q_{2 S - 1} \frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}} = - M I u_{g} (t_{k + 1})

(37)

Equation (36) is substituted into equation (37). $\frac{d^{2 S - 1} x (t_{k + 1})}{d t^{2 S - 1}}$ can be obtained based on $x (t_{k}), \frac{d x (t_{k})}{d t}, \frac{d^{2} x (t_{k})}{d t^{2}}, \dots, \frac{d^{2 S - 2} x (t_{k})}{d t^{2 S - 2}}$ . Then $x (t_{k + 1}), \frac{d x (t_{k + 1})}{d t}, \frac{d^{2} x (t_{k + 1})}{d t^{2}}, \dots, \frac{d^{2 S - 2} x (t_{k + 1})}{d t^{2 S - 2}}$ can be obtained from equation (37). An improved Newmark-β method can be realized based on the equivalent complex damping model.

Numerical examples

The five-story steel frame model is taken as example. The mass and stiffness distribution are shown in Figure 1. The natural frequencies are shown in Table 1. The damping ratio of steel is 0.02.²³ The loss factor is twice the damping ratio, namely, the loss factor is 0.04.²⁴

Figure 1.

The mass and stiffness distribution of steel frame.

Table 1.

Natural frequencies of steel frame.

Mode	1	2	3	4	5
Natural frequency (rad/s)	2.9706	7.1805	11.0994	14.1601	16.0076

The mass matrix is

M = 10^{3} kg \times [\begin{array}{l} 1.0 & 0 & 0 & 0 & 0 \\ 0 & 1.5 & 0 & 0 & 0 \\ 0 & 0 & 2.0 & 0 & 0 \\ 0 & 0 & 0 & 2.5 & 0 \\ 0 & 0 & 0 & 0 & 3.0 \end{array}]

(38)

The stiffness matrix is

K = 10^{5} N / m \times [\begin{array}{l} 1.0 & - 1.0 & 0 & 0 & 0 \\ - 1.0 & 2.2 & - 1.2 & 0 & 0 \\ 0 & - 1.2 & 2.7 & - 1.5 & 0 \\ 0 & 0 & - 1.5 & 3.3 & - 1.8 \\ 0 & 0 & 0 & - 1.8 & 3.8 \end{array}]

(39)

The damping matrix of the complex damping model is

v K = 10^{4} N / m \times [\begin{array}{l} 0.40 & - 0.40 & 0 & 0 & 0 \\ - 0.40 & 0.88 & - 0.48 & 0 & 0 \\ 0 & - 0.48 & 1.08 & - 0.60 & 0 \\ 0 & 0 & - 0.60 & 1.32 & - 0.72 \\ 0 & 0 & 0 & - 0.72 & 1.52 \end{array}]

(40)

According to the mass participation coefficient method,²⁵ only considering the first three vibration modes of the five-story steel frame mode can meet the requirements of calculation accuracy. Based on the natural frequencies of the first three vibration modes, the parameters of equivalent complex damping model can be determined from equation (14). The proposed improved Gaussian precise integration method (IGEC) and the improved Newmark-β method (INEC) based on equivalent complex damping model are adopted. Besides, the improved time-domain method (IT)^3,13 and the Gaussian precise integration method (GC) based on traditional complex damping model²⁶ are also adopted. IT is a simplified calculation method, which subjectively removes the divergent term in the general solution and is equal to the frequency-domain results of complex damping model. IT can be regarded as a reference method. GC is the time-domain method, which does not remove the divergent term in the general solution. The external excitation is selected as earthquake wave. The information of earthquake waves is shown in Table 2.

Table 2.

Information of earthquake waves.

Earthquake wave	Peak acceleration (g)	Main frequency (rad/s)	Time (s)
Morgan hill	0.1419	2.4544	38
Coalinga Cantua	0.2885	10.2777	35
Coalinga Parkfield	0.0511	16.4136	30

Due to the three earthquake waves, IGEC, INEC, IT, and GC are used to calculate the time-history displacements of frame model, and the results are shown in Figures 2–4. With the increase of earthquake action time, the calculated results based on GC are divergent. Figures 2–4 show that GC is only applied for the cases with short time earthquake action time. Due to Morgan Hill wave, when the earthquake action time is greater than 32 s, the divergent phenomenon occurs in the time-history displacements as shown in Figure 2. Due to Coalinga Cantua wave, when the earthquake action time is greater than 30 s, the divergent phenomenon occurs in the time-history displacements as shown in Figure 3. Due to Coalinga Parkfield wave, when the earthquake action time is greater than 24 s, the divergent phenomenon occurs in the time-history displacements as shown in Figure 4. Therefore, the complex damping model cannot be applied to the time-domain analysis.

Figure 2.

The top time-history displacements with different methods due to Morgan Hill wave. (a) Acceleration of earthquake wave and (b) top time-history displacements.

Figure 3.

The top time-history displacements with different methods due to Coalinga Cantua wave. (a) Acceleration of earthquake wave and (b) top time-history displacements.

Figure 4.

The top time-history displacements with different methods due to Coalinga Parkfield wave. (a) Acceleration of earthquake wave and (b) top time-history displacements.

Compared with GC, the time-history displacements based on IGEC, INEC, and IT are convergent. The comparative analysis of the calculated results is shown in Table 3. The calculated results based on IGEC and INEC are consistent with the ones based on IT. The relative errors are all smaller than 5%. The maximum relative error based on IGEC is 0.18%. The maximum relative error based on INEC is 4.63%. So the correctness of IGEC and INEC is proved. Besides, IGEC and INEC are different from IT. The divergent terms in the general solution is subjectively abandoned in the calculated process of IT. The calculated results based on IGEC and INEC are objective.

Table 3.

Comparative analysis of the calculated results with different methods.

Methods	IGEC	INEC	IT
Earthquake wave	Morgan Hill
Peak time (s)	18.25	16.95	18.26
Peak displacement (mm)	220.88	210.54	220.75
Relative error (%)	0.06	4.63	—
Earthquake wave	Coalinga Cantua
Peak time (s)	15.91	15.92	15.93
Peak displacement (mm)	195.79	187.67	195.76
Relative error (%)	0.02	4.13	—
Earthquake wave	Coalinga Parkfield
Peak time (s)	9.01	9.01	9.01
Peak displacement (mm)	33.10	32.25	33.16
Relative error (%)	0.18	2.74	—

Compared with INEC, the relative errors based on IGEC are smaller. IGEC has stronger equivalence and higher accuracy. Especially, due to Morgan Hill wave, the time of peak displacement based on INEC is obviously different from the one based on IGEC and IT. The reason is that IGEC assumes that the acceleration within a time step changes nonlinearly and INEC assumes that the acceleration within a time step changes linearly. The vibration responses change nonlinearly due to earthquake waves. Therefore, the calculation error of improved Gaussian precise integration method is less and the corresponding accuracy is higher. However, the calculated process of INEC is more simple, which does not involve the state space and complex integral calculation. Namely, the computational efficiency of INEC is higher.

Conclusions

Based on the constitutive relationship of viscoelastic materials, some main vibration modes are considered and an equivalent complex damping model is proposed. Then the corresponding Gaussian precise integration method and Newmark-β method are obtained. Numerical examples validate the effectiveness of the numerical methods. Some conclusions are as follows.

(1) With the increase of action time of external excitation, the calculated results of numerical method for complex damped systems are divergent. Compared with this method, the calculated results of the proposed numerical method for equivalent complex damped systems are stable and convergent, and can effectively avoid divergence in the time-domain numerical solution of the complex damping motion equation.

(2) Compared with the improved Newmark-β method, the improved Gaussian precise integration method assumes that the acceleration within a time step changes nonlinearly, which can better fit nonlinear vibration response. So the numerical results of the improved Gaussian precise integration method have higher accuracy. However, the state space equation will increase the dimension of the matrix. The calculated process of the improved Gaussian precise integration method is more complex.

In this paper, the equivalent complex damping model can overcome the non-causal shortcoming and consider more natural frequencies. However, when the high modes are introduced, the calculated process will be complex. In the future work, a more simplified equivalent damping model that can consider all vibration frequencies within the interesting frequency range needs further to be analyzed, and some complex and high-precision numerical methods will be studied.

Footnotes

Acknowledgments

The ground-motion data was obtained from the PEER Ground Motion Database ().

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The work was supported by the National Natural Science Foundation of China (Grant No. 52208322), Postdoctoral Research Foundation of China (Grant No. 2022M712905), Henan Province Science and Technology Research Funding Project (No. 222102230010), and Key Research Project of Henan Higher Education Institutions (Grant No. 22A560005).

ORCID iDs

Shuxia Wang

Panxu Sun

Data availability statement

The data that support the findings of this study are available from the corresponding author upon reasonable request.*

References

Yuan

Liang

. Dynamic properties of composite panel sandwiched by a perforated fiber reinforced damping film. Polym Compos 2022; 43(2): 753–767.

Chen

Zhai

Han

. Vibration and damping analysis of the bladed disk with damping hard coating on blades. Aero Sci Technol 2016; 58: 248–257.

Wang

Zhang

, et al. Dynamic characteristics of blade with viscoelastic damping block based on complex eigenvalue method. Shock Vib 2018; 2018: 1–16.

Miklestad

. The concept of complex damping. J Appl Mech 1952; 19(3): 284–286.

Bert

. Material damping: an introductory review of mathematical models, measures and experimental techniques. J Sound Vib 1973; 29(2): 129–153.

Sun

Yang

Deng

. Complex mode superposition method of nonproportionally damped linear systems with hysteretic damping. J Vib Control 2021; 27(13-14): 1453–1465.

Inaudi

Kelly

. Linear hysteretic damping and the Hilbert transform. J Eng Mech 1995; 121(5): 626–632.

Wang

Fan

, et al. A time-domain method for free vibration responses of an equivalent viscous damped system based on a complex damping model. J Low Freq Noise, V A 2023; 42: 146134842311575.

Nakamura

. Practical causal hysteretic damping. Earthq Eng Struct Dynam 2007; 36(5): 597–617.

10.

Makris

. Causal hysteretic element. J Eng Mech 1997; 123(11): 1209–1214.

11.

Nakamura

. Improved methods to transform frequency-dependent complex stiffness to time domain. Earthq Eng Struct Dynam 2006; 35(8): 1037–1050.

12.

Barkanov

Hufenbach

Kroll

. Transient response analysis of systems with different damping models. Comput Methods Appl Math 2003; 192(1): 33–46.

13.

Pan

. The direct integration method with virtual Initial conditions on the free and forced vibration of a system with hysteretic damping. Appl Sci-basel 2019; 9(18): 3707.

14.

Chopra

. Dynamics of structures: theory and applications to earthquake engineering. New Jersey: Prentice-Hall, 2007.

15.

Jangid

. Response of SDOF system to non-stationary earthquake excitation. Earthq Eng Struct Dynam 2004; 33(15): 1417–1428.

16.

Jangid

. Equivalent linear stochastic seismic response of isolated bridges. J Sound Vib 2008; 309(3–5): 805–822.

17.

Henwood

. Approximating the hysteretic damping matrix by a viscous matrix for modelling in the time domain. J Sound Vib 2002; 254(3): 575–593.

18.

Sun

Yang

Kang

. Time-domain analysis for dynamic responses of non-classically damped composite structures. Compos Struct 2020; 251: 112554.

19.

Reggio

De Angelis

. Modelling and identification of structures with rate-independent linear damping. Meccanica 2015; 50(3): 617–632.

20.

Lin

Rao

. Vibration analysis of a multiple-layered viscoelastic structure using the Biot damping model. AIAA J 2010; 48(3): 624–634.

21.

Zhong

Cai

. Precise integration method for LQG optimal measurement feedback control problem. Appl Math Mech 2000; 21(12): 1417–1422.

22.

Foss

. Coordinates which uncouple the equations of motion of damped linear dynamic systems. J Appl Mech 1958; 25: 361–364.

23.

Huang

Qian

Zhou

, et al. An approach to equivalent damping ratio of vertically mixed structures based on response error minimization. Soil Dynam Earthq Eng 2015; 72: 119–128.

24.

Lin

Zhu

. On the relationship between viscous and hysteretic damping models and the importance of correct interpretation for system identification. J Sound Vib 2009; 325(1-2): 14–33.

25.

Wilson

. Three dimensional static and dynamic analysis of structures: a physical approach with emphasis on earthquake engineering. Berkeley, California, USA: Computers and Structures, 1999.

26.

Pan

Wang

. Gauss precise time-integration of complex damping vibration systems. Eng Mech 2012; 29(2): 16–20, (In Chinese).

A time-domain integration method with equivalent complex damping model based on viscoelastic material constitutive relation

Abstract

Keywords

Introduction

An equivalent complex damping model based on viscoelastic material constitutive relation

Construction of equivalent complex damped motion equation

Initial values of motion equations

An improved Gaussian precise integration method

State space equation

Gaussian precise integration method

An improved Newmark-β method

Numerical examples

Conclusions

Footnotes

Acknowledgments

Declaration of conflicting interests

Funding

ORCID iDs

Data availability statement

References