Problem Formulation

1.1 Dynamic Analysis

Consider the equations of motion for a linear system subjected to dynamic forces

Mr + Cr + Kr = R

where M is the mass matrix, C is the damping matrix, K is the stiffness matrix, r is the unknown displacement vector, and R is the load vector.

Considering mode superposition, we use the following transformation from the nodal displacements to the generalized displacements



where P is the number of mode shapes considered (in general P<<m, where m is the number of degrees of freedom), Z is a vector of generalized displacements, and Ф is the matrix of eigenvectors (mode shapes). The eigenvectors Фк and eigenvalues Xk = юк2 (юк are the circular frequencies) are obtained by solving the eigenproblem

In the presentation that follows we assume damping such that classical modal analysis can be used. Substituting Eq. (2) into Eq. (1) and pre-multiplying the resulting equations by Фт, we obtain

In these equations the right hand side vector in normalized coordinates is P = ®TR, and the mass matrix is an identity matrix І=ФТМ Ф. The damping matrix and the stiffness matrix are diagonal low – order matrices

л=фтс ф a2 = фтк ф (5)

The elements on the diagonals of these matrices are given by 2ak^k, <х>1, respectively, k = 1, ..p, ifk being the damping ratios. Thus, Eq. (4) consists of the p uncoupled equations

In many problems (e. g. earthquake loading) the load vector R, and therefore the right hand side terms


are given as discrete values at each time step, and not by analytical functions.

In summary, computation of the dynamic response by modal analysis involves the following steps.

a. Determine the matrices K, M, and C.

b. Determine the p requested eigenpairs Xk, ®k by solving the eigenproblem of Eq. (3).

c. Compute the modal coordinates Zk by solving Eqs. (6).

d. Compute the nodal displacements r by Eq. (2).

e. Calculate the element forces using the element stiffness properties.

1.2 Displacement Derivatives

The derivative expressions of the displacement vector r with respect to a design variable X, д r/дX, are given by differentiating Eq. (2)

The derivatives d Ф/Jd Xj can be evaluated efficiently by finite-differences using the CA approach, as will be shown later. Assuming that the damping ratios dfk are independent of the design variables

(which is typical, for example, in civil engineering structures), we calculate d ZfidXj by differentiation of Eq. (6)


qk = dZk / dX, qk = dZk / dX, qk = dZk / dX,

and substituting Eqs. (7), (10) into Eq. (9) yields

Note that the left hand sides of Eqs. (6) and (11) are similar, whereas the right hand sides are different. For identical initial conditions (e. g. qk = qk = 0 for t = 0) this similarity can be used to reduce the number of differential equations that must be solved during the solution process.

In summary, given the eigenpairs and the response for a certain design and time, evaluation of the displacement derivatives involves the following steps.

a. Evaluate the derivatives of the eigenpairs (д ФjjdXj and d XjjdXj ).

b. Compute the right side of Eq. (11).

c. Compute the derivatives qk = d ZjjdXj by solving Eq. (11).

d. Evaluate the displacement derivatives д r/dXj by Eq. (8).

Assuming a problem with p mode shapes and n design variables, the main computational effort is involved in the following two steps:

a. Solution of the pn differential equations (11).

b. Evaluation of pn derivatives of the eigenpairs ( д Ф jJ d Xj and д ak/ d Xj).

A procedure intended to reduce the number of differential equations to be solved during the solution process is proposed below. Efficient evaluation of the derivatives of the eigenpairs, using finite – difference and the CA approach, is presented later.

2. Reducing the Number of Differential Equations

Due to the linearity of Eq. (11), we can use superposition and divide it into the following 3 equations with identical initial conditions

and assuming that the load vector can be expressed in the form R(X, t)=R(X) g(t), then Eqs. (15), (16) describe similar functions in time with different amplitudes. For zero initial conditions (or, if we neglect the influence of the homogeneous solution), the ratio between the two displacement functions of Eqs. (6) and (12) is equal to the ratio between the right-hand side terms. Thus, given the solutions Zk of Eq. (6) for all p modes, the solutions q(1) of Eq. (12) for i = 1 can be determined directly by

To find q1’2’1,q(3), Eq. (12) must be solved for i=2 and i=3. ForX we have to solve the two equations

qf + 2aktkq (2) + Co q ™ = -2^Z k

q? + ^qk3) + – kqk3) = ~gk

Given the solutions of Eqs. (18), (19) with respect to X1, it is observed that the solutions for any other variable Xj can be determined directly by



qk2,( X)=j^q?( X)


q(3)(x -) =|X – qf( xf)



In the particular case where Ok and R are orthogonal we obtain Pk=OkTR = 0. From Eq. (6) we have Zk = Zk = Zk = 0 , and from Eqs. (17) – (19) we find q® = qf’1 = qf’1 = 0 .

In summary, assuming a problem with p considered mode shapes and n design variables, the number of times that the differential equations (11) must be solved in order to perform sensitivity analysis is usually pn. Considering the procedure presented in this section and assuming that the solution of the analysis problem [Eq. (6)] is known, the number of times that the differential equations must be solved in order to perform sensitivity analysis is only 2p [Eqs. (18), (19)]. Thus, the ratio between the two numbers is pn/2p=n/2, which means a significant reduction in the computational cost. For example, for a problem with 10 design variables, the procedure presented requires about 20% of the effort involved in complete sensitivity analysis.