广义特征值问题GEVP高层建筑结构动力学算例和Matlab程序.doc

资料内容仅供您学习参考，如有不当或者侵权，请联系改正或者删除。

m4
x4(t)
k4 / 2
k3 / 2
k4 / 2
m3
x3(t)
k3 / 2
m2
k / 2
k2 / 2
x2(t)
x1(t)
k2 / 2
k1 / 2
m1
k1 / 2
Answer:
a)
Determine the governing equations for the system.
From the analysis of free-body diagrams, we write
m1&x&1 + (k1 + k2)x1 - k2x2 = 0
m2&x&2 - k2x1 + (k2 + k3)x2 - k3x3 = 0
m3&x&3 - k3x2 + (k3 + k4)x3 - k4x4 = 0
,
m4&x&4 - k4x3 + k4x4 = 0
which is, in matrix form,
ém1
0
0
0
0 ù
ék1 + k2
- k2
k2 + k3
- k3
0
0 ù
ê
ú
0
ê
ú
0
m2
- k2
- k3
k3 + k - k ú
- k4
0
ê
ú
ê
ú
&x&+
x = 0
(1)
ê
ê
ë
ú
ú
û
ê
ê
ë
0
0
0
m3
0
0
0
0
4
4
ú
û
0
m4
0
k4
1

b)
Determine the modes of vibration. Assume mi = 4000 kg and ki = 5000 N/m, for i = 1, . . . , 4.
For the assumed values of mi and ki, we have, from Eq.(1),
é4000
0
4000
0
0
0
0 ù
é10000 - 5000
0
0 ù
ê
ú
ê
ú
0
0
- 5000 10000 - 5000
0
ê
ú
ê
ê
ê
ë
ú
&x&+
x = 0
(2)
(3)
ê
ê
ë
ú
ú
û
- 5000 10000 - 5000
ú
ú
û
0
0
4000
0
0
0
0
0
4000
0
- 5000
5000
Simplified, we obtain
é4 0 0 0ù
é10 - 5
0
0 ù
ê
ê
ê
ê
ë
ú
ú
ú
ú
û
ê
ú
0 4 0 0
0 0 4 0
0 0 0 4
- 5 10 - 5
0
ê
ê
ê
ë
ú
&x&+
x = 0
- 5 10 - 5ú
ú
0
0
0
- 5
5
û
Solving the above GEVP with a MATLAB script program (Appendix A), we obtain the 4 modes of
the 4-story building
Λ = diag(0.1508,1.2500, 2.9341,4.4151) , ω = diag(0.3883,1.1180,1.7129,2.1012)
(4)
(5)
é - 0.1140 0.2887 0.3283 - 0.2143ù
ê
ú
ú
ú
ú
û
- 0.2143 0.2887 - 0.1140 0.3283
- 0.2887 0.0000 - 0.2887 - 0.2887
- 0.3283 - 0.2887 0.2143 0.1140
U = [u1 u2 u3 u4] = êê
ê
ë
c)
Plot of the mode shapes.
The mode shapes can be plotted by the MATLAB script program (Appendix A) as well
2

d)
Determine the free response of the system due to the initial displacement.
Under arbitrary initial conditions, the response of the system can be expressed as
4
å
x = ciui sin(wit +y i )
i=1
Due to the initial displacement, we have
x = c1u1 sin(w1t + p ) + c2u2 sin(w2t + p ) + c3u3 sin(w3t + p ) + c4u4 sin(w4t + p )
(6)
(7)
2
2
2
2
Solving a set of linear simultaneous equations
[u1 u2 u3 u4](c1 c2 c3 c4)T = x0
we have
(c1 c2 c3 c4)T = (- 0.0414 0.0508 0.0130 - 0.0063)
Thus, the free response of the system due to the intial displacement is as expressed inn Eq.(6), where
wi and ci,(i = 1,2,3,4) are given by Eqs.(4) and (7) respectively.
e)
By definition, a node of a mode shape is simply the location of zero displacement. Determine the
node locations for modes 2 to 4. Under what conditions may the node locations be important?
Please give at least 2 examples.
Linear interpolation of the curves of modes 2 to 4 gives the three nodes locations
mode
2
3
4
Node location ( ith story)
3.0
1.7423, 3.5740
1.3949, 2.5321, 3.7169
If we want to excite the structure to measure the vibration and do some modal analyses, it’s important to
select the proper nodes for excitation. If the point for excitation happens to be at or near a node of a mode,
the specific mode will not be excited and can not be measured.
A second example to illustrate the importance of the locations of nodes is related to joint two large machine
parts or components together. If the points to link (assembly) the two elements together coincide with the
nodes of the principle modes, the bolting force will vary little due to vibration of the structures.
3

Answer:
a)
Determine the third normalized eigenvector of the system
Assume the third normalized eigenvector of the system to be u3 = (u31,u32,u33)T , then, from the
orthogonality of the eigenvectors, we have
u
u
u
T
1
Mu3 = 0
Mu3 = 0
Mu3 = 1
T
2
T
3
Solving the above equations (Appendix B), we have
u3 = ( 0.6591 - 0.4950 0.1588)T
b)
If the stiffness matrix of the system is known, how would you determine the natural frequencies
of the system without having to solve a generalized eigenvalue problem? What is the benefit of
determining the natural frequencies using this particular approach?
From equation
Ku = w Mu
2
we know
wi = (Kui )1 ,i = 1,2,L
(Mui )1
By this approach, it avoids the time-consuming computation to solve a generalized eigenvalue problem
and can get the result faster and more economically.
3. In class, we considered the following symmetric generalized eigenvalue problem (GEVP)
[K]u = l [M]u
and proved the orthogonality conditions of the mode shapes with respect to the mass and stiffness
matrices. Often it is more desirable to solve a standard eigenvalue problem (SEVP) of the form
[A]u = l u
where [A] is generally nonsymmetric. The rationale for solving the SEVP is two-fold. First, many
4

numerically efficient algorithms have been developed over the years specifically for the solution of a
SEVP. Second, the behavior of a SEVP is much better understood, and many theorems have been
formulated for a SEVP as opposed to a GEVP. Prove the orthogonality conditions of the mode shapes
for the SEVP. [Hint: you may want to consider the adjoint eigenvalue problem.]
Prove:
If [A] is symmetric, we can prove the orthogonality conditions of the mode shapes for the SEVP.
However, if [A] is nonsymmetric, we should prove the biorthogonality conditions instead. Actually the
orthogonality condition can be derived easily if the biorthogonality conditions are proved.
a)
biorthogonality conditions
For the eigenvalue problem
Ax j = lx j
(3-1)
For the adjoint eigenvalue problem
A yi = li yi
T
(3-2)
premultiplying
yi
T
on both sides of Eqs.(3-1) gives
yi Ax j = lyi x j
T T
(3-3)
Transpose of Eq.(3-2) and postmultiplying x j on both sides gives
yi Ax j = li yi x j
T T
(3-4)
Subtracting Eq.(3-4) from Eq.(3-3) gives
(li - l j )yi
If i ¹ j , li ¹ l j , then
yi
x j = 0 (3-5)
T
x j = 0
T
which is the biorthogonality conditions.
b)
orthogonality conditions
When [A] is symmetry, then A = A and yi = xi , where Axi = lixi . Then Eq.(3-5) gives
T
xi
x j = 0 (3-6)
T
which is the orthogonality conditions.
Appendix A
%
%
[V,D] = EIG(A,B) produces a diagonal matrix D of generalized
eigenvalues and a full matrix V whose columns are the
5

%
corresponding eigenvectors so that A*V = B*V*D.
K=[10 -5 0 0; -5 10 -5 0; 0 -5 10 -5; 0 0 -5 5 ]
M=4.*eye(4)
[V,D]=eig(K,M)
%
clf;
plot(V(:,:),[1 2 3 4]','o-');
legend('Mode1','Mode2','Mode3','Mode4');
K*V
M*V*D
%
x0=[0.025 0.020 0.01 0.001]';
C=inv(V)*x0
Appendix B
function hw1new_2
u30=[1, 1, 1]';
[u3,fval,exitflg]=fminsearch(@fu3,u30,[]);
u3,fval
u3=u3'
function f=fu3(x)
u1=[0.2754946 0.3994672 0.4490562]';
u2=[0.6916979 0.2974301 -0.3389320]';
m=[1 0 0; 0 2 0; 0 0 3];
f=abs(u1'*m*x);
f=f+abs(u2'*m*x);
f=f+abs(x'*m*x-1);
6