平面框架结构的数值分析与有限元分析实例
姓名: 张建
学号: 11011365
联系方式: 159********
指导教师: 王骑
2011年12月3日
题目:平面框架结构的数值分析与有限元分析
如下图所示平面框架结构,所受荷载如图所示,试分别用MATLAB 数值分析软件和ANSYS有限元分析软件分析该结构各杆件的位移和内力。结构中各个截面的参数都为：E=210GPa,I=32×10?5m4，A=1.2×10?2m2。在ANSYS平台上，平面梁单元采用Beam3单元，不考虑剪切影响。
分析：
（1）第一种方法：采用MATLAB进行数值分析
一、步骤
1、对结点和单元进行编号，选定整体坐标系和局部坐标系，如下图所示。
2、求出各单元在整体坐标系中的单元刚度矩阵，通过调用MATLAB的PlaneFrameElementStiffness（E，A，I，L，theta）函数实现。
3、将各单刚子块对号入座，形成结构原始刚度矩阵，通过调用MATLAB 的PlaneFrameAssemble（K，k，i，j）函数实现。
4、计算非结点荷载作用下的各单元固端力、等效结点荷载及综合结点荷载，根据（《结构力学》上册-李廉锟主编）表10-3等直杆单元的固端力查表，并在MATLAB中编程计算所需数据。
5、引入支承条件，修改结构的原始刚度方程。结点1和10为固定端，所以结点1和10处的线位移和角位移均为0，在原始刚度矩阵中删去与零位移对应的行和列，同时在结点位移列向量和结点外力列向量中删去相应的行，在MATLAB中通过命令k=K（4:27，4:27）便可得到修改后的结构的原始刚度矩阵。
6、解结构的刚度方程，求出结点位移，在MATLAB中通过命令u=k\f 来求得所需数据。
7、计算各单元杆端力，在MATLAB中通过命令PlaneFrameElement
Forces(E,A,I,L,theta,u)+各单元在局部坐标系下的固端力来求解。
二、计算结果数据整理
MATLAB计算结构结点位移
MATLAB计算结构单元杆端力
三、MATLAB编程代码如下
addpath C:\temp\toolbox\svm
which PlaneFrameElementStiffness.m
E=210E6
I=32E-5
A=1.2E-2
L1=2
L2=2
L3=1
L4=3
L5=3*sqrt(3)
L6=1
L7=3
L8=3
L9=2
L10=2
k1=PlaneFrameElementStiffness(E,A,I,L1,90)
k2=PlaneFrameElementStiffness(E,A,I,L2,90)
k3=PlaneFrameElementStiffness(E,A,I,L3,90)
k4=PlaneFrameElementStiffness(E,A,I,L4,60)
k5=PlaneFrameElementStiffness(E,A,I,L5,-30)
k6=PlaneFrameElementStiffness(E,A,I,L6,-90)
k7=PlaneFrameElementStiffness(E,A,I,L7,0)
k8=PlaneFrameElementStiffness(E,A,I,L8,0)
k9=PlaneFrameElementStiffness(E,A,I,L9,-90)
k10=PlaneFrameElementStiffness(E,A,I,L10,-90) %整体坐标系下的单元刚度矩阵
K=zeros(30,30)
K=PlaneFrameAssemble(K,k1,1,2)
K=PlaneFrameAssemble(K,k2,2,3)
K=PlaneFrameAssemble(K,k3,3,4)
K=PlaneFrameAssemble(K,k4,4,5)
K=PlaneFrameAssemble(K,k5,5,6)
K=PlaneFrameAssemble(K,k6,6,8)
K=PlaneFrameAssemble(K,k7,3,7)
K=PlaneFrameAssemble(K,k8,7,8)
K=PlaneFrameAssemble(K,k9,8,9)
K=PlaneFrameAssemble(K,k10,9,10) %集成原始刚度矩阵
F1=0
F10=0
F2n2=0 %2单元2结点处的轴力，其他同理
a=2,b=0,l=2,q=10
F2s2=-q*a*(2*l^3-2*l*a^2+a^3)/2/(l^3) %2单元2结点处的剪力，其他同理
M2f2=-q*a^2*(6*l^2-8*l*a+3*a^2)/12/(l^2) %2单元2结点处的弯矩，其他同理
F3n2=0
F3s2=-q*a^3*(2*l-a)/2/(l^3)
M3f2=q*a^3*(4*l-3*a)/12/(l^2)
F2=[F2n2;F2s2;M2f2;F3n2;F3s2;M3f2]
a=1,b=0,l=1,q=10
F3n3=0,F3s3=-q*a*(2*l^3-2*l*a^2+a^3)/2/(l^3),M3f3=-q*a^2*(6*l^2-8*l*a+3*a^2)/12/(l^2) F4n3=0
F4s3=-q*a^3*(2*l-a)/2/(l^3)
M4f3=q*a^3*(4*l-3*a)/12/(l^2)
F3=[F3n3;F3s3;M3f3;F4n3;F4s3;M4f3]
a=3,b=0,l=3,q=5
F4n4=0,F4s4=-q*a*(2*l^3-2*l*a^2+a^3)/2/(l^3),M4f4=-q*a^2*(6*l^2-8*l*a+3*a^2)/12/(l^2)
F5n4=0,F5s4=-q*a^3*(2*l-a)/2/(l^3),M5f4=q*a^3*(4*l-3*a)/12/(l^2)
F4=[F4n4;F4s4;M4f4;F5n4;F5s4;M5f4]
a=3*sqrt(3),b=0,l=3*sqrt(3),q=15
F5n5=0,F5s5=q*a*(2*l^3-2*l*a^2+a^3)/2/(l^3),M5f5=q*a^2*(6*l^2-8*l*a+3*a^2)/12/(l^2)
F6n5=0,F6s5=q*a^3*(2*l-a)/2/(l^3),M6f5=-q*a^3*(4*l-3*a)/12/(l^2)
F5=[F5n5;F5s5;M5f5;F6n5;F6s5;M6f5] %F5n5为5单元5结点处的轴力，其他同理
a=1,b=0,l=1,q=30
F6n6=0,F6s6=q*a*(2*l^3-2*l*a^2+a^3)/2/(l^3),M6f6=q*a^2*(6*l^2-8*l*a+3*a^2)/12/(l^2)
F8n6=0,F8s6=q*a^3*(2*l-a)/2/(l^3),M8f6=-q*a^3*(4*l-3*a)/12/(l^2)
F6=[F6n6;F6s6;M6f6;F8n6;F8s6;M8f6] %6单元在局部坐标系中的固端力，其他同理
a=2,b=0,l=2,q=30
F8n9=0,F8s9=q*a*(2*l^3-2*l*a^2+a^3)/2/(l^3),M8f9=q*a^2*(6*l^2-8*l*a+3*a^2)/12/(l^2)
F9n9=0,F9s9=q*a^3*(2*l-a)/2/(l^3),M9f9=-q*a^3*(4*l-3*a)/12/(l^2)
F9=[F8n9;F8s9;M8f9;F9n9;F9s9;M9f9]
a=3,b=0,l=3,F=50
F3n7=0,F3s7=0,M3f7=0
F7n7=0,F7s7=0,M7f7=0
F7=[F3n7;F3s7;M3f7;F7n7;F7s7;M7f7]
a=0,b=3,l=3,F=50
F7n8=0,F7s8=0,M7f8=0
F8n8=0,F8s8=0,M8f8=0
F8=[F7n8;F7s8;M7f8;F8n8;F8s8;M8f8]
alpha1=pi/2,alpha2=pi/2,alpha3=pi/2,alpha4=pi/3,alpha5=-pi/6,alpha6=-pi/2,alpha7=0,alpha8=0,alpha9=-pi/2, alpha10=-pi/2
T1=[cos(alpha1) sin(alpha1) 0 0 0 0;-sin(alpha1) cos(alpha1) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha1) sin(alpha1) 0;
0 0 0 -sin(alpha1) cos(alpha1) 0;0 0 0 0 0 1]
T2=[cos(alpha2) sin(alpha2) 0 0 0 0;-sin(alpha2) cos(alpha2) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha2) sin(alpha2) 0;
0 0 0 -sin(alpha2) cos(alpha2) 0;0 0 0 0 0 1]
T3=[cos(alpha3) sin(alpha3) 0 0 0 0;-sin(alpha3) cos(alpha3) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha3) sin(alpha3) 0;
0 0 0 -sin(alpha3) cos(alpha3) 0;0 0 0 0 0 1]
T4=[cos(alpha4) sin(alpha4) 0 0 0 0;-sin(alpha4) cos(alpha4) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha4) sin(alpha4) 0;
0 0 0 -sin(alpha4) cos(alpha4) 0;0 0 0 0 0 1]
T5=[cos(alpha5) sin(alpha5) 0 0 0 0;-sin(alpha5) cos(alpha5) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha5) sin(alpha5) 0;
0 0 0 -sin(alpha5) cos(alpha5) 0;0 0 0 0 0 1]
T6=[cos(alpha6) sin(alpha6) 0 0 0 0;-sin(alpha6) cos(alpha6) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha6) sin(alpha6) 0;
0 0 0 -sin(alpha6) cos(alpha6) 0;0 0 0 0 0 1]
T7=[cos(alpha7) sin(alpha7) 0 0 0 0;-sin(alpha7) cos(alpha7) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha7) sin(alpha7) 0;
0 0 0 -sin(alpha7) cos(alpha7) 0;0 0 0 0 0 1]
T8=[cos(alpha8) sin(alpha8) 0 0 0 0;-sin(alpha8) cos(alpha8) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha8) sin(alpha8) 0;
0 0 0 -sin(alpha8) cos(alpha8) 0;0 0 0 0 0 1]
T9=[cos(alpha9) sin(alpha9) 0 0 0 0;-sin(alpha9) cos(alpha9) 0 0 0 0;0 0 1 0 0 0;0 0 0 cos(alpha9) sin(alpha9) 0;
0 0 0 -sin(alpha9) cos(alpha9) 0;0 0 0 0 0 1]
T10=[cos(alpha10) sin(alpha10) 0 0 0 0;-sin(alpha10) cos(alpha10) 0 0 0 0;0 0 1 0 0 0;
0 0 0 cos(alpha10) sin(alpha10) 0;0 0 0 -sin(alpha10) cos(alpha10) 0;0 0 0 0 0 1]
T1z=T1'
T2z=T2'
T3z=T3'
T4z=T4'
T5z=T5'
T6z=T6'
T7z=T7'
T8z=T8'
T9z=T9'
T10z=T10'
Fg1=0 %Fg1为1单元在整体坐标系下固端力
Fg10=0
Fg2=T2z*F2 %Fg2为2单元在整体坐标系下固端力,其他同理,T2z为T2的转置矩阵,T为坐标转换矩阵Fg3=T3z*F3
Fg4=T4z*F4
Fg5=T5z*F5
Fg6=T6z*F6
Fg7=T7z*F7
Fg8=T8z*F8
Fg9=T9z*F9
F2g1=[0;0;0]
F2g2=Fg2(1:3)
F3g2=Fg2(4:6)
F3g3=Fg3(1:3) %整体坐标系下3单元在结点3处的固端力
F4g3=Fg3(4:6)
F4g4=Fg4(1:3)
F5g4=Fg4(4:6)
F5g5=Fg5(1:3)
F6g5=Fg5(4:6)
F6g6=Fg6(1:3)
F8g6=Fg6(4:6)
F3g7=[0;0;0]
F7g7=[0;0;0]
F7g8=[0;0;0]
F8g8=[0;0;0]
F8g9=Fg9(1:3)
F9g9=Fg9(4:6)F9g10=[0;0;0]
Fe2=-F2g1-F2g2
Fe3=-F3g2-F3g3-F3g7 %结点3处的等效结点荷载
Fe4=-F4g3-F4g4,Fe5=-F5g4-F5g5,Fe6=-F6g5-F6g6,Fe7=-F7g7-F7g8,Fe8=-F8g6-F8g8-F8g9,Fe9=-F9g9-F9g10 Fz2=[0;0;0]+Fe2 %综合结点荷载
Fz3=[0;0;0]+Fe3,Fz4=[0;0;0]+Fe4,Fz5=[0;0;0]+Fe5,Fz6=[0;0;0]+Fe6,Fz7=[0;-50;0]+Fe7,Fz8=[0;0;0]+Fe8,
Fz9=[0;0;0]+Fe9
F=[Fz2;Fz3;Fz4;Fz5;Fz6;Fz7;Fz8;Fz9]
k=K(4:27,4:27) %引入支承条件，提取刚度矩阵
u=k\F %反斜线符号“\”用于高斯消元法
U=[0;0;0;u;0;0;0]
F=K*U
u1=[U(1);U(2);U(3);U(4);U(5);U(6)] %1单元结点位移矢量，其他同理
u2=[U(4);U(5);U(6);U(7);U(8);U(9)]
u3=[U(7);U(8);U(9);U(10);U(11);U(12)]
u4=[U(10);U(11);U(12);U(13);U(14);U(15)]
u5=[U(13);U(14);U(15);U(16);U(17);U(18)]
u6=[U(16);U(17);U(18);U(22);U(23);U(24)]
u7=[U(7);U(8);U(9);U(19);U(20);U(21)]
u8=[U(19);U(20);U(21);U(22);U(23);U(24)]
u9=[U(22);U(23);U(24);U(25);U(26);U(27)]
u10=[U(25);U(26);U(27);U(28);U(29);U(30)]
f1=PlaneFrameElementForces(E,A,I,L1,90,u1) %1单元结点力矢量，其他同理
f2=PlaneFrameElementForces(E,A,I,L2,90,u2),f3=PlaneFrameElementForces(E,A,I,L3,90,u3), f4=PlaneFrameElementForces(E,A,I,L4,60,u4),
f5=PlaneFrameElementForces(E,A,I,L5,-30,u5),f6=PlaneFrameElementForces(E,A,I,L6,-90,u6), f7=PlaneFrameElementForces(E,A,I,L7,0,u7),
f8=PlaneFrameElementForces(E,A,I,L8,0,u8),f9=PlaneFrameElementForces(E,A,I,L9,-90,u9), f10=PlaneFrameElementForces(E,A,I,L10,-90,u10)
Fgd1=PlaneFrameElementForces(E,A,I,L1,90,u1)+F1
Fgd2=PlaneFrameElementForces(E,A,I,L2,90,u2)+F2
Fgd3=PlaneFrameElementForces(E,A,I,L3,90,u3)+F3
Fgd4=PlaneFrameElementForces(E,A,I,L4,60,u4)+F4
Fgd5=PlaneFrameElementForces(E,A,I,L5,-30,u5)+F5
Fgd6=PlaneFrameElementForces(E,A,I,L6,-90,u6)+F6
Fgd7=PlaneFrameElementForces(E,A,I,L7,0,u7)+F7
Fgd8=PlaneFrameElementForces(E,A,I,L8,0,u8)+F8
Fgd9=PlaneFrameElementForces(E,A,I,L9,-90,u9)+F9
Fgd10=PlaneFrameElementForces(E,A,I,L10,-90,u10)+F10
备注：所调用的几个M文件源代码如下：
function y = PlaneFrameElementStiffness(E,A,I,L,theta)
%PlaneFrameElementStiffness This function returns the element
% stiffness matrix for a plane frame
% element with modulus of elasticity E,
% cross-sectional area A, moment of
% inertia I, length L, and angle
% theta (in degrees).
% The size of the element stiffness
% matrix is 6 x 6.
x = theta*pi/180;
C = cos(x);
S = sin(x);
w1 = A*C*C + 12*I*S*S/(L*L);
w2 = A*S*S + 12*I*C*C/(L*L);
w3 = (A-12*I/(L*L))*C*S;
w4 = 6*I*S/L;
w5 = 6*I*C/L;
y = E/L*[w1 w3 -w4 -w1 -w3 -w4 ; w3 w2 w5 -w3 -w2 w5 ;
-w4 w5 4*I w4 -w5 2*I ; -w1 -w3 w4 w1 w3 w4 ;
-w3 -w2 -w5 w3 w2 -w5 ; -w4 w5 2*I w4 -w5 4*I];
function y = PlaneFrameAssemble(K,k,i,j)
%PlaneFrameAssemble This function assembles the element stiffness % matrix k of the plane frame element with nodes % i and j into the global stiffness matrix K.
% This function returns the global stiffness
% matrix K after the element stiffness matrix
% k is assembled.
K(3*i-2,3*i-2) = K(3*i-2,3*i-2) + k(1,1);
K(3*i-2,3*i-1) = K(3*i-2,3*i-1) + k(1,2);
K(3*i-2,3*i) = K(3*i-2,3*i) + k(1,3);
K(3*i-2,3*j-2) = K(3*i-2,3*j-2) + k(1,4);
K(3*i-2,3*j-1) = K(3*i-2,3*j-1) + k(1,5);
K(3*i-2,3*j) = K(3*i-2,3*j) + k(1,6);
K(3*i-1,3*i-2) = K(3*i-1,3*i-2) + k(2,1);
K(3*i-1,3*i-1) = K(3*i-1,3*i-1) + k(2,2);
K(3*i-1,3*i) = K(3*i-1,3*i) + k(2,3);
K(3*i-1,3*j-2) = K(3*i-1,3*j-2) + k(2,4);
K(3*i-1,3*j-1) = K(3*i-1,3*j-1) + k(2,5);
K(3*i-1,3*j) = K(3*i-1,3*j) + k(2,6);
K(3*i,3*i-2) = K(3*i,3*i-2) + k(3,1);
K(3*i,3*i-1) = K(3*i,3*i-1) + k(3,2);
K(3*i,3*i) = K(3*i,3*i) + k(3,3);
K(3*i,3*j-2) = K(3*i,3*j-2) + k(3,4);
K(3*i,3*j-1) = K(3*i,3*j-1) + k(3,5);
K(3*i,3*j) = K(3*i,3*j) + k(3,6);
K(3*j-2,3*i-2) = K(3*j-2,3*i-2) + k(4,1);
K(3*j-2,3*i-1) = K(3*j-2,3*i-1) + k(4,2);
K(3*j-2,3*i) = K(3*j-2,3*i) + k(4,3);
K(3*j-2,3*j-2) = K(3*j-2,3*j-2) + k(4,4);
K(3*j-2,3*j-1) = K(3*j-2,3*j-1) + k(4,5);
K(3*j-2,3*j) = K(3*j-2,3*j) + k(4,6);
K(3*j-1,3*i-2) = K(3*j-1,3*i-2) + k(5,1);
K(3*j-1,3*i-1) = K(3*j-1,3*i-1) + k(5,2);
K(3*j-1,3*i) = K(3*j-1,3*i) + k(5,3);
K(3*j-1,3*j-2) = K(3*j-1,3*j-2) + k(5,4);
K(3*j-1,3*j-1) = K(3*j-1,3*j-1) + k(5,5);
K(3*j-1,3*j) = K(3*j-1,3*j) + k(5,6);
K(3*j,3*i-2) = K(3*j,3*i-2) + k(6,1);
K(3*j,3*i-1) = K(3*j,3*i-1) + k(6,2);
K(3*j,3*i) = K(3*j,3*i) + k(6,3);
K(3*j,3*j-2) = K(3*j,3*j-2) + k(6,4);
K(3*j,3*j-1) = K(3*j,3*j-1) + k(6,5);
K(3*j,3*j) = K(3*j,3*j) + k(6,6);
y = K;
function y = PlaneFrameElementForces(E,A,I,L,theta,u)
%PlaneFrameElementForces This function returns the element force
% vector given the modulus of elasticity E, % the cross-sectional area A, the moment of % inertia I, the length L, the angle theta
% (in degrees), and the element nodal
% displacement vector u.
x = theta * pi/180;
C = cos(x);
S = sin(x);
w1 = E*A/L;
w2 = 12*E*I/(L*L*L);
w3 = 6*E*I/(L*L);
w4 = 4*E*I/L;
w5 = 2*E*I/L;
kprime = [w1 0 0 -w1 0 0 ; 0 w2 w3 0 -w2 w3 ;
0 w3 w4 0 -w3 w5 ; -w1 0 0 w1 0 0 ;
0 -w2 -w3 0 w2 -w3 ; 0 w3 w5 0 -w3 w4];
T = [C S 0 0 0 0 ; -S C 0 0 0 0 ; 0 0 1 0 0 0 ;
0 0 0 C S 0 ; 0 0 0 -S C 0 ; 0 0 0 0 0 1];
y = kprime*T* u;
（2）第二种方法：采用ANSYS进行有限元分析
一、步骤（GUI操作方式）
1、设置计算类型
ANSYS Main Menu:Preferences…→Structural→OK
2、选择单元类型
ANSYS Main Menu: Preprocessor→Element Type→Add/Edit/Delete…
→Add…→beam：2D elastic 3→OK (返回到Element Types窗口)→Close 3、定义材料参数
ANSYS Main Menu: Preprocessor→Material Props→Material Models→Structural→Linear→Elastic→Isotropic: EX:3e11 (弹性模量)→OK→鼠标点击该窗口右上角的“×”来关闭该窗口
4、定义实常数以确定平面问题的厚度
ANSYS Main Menu: Preprocessor→Real Constants…→Add/Edit/ Delete→Add→Type 1 Beam3→OK→Real Constant Set No: 1 (第1号实常数), Cross-sectional area:1.2e-2 (梁的横截面积),Area moment of inertia:3.2e-4(惯性矩)→OK→Close
5、生成几何模型
生成节点
ANSYS Main Menu: Preprocessor→Modeling→Creat→Keypoints→In Active CS→Keypoint number 1→X:0,Y:0,Z:0→Apply→Keypoint number 2→X:0,Y:2,Z:0→Apply→Keypoint number 3→X:0,Y:4,Z:0→Apply→Keypoint number 4→X:0,Y:5,Z:0→Apply→Keypoint number 5→X:1.5,
Y:1.5*(3**0.5)+5,Z:0→Apply→Keypoint number 6→X:6,Y:5,Z:0→Apply
→Keypoint number 7→X:3,Y:4,Z:0→Apply→Keypoint number 8→X:6, Y:4,Z:0→Apply→Keypoint number 9→X:6,Y:2,Z:0→Apply→Keypoint number 10→X:6,Y:0,Z:0→OK
生成单元
ANSYS Main Menu: Preprocessor→Modeling→Create→Lines→Straight Line→选择节点1、2(生成单元1)→Apply→选择节点2、3(生成单元2)→Apply→选择节点3、4(生成单元3)→Apply→选择节点4、5(生成单元4)→Apply→选择节点5、6(生成单元5)→Apply→选择节点6、8(生成单元6)→Apply→选择节点3、7(生成单元7)→Apply→选择节点7、8(生成单元8)→Apply→选择节点8、9(生成单元9)→Apply→选择节点9、10(生成单元10)→OK
划分单元
ANSYS Main Menu: Preprocessor→Meshing→Mesh Tool…→Lines Set…→Pick All…→Element edge length:0.5→Mesh…→Pick All→OK
6、模型施加约束和外载
显示单元编号(划分单元后单元自动重新编号)
ANSYS Utility Menu:PlotCtrls→Numbering…→Elem/Attrib numbering：Element numbers→OK
施加均布荷载
ANSYS Main Menu: Solution→Define Loads→Apply→Structural→Pressure→On Beams→选取单元5、6、7、8、9、10→apply→VALI：
-10000→VALJ：-10000→apply→选取单元11、12、13、14、15、16→apply→VALI：-5000→VALJ：-5000→apply→选取单元17、18、19、20、21、22、23、24、25、26、27→apply→VALI：15000→VALJ：15000→apply→选取单元28、29、42、43、44、45→apply→VALI：30000→VALJ：30000→OK
显示节点编号
ANSYS Utility Menu:PlotCtrls→Numbering…→Elem/Attrib numbering：No numbers，Node numbers：On→OK
施加集中荷载
ANSYS Main Menu: Solution→Define Loads→Apply→Structural→Force/Moment→On Nodes→选择节点31→apply→Direction of force: FY →VALUE：-50000→OK
左、右下角节点加约束
ANSYS Main Menu: Solution→Define Loads→Apply→Structural→Displacement→On Nodes→选取节点1和节点46→Apply→Lab:ALL DOF→OK
7、分析计算
ANSYS Main Menu: Solution→Solve→Current LS→OK→Should The Solve Command be Executed? Y→Close (Solution is done! )→关闭文字
窗口
8、结果显示
ANSYS Main Menu: General Postproc→Plot Results→Deformed
Shape …→Def + Undeformed→OK (返回到Plot Results)
得变形图如下：
ANSYS Main Menu: General Postproc→Element Table→Define Table…→Add…→Lab：NI→Item，Comp：By sequence num→SMISC：SMISC1→Apply→Lab：NJ→Item，Comp：By sequence num→SMISC：SMISC7→Apply→Lab：QI→Item，Comp：By sequence num→SMISC：SMISC2→Apply→Lab：QJ→Item，Comp：By sequence num→SMISC：SMISC8→Apply→Lab：MI→Item，Comp：By sequence num→SMISC：SMISC6→Apply→Lab：MJ→Item，Comp：By sequence num→SMISC：SMISC12→OK→Close
ANSYS Main Menu: General Postproc→Plot Results→Contour Plot
→Line Elem Res…→LabI：NI，LabJ：NJ→OK(返回到Plot Results) 得轴力图如下：
ANSYS Main Menu: General Postproc→Plot Results→Contour Plot →Line Elem Res…→LabI：QI，LabJ：QJ→OK(返回到Plot Results) 得剪力图如下：
ANSYS Main Menu: General Postproc→Plot Results→Contour Plot →Line Elem Res…→LabI：MI，LabJ：MJ→Fact：-1→OK(返回到Plot Results)
得弯矩图如下：
9、数据提取
结点位移提取
ANSYS Main Menu: General Postproc→List Results→Nodal Solution …→Nodal Solution→DOF Solution→X-Component of displacement→Apply→提取结点X方向位移→关闭窗口→Y-Component of displacement→Apply→提取结点Y方向位移→关闭窗口→
Z-Component of rotation→OK→提取结点转角位移→关闭窗口单元杆端力提取
ANSYS Main Menu: General Postproc→List Results→Elem Table Data…→Lab1-9：NI，NJ，QI，QJ，MI，MJ →OK→提取单元杆端力→关闭窗口
二、计算结果数据整理
ANSYS计算结构结点位移
ANSYS计算结构单元杆端力
（3）MATLAB和ANYSY计算结果对比
通过二者计算结果数据的对比，发现二者的相对误差D=(前者-后者)/前者，均不超过5%，且最大结点位移相对误差为0.0039%（产生在4单元u i处），最大单元杆端力相对误差为0.0041%(产生在1、2单元杆端轴力处)，可以认为二者计算结果是相等的。
MATLAB和ANYSY计算结果结点位移相对误差
MATLAB和ANYSY计算结果单元杆端力相对误差
（4）总结
通过本次作业，本人对矩阵位移法有了更清楚的认识和掌握，对于MATLAB和ANSYS，从入门到它们在作业中的运用，有了较大的进步，为以后的学习打好了基础，同时在作业过程中，也遇到了不少的问题，比如在用MATLAB编程计算单元刚度矩阵和整体刚度矩阵的时，由于数据的繁琐，编程遇到了较大的麻烦，经过反复的思考和对程序的调试，终于柳暗花明又一村，同时发现也可对MATLAB的强大
矩阵运算功能的函数进行调用，为此带来了很大的方便；在用ANSYS 分析的过程中，一开始就忘了将杆单元进行划分，所以计算的结果差强人意，而且开始采用的是Nodes来建立节点，发现建好模型后，无法进行单元的进一步划分，后面改变策略，采用Keypoints来建立节点，就取得良好的效果，最后在模型计算结果数据的提取中有遇到了一点问题，不过通过自己不断的查阅资料，将问题一一解开，最后将两种软件的结果数据进行对比，发现几乎完全相等，最终圆满的完成了本次作业。