Home | | Structural Dynamics and Earthquake Engineering | Numerical solution methods for natural frequencies and mode shapes

# Numerical solution methods for natural frequencies and mode shapes

Numerical solution methods for natural frequencies and mode shapes in relation to structural dynamics during earthquakes: Abstract: This chapter discusses basic solution schemes as well as approximate methods for finding natural frequencies and mode shapes. The methods include VianelloŌĆ'Stoodala power method, transfer matrix method, Jacobi method, Holzer method, RayleighŌĆÖs approximation and DunkerleyŌĆÖs approximation. Some of the approximate methods will lead to upper bound solutions and some lower bound solutions. A relevant program in MATHEMATICA is also given.

Numerical solution methods for natural frequencies and mode shapes in relation to structural dynamics during earthquakes

Abstract: This chapter discusses basic solution schemes as well as approximate methods for finding natural frequencies and mode shapes. The methods include VianelloŌĆ'Stoodala power method, transfer matrix method, Jacobi method, Holzer method, RayleighŌĆÖs approximation and DunkerleyŌĆÖs approximation. Some of the approximate methods will lead to upper bound solutions and some lower bound solutions. A relevant program in MATHEMATICA is also given.

Key words: banded matrix, sweeping technique, deflation, transfer matrix, Holzer method.

Introduction

The important step in the dynamic analysis of a multiple-degrees-of-freedom (MDOF) system is quite often the solution of the eigenvalue problem, or the determination of the system natural frequencies and the corresponding normal vibration modes. This is particularly true if a mode superposition analysis is to be conducted. Several procedures for solving the eigenvalue problem have been discussed in many books. Predicting or finding the roots of the characteristic polynomial as discussed in Chapter 10 is satisfactory only for systems having few degrees of freedom. For large MDOF systems, extracting the roots of the characteristic polynomial requires computational effort and is quite often an indeterminable task. This chapter discusses the basic solution schemes as well as approximate methods for finding frequencies.

General solution methods for eigen problems

In structural dynamics, the basic eigen problem for an MDOF system having ŌĆśnŌĆÖ degrees of freedom is represented as

[k]{Žå} = ╬╗[m]{Žå}                        ŌĆ”  . . . . 11.1

Let [k] be the stiffness matrix of order ŌĆśnŌĆÖ and [m], the mass matrix, also of order ŌĆśnŌĆÖ. For most structural systems, [k] is normally banded matrix and [m] is a diagonal matrix for a lumped mass formulation (without rotary inertia) coupling or a narrowly banded matrix for a consistent mass formulation.

There are ŌĆśnŌĆÖ eigenvalues and ŌĆśnŌĆÖ eigenvectors satisfying the above equation. The rth eigen pair is determined by (╬╗r,{Žå}r). In dynamic problems, the eigenvalues are the square of the natural frequencies Žē n2 such that

0 < ╬╗1 < ╬╗2 <ŌĆ”< ╬╗n           ŌĆ”ŌĆ”ŌĆ”.  .. . 11.2

The dynamic response of MDOF systems having a large number of degrees of freedom is generally confined to a relatively small subset of the lowest vibration modes of the system. Therefore for such systems, only ŌĆśpŌĆÖ eigen pairs need to be solved for where p << n. The solution of p eigen values and the corresponding eigenvectors can be written as The majority of the eigen problem solution techniques can be classified as

ŌĆó       vector iteration methods;

ŌĆó       transformation methods;

ŌĆó       polynomial iteration methods.

Clearly all the methods are iterative in nature because the solution of the eigen problem as defined in Eq. 11.1 is tantamount to solving the characteristic polynomial of order ŌĆśnŌĆÖ. Since explicit formulas for the determination of roots to the characteristic polynomial having an order higher than 4 do not exist, an iterative solution is mandatory.

The main essence of each method is very distinctive. The vector iteration methods are based on the property Let [K] and [M] be modal stiffness and modal mass respectively. The polynomial iteration is based on the property that the characteristic polynomial is a function of ╬╗r, and The characteristic polynomial is of order ŌĆśnŌĆÖ. Solution of characteristic polynomial has been discussed in Previous Pages.

Vector iteration technique

1  Vianello and Stoodala method (power method)

This iterative method can be applied to extract the highest eigenvalue of either symmetric or unsymmetric matrix of any order. Consider the homogeneous equation written in the form of Assume any vector {x}0 and multiplying with [A] matrix gives {y}1 can be written in terms of ╬╗1{x}1 by taking highest element (associated vector) outside and this can be used in the next iteration as It can be seen that the eigenvalue ╬╗ as well as eigenvector will converge. The iteration can be stopped when |╬╗n+1 ŌĆ' ╬╗n| < ╬Ą (tolerance) and ╬╗ is the highest eigenvalue in the absolute sense.

Consider the problem in dynamics [D]1 is called the dynamic matrix. Application of the power method will converge to (1/╬╗)max or ╬╗min = ╬╗1 ┬Ę ╬╗1 is called the least dominant eigenvalue and {ŽĢ}1 is called least dominant eigenvector. However the method may be modified to calculate eigenvalues and the corresponding eigenvectors for higher modes by matrix deflation or deflation of the iteration vectors.

Assume using the power method we get the (1/╬╗)max or ╬╗min = ╬╗1 and {Žå}1 pair. The basic premise for vector deflation is that, for an iteration vector to converge to a required eigenvector, the iteration vector must be orthogonal to the eigenvector. Therefore, for this case at hand, this can be interpreted as meaning that if the iterative vector is orthogonalized to the eigenvectors already calculated (for example {Žå}1) the vector is precluded and occurs instead to another (higher) eigenvector. More succinctly, an eigen pair other than (╬╗1, {Žå}1) becomes the least dominant eigen pair.

2  Method 1 sweeping technique

To find the second eigenvalue assume a trial vector that is a linear combination of all eigenvectors. In the trial vector assumed in Eq. 11.13, let us sweep out the effect of mode 1 as If {Žå}1 is the normalized eigenvector, as proved in Chapter 10, {Žå }1T [ m]{Žå}1 = 1 because of the normalization principle. In that case dynamic matrix is written as [D]2 = [D]1[S]1 where [S]1 is given as in Eq. 11.18 as In Eq. 11.19, [S]1 is called sweeping matrix for the first mode. The purpose of [S]1 is to eliminate or sweep out the effect of mode 1. Similarly [S]n is to sweep out the effect of modes for 1 to n and allow the mode (n + 1) to become least dominant.

3  Method 2 ŌĆ' deflation method

Assume we find (╬╗1, {Žå}1) pair to obtain some eigenvalue using the ortho-normal principle. The power method will converge to the second eigenvalue and the corresponding normalized eigenvector {ŽĢ}2 can be obtained.

To extract the third eigenvalue, Using the power method, the iteration will converge to third eigenvalue and the corresponding normalized eigenvector can be obtained.

Let us use the same example and assume we obtain the highest eigenvalue as 2.8453 and the corresponding normalized eigenvector as given in method 1. Substituting the corresponding values in Eq. 11.25, we get  Using the first normalized eigenvector we obtain the equation

x1 = ŌĆ'1.611x2 ŌĆ' 1.656x(a)

Using the second normalized eigenvector we obtain the equation The procedure can be programmed very easily in the EXCEL package.

JacobiŌĆÖs method

The matrix iteration method discussed in Section 11.3 produces the eigenvalues and eigenvectors of matrix [D]1 one at a time. JacobiŌĆÖs method is also an iterative method but produces all the eigenvalues and eigenvectors of the matrix [D]1 simultaneously. [D]1 is a real symmetrix matrix and has only real eigenvalues. There is an orthogonal matrix [R] such that [R]T[D]1[R] is a diagonal matrix. The diagonal elements are the eigenvalues and the columns of [R] is generated as a product of several rotation matrices as

[R] = [R1][R2][R3]ŌĆ”  ---- 11.36

Consider the highest off-diagonal term of matrix [D]1. Let it be dij. Find Again find the highest off-diagonal term and form [R2] matrix. While making this off-diagonal term as zero, it introduces non-zero contributions to formerly zero positive. However, successive matrices of the form

[R2]T[R1]T[D]1[R1][R2]                    - - - - - - 11.41

[R3]T[R2]T[R1]T[D]1[R1][R2][R3]                 - - - - - - 11.42

converges to the required diagonal form. Find matrix [R] such that

[R] = [R1][R2][R3]ŌĆ”[Rn]               - - - - - - 11.43

is the required eigenvector.

Consider Example 11.1 with mass matrix as [I]. The dynamic matrix is written as  Transfer matrix method to find the fundamental frequency of a multi-storeyed building (shear frame)

Consider a shear frame shown in Fig. 11.1 consisting of ŌĆśnŌĆÖ storeys. Consider the free body diagram of the ith storey as shown in Fig. 11.2. Here ŌĆśmŌĆÖ is the mass, ŌĆśkŌĆÖ is the stiffness, ŌĆśVŌĆÖ is the shear, ŌĆśŽēŌĆÖ is the natural frequency and ŌĆśvŌĆÖ is the displacement. Considering inertia force, the shear in the ith storey can be written in terms of shear of (i+1)th storey. Initially the displacement at the top storey level and the natural frequency are both assumed. By using the transfer matrix method, it is possible to find the displacement at the base as well as shear in the base storey. If the support displacement is not zero, a new value for the natural frequency is assumed and the procedure is repeated till we get the value of the base displacement as zero.

Example 11.2

Find the natural frequency of the three storeyed shear building as shown in Fig. 11.3, given the mass and the stiffness.

Solution

Since ŌĆśMATHEMATICAŌĆÖ can solve the problem using symbolic processing, we will write in symbolic form. Transfer matrix at top frame ( assume p = Žē n2 )

One we get natural frequency one can also get the mode shape as shown below. Holzer method for torsional vibrations

The Holzer method falls under the determinant search technique. The method can be applied to rectilinear or angular motions for damped as well as undamped systems. The method is best suited for systems where the components are arranged along the basic axis. Let us consider the torsional vibration for shafts as shown in Fig. 11.4.

The torsional moment equilibrium at node 1 can be written as Approximate methods for finding the natural frequencies

1  RayleighŌĆÖs quotient

In many practical situations involving MDOF systems, only the accurate estimation of the fundamental frequency is required. In such cases, laborious calculations to extract all the normal vibration modes of the system are not warranted and the approximate methods are desirable. This section discusses two approximate methods for estimating the fundamental frequency of MDOF systems.

The first method, RayleighŌĆÖs method, is an upper bound method based on energy principles and stiffness approach. The second method, DunkerleyŌĆÖs approximation, is based on the flexibility of the system eigenvalue problem and therefore provides lower bound estimation of the fundamental frequency. Thus the upper bound estimation of the fundamental frequency provided by RayleighŌĆÖs method can be complemented by the lower bound estimation afforded by DunkerleyŌĆÖs approximations to envelope true fundamental frequency.

Consider that an undamped single-degree-of-freedom (SDOF) mass spring system is in free harmonic motion given by 2  RayleighŌĆÖs quotient method to MDOF

In the above we have discussed RayleighŌĆÖs method to determine the fundamental frequency of an SDOF system. Application of the RayleighŌĆÖs method to determine the fundamental frequency of an MDOF system is presented in this section.

Consider the eigenvalue problem of the MDOF system represented by equation In Eq. 11.70b, the denominator is related to the kinetic energy for the rth mode and the numerator is related to the potential energy, or the strain energy of the rth mode. If the modal vector {Žå}r is replaced with any arbitrary vector {A}, Eq. 11.70b is written as where R({A}) is a scalar quantity referred to as RayleighŌĆÖs quotient. It is evolved from Eq. 7.51 that RayleighŌĆÖs quotient is dependent upon the known matrix [m] and [k] and the unknown arbitrary vector {A}. Obviously if {A}

coincides with one of the systems normal modes then ╬╗R is the corresponding eigenvalue or normal frequency of the system. A very important property of RayleighŌĆÖs quotient is and it also follows that for any vector {A}, if [K] is positive definite, R({A}) > 0. If [K] is positive semi-definite, Equation 11.65 thus indicates that the RayleighŌĆÖs quotient is never lower than the fundamental eigenvalue, and furthermore the minimum value the RayleighŌĆÖs quotient can assume is that of the fundamental eigenvalue itself. Therefore, RayleighŌĆÖs quotient is very good technique to estimate the fundamental frequency of MDOF systems. A reasonable estimate for the vector {A} corresponding to the fundamental mode is the vector of static displacement resulting from subjecting the masses in the system to forces proportional to their weights. Many seismic design code present expressions to estimate the fundamental frequency of high-rise building based on this concept. The natural frequency thus obtained is called the Rayleigh frequency ŽēR expressed as The accuracy of the Rayleigh frequency ŽēR depends entirely on the displacement vector {A} used to represent the vibration mode shape. In principle, any vector {A} may be selected which satisfies the geometric boundary conditions. However, any vector other than the true modal vector requires the action of additional external constraints to maintain equilibrium, which would in turn stiffen the structure, resulting in increased computed frequency. Therefore the true vibration mode will yield the lowest frequency

obtained by RayleighŌĆÖs method. Hence the approximation yielding the lowest frequency for a particular case is the best result.

Steps

1.     Estimate the fundamental mode of vibration. This may be done either by assuming the displacement of the modes directly or computing the displacement from the associated forces. If the frequency is computed for several displacement-assumed configurations, the smallest of the computed values will be close to the exact value of Žēn and the associated configuration is closest to the actual configuration.

Example 11.6

Determine the fundamental frequency of the shear frame shown in Fig. 11.3 by the improved Rayleigh method.

Solution

R00 method

Assume the mode shape as shown in Fig. 11.9. Maximum potential energy: R01 method

By calculating the shear in each floor let us improve the mode shape. The shear in each floor can be calculated as (see Fig. 11.10) Equating this improved potential energy to previously calculated kinetic energy

0.009 06 Žē n4  = 2.25Žē n2

Or

Žē n2  = 234.375; Žē n  = 15.309

R11 method

Improve the kinetic energy. Velocity at each storey level:

v1  = 0.0025 Žē n3 ; v 2  = 0.004 58 Žē n3 ; v3  = 0.006 25 Žē n3

Tmax(improved)  = Žē n6 /2[2 ├- 0.00252  + 1.5 ├- 0.004 582  + 1 ├- 0.006 252 ]

= 4.15135 ├- 10 Ōł'5 Žē n6

Equating maximum kinetic energy to maximum potential energy, we get

0.009 06 Žē n4  = 4.15135 ├- 10 Ōł'5 Žē n6 DunkerleyŌĆÖs approximation

It is another approximate method for estimating the fundamental frequency for MDOF systems. The method yields accurate results for systems for which

damping is negligible and the natural frequencies are well separated. DunkerleyŌĆÖs equation provides a ŌĆślower boundŌĆÖ estimate with fundamental frequency and is therefore complementary with the Rayleigh method that provides an ŌĆśupper boundŌĆÖ estimate with fundamental frequency.

To derive DunkerleyŌĆÖs equation consider the equation where [a] is the flexibility matrix.

The frequency equation is obtained by expanding the determinant of the characteristic matrix in Eq. 11.101.

Let us consider a two-degrees-of-freedom system with lumped mass diagonal matrix. Thus the resulting characteristic determinant becomes Expanding Eq. 11.102 results in the system frequency equation, i.e. second order equation in ╬╗ = 1/Žē2 given by,

If the roots are 1/ Žē12 , 1/ Žē 22

The relationship represented by Eq. 11.104 also holds true for systems having ŌĆśnŌĆÖ degrees of freedom. Extending this ŌĆśnŌĆÖ degrees of freedom system DunkerleyŌĆÖs approximation to the fundamental frequency is made on the assumption that if the fundamental frequency ŌĆśŽēŌĆÖ is much lower than the higher harmonics (Žē2,ŌĆ” Žēn) then the terms on the left hand side 1/ Žē 22 ŌĆ” 1/ Žē n2 can be calculated. The elimination of these terms yields an estimate of 1/ Žē12 which is higher than the true value thereby making the estimate of ŌĆśŽē1ŌĆÖ lower than the exact value of fundamental frequency. The DunkerleyŌĆÖs lower bound estimate of ŌĆśŽē1ŌĆÖ is approximated to In Eq. 11.106 the term ŌĆśaiimiŌĆÖ represents the contribution of each mass to 1/Žē12 in the absence of all other masses. Thus where Žē ii2 is the natural frequency of an SDOF system with mass ŌĆśmiŌĆÖ acting alone at state i. Hence DunkerleyŌĆÖs equation is given by Summary

In this chapter the sweeping technique combined with power method and transfer matrix methods have been discussed to find the natural frequencies of n-degrees-of-freedom system. In addition, RayleighŌĆÖs coefficient method and DunkerleyŌĆÖs approximate methods are also discussed to find the approximate fundamental frequency of an n-degrees-of-freedom system.

Study Material, Lecturing Notes, Assignment, Reference, Wiki description explanation, brief detail
Civil : Structural dynamics of earthquake engineering : Numerical solution methods for natural frequencies and mode shapes |