一个变截面梁的有限元求解
本帖最后由 whutjp 于 2012-4-18 08:59 编辑这是我的matlab程序：
clc
clear
E=6.02e11; %弹性模量
L=[0.02 0.18 0.015 0.335 0.39 0.226 0.16 0.64 0.089 0.03 0.66 0.66 0.701 0.2 0.2...
0.086 0.175 0.2 0.953 0.4 0.4 0.088 0.424 0.04 0.02 0.06 0.08 0.055 0.045 0.127...
0.07 0.1 0.049 0.135 0.135 0.0375 0.1975 0.1975 0.0725 0.135 0.135 0.005];%单元长度
D=[0.225 0.24 0.225 0.2738 0.2905 0.31 0.31 0.31 0.31 0.31 0.3 0.3 0.3 0.308 0.3 0.3 0.3 0.3 0.3...
0.285 0.285 0.285 0.455 0.383 0.383 0.65 0.65 0.34 0.34 0.35 0.34 0.71 0.34 0.35 0.35...
0.42 0.42 0.42 0.35 0.35 0.35 0.35]; %单元外径
d=;%单元内径
F=;%节点外力载荷
M=;%节点外加弯矩载荷
s_zc=;%支承所在节点
q=7.85.*9.8e3.*pi.*(D.^2-d.^2)/4; %单元均布载荷
I=pi.*(D.^4-d.^4)./64; %单元截面惯性矩
n=size(L,2);%获取所划分的单元长度
for i=1:1:n
k(:,:,i)=BeamElementStiffness(E,I(i),L(i));
end %每个单元的刚度矩阵
K=zeros(2*n+2,2*n+2); %定义整体刚度矩阵
for i=1:1:n
K=BeamAssemble(K,k(:,:,i),i,i+1);
end %将单元刚度矩阵整合入整体刚度矩阵
K1=SparseMatrixChange(K,s_zc);%应用波阵法进行处理
WaiLi=BeamWaiLi(L,D,d,F,M);
n1=size(s_zc,2); %获取支承个数
WaiLi1=WaiLi';
for i=1:1:n1
WaiLi1(s_zc(i),:)=[];
end %应用波阵法进行处理
u=K1\WaiLi1; %解方程
运行后提示：Warning: Matrix is singular to working precision. 然后出来的结果都是NaN啊
整体思路就是：将变截面的轴系分为很多段，轴系的重量当均布载荷处理为节点的应力和弯矩！最后求解的时候我已经加了边界条件（轴承处位移为0）了，怎么还是无法求解呢？
最好把你的有限元推导过程列出来，这样方便大家查找错误。 tonnyw 发表于 2012-4-17 23:55 static/image/common/back.gif
最好把你的有限元推导过程列出来，这样方便大家查找错误。
我已经把推导过程也上传了！我完全是按照这个推导过程来的，怎么我的就老解不出来呢？ How about starting with a simple Euler-Bernoulli beam so we can understand things better? Once you are confident with the Euler-Bernoulli, you can go ahead with Timoshenko beam. tonnyw 发表于 2012-4-18 23:40 static/image/common/back.gif
How about starting with a simple Euler-Bernoulli beam so we can understand things better? Once you a ...
Thanks for your advice, let me try this Timoshenko beam.
页:
[1]