资源描述：
《Euler-Bernoulli beam theroy(欧拉伯努利梁理论与有限元计算).doc》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。
1、Euler-Bernoullibeam一、理论部分Euler-Bernoullibeam假设(1.1)由(1.1)式可得(1.2)虚位移原理(1.3)其中(1.4)令，为单元长度，则上式成为(1.5)单元节点位移取为(1.6)令(1.7)其中形函数(1.8)(1.9)分别对式(1.7)、(1.8)求一阶和两阶导数得(1.10)(1.11)由，可得(1.12)(1.13)将(1.6)、(1.11)式代入(1.4)式，可得(1.14)刚度矩阵(1.15)(1.16)(1.17)二、算例本文所举算例为悬臂梁端部受
2、竖向载荷，挠度理论曲线为具体数据如下：，，，，，划分10个单元，11个结点。由上图可以看出，数值模拟结构跟理论解吻合非常好。源代码(Matlab)%program_maininitial_statistic;le=pretreatment(ND,NE,lamda,x,y);[Dis,K]=Euler_Bernoulli_beam(ND,NE,lamda,EA,EI,le,P,Node_c,Dis_c);post(Dis,x,y,ND,NE);%initialstatistic;ND=11;NE=10;EA=
3、2.1e8*ones(NE,1);EI=4.375e4*ones(NE,1);x=(0:0.1:1)';y=zeros(11,1);lamda=[1,2;2,3;3,4;4,5;5,6;6,7;7,8;8,9;9,10;10,11;];P=zeros(3*ND,1);P(32)=-1.0e4;Node_c=[1,2,3]';Dis_c=[0,0,0]';functionle=pretreatment(ND,NE,lamda,x,y)%theprogramofpretreatment;fori=1:NEle(
4、i)=sqrt((x(lamda(i,1))-x(lamda(i,2)))^2+(y(lamda(i,1))-y(lamda(i,2)))^2);endendfunction[Dis,K]=Euler_Bernoulli_beam(ND,NE,lamda,EA,EI,le,P,Node_c,Dis_c)%solvetheEulerBernoullibeamproblembythefiniteelementmethod;%ND,NE,lamdatheinformationofNode,Element;%EA,
5、EI,lethepropertiesofmaterial;K=Stiffness(ND,NE,lamda,EA,EI,le);[K,P]=Import_c(K,P,Node_c,Dis_c);Dis=inv(K)*P;EndfunctionK=Stiffness(ND,NE,lamda,EA,EI,le)%集成总体刚度矩阵；%ND节点数，NE单元数，lamda单元节点定位向量；K=zeros(3*ND);Ke=zeros(6);forie=1:NEKe=Estiffness(EA(ie),EI(ie),le
6、(ie),2);l1=lamda(ie,1);l2=lamda(ie,2);K(3*l1-2:3*l1,3*l1-2:3*l1)=K(3*l1-2:3*l1,3*l1-2:3*l1)+Ke(1:3,1:3);K(3*l1-2:3*l1,3*l2-2:3*l2)=K(3*l1-2:3*l1,3*l2-2:3*l2)+Ke(1:3,4:6);K(3*l2-2:3*l2,3*l1-2:3*l1)=K(3*l2-2:3*l2,3*l1-2:3*l1)+Ke(4:6,1:3);K(3*l2-2:3*l2,3*l2-2
7、:3*l2)=K(3*l2-2:3*l2,3*l2-2:3*l2)+Ke(4:6,4:6);endendfunctionKe=Estiffness(EA,EI,le,N)%求单元刚度矩阵；%N为高斯点个数；[x,Hw]=Ngauss(N);Ke=zeros(6);fori=1:NBu=[-0.5,0,0,0.5,0,0]';Bv=[0,1.5*x(i),(0.75*x(i)-0.25)*le,0,-1.5*x(i),(0.75*x(i)+0.25)*le]';Ke=Ke+(Bu*Bu'*2/le*EA+Bv
8、*Bv'*(2/le)^3*EI)*Hw(i);endendfunction[x,Hw]=Ngauss(N)%给定高斯点个数返回高斯点出坐标和权系数值；switchNcase1x=0;Hw=2;case2x=[-0.577350269189626,0.577350269189626]';Hw=[1,1]';case3x=[-0.774596669241483,0,0.774596669241483]';Hw=[0
