平面桁架的有限元法课件

平面桁架的有限元法一、建立有限元模型平面桁架是现成的有限元模型，只有节点力二、单元分析位移假设，几何方程，物理方程，单元刚度矩阵三、总刚度矩阵、总载荷列阵、引入支撑每个节点建立平衡方程集成技术（叠加法、对号入座）四、求解、结果分析，后处理平面桁架的有限元法一、建立有限元模型1一、建立有限元模型平面桁架是现成的有限元模型，只有节点力桁架原型=有限元模型一、建立有限元模型桁架原型=有限元模型2二、单元分析平面二力杆单元注意：一个单元的节点数目（2）；一个节点的自由度数目（2）；有整体坐标系和局部坐标系之分二、单元分析注意：一个单元的节点数目（2）；3二、单元分析（先在局部坐标系中进行）问题：当节点位移一定，相应的节点力为多大？节点力与节点位移之间的比例常数为刚度系数目标：建立单元刚度矩阵。方法：结构力学的位移法加虚功原理，皆用节点位移表示。二、单元分析（先在局部坐标系中进行）节点力与节点位移之间的比41、位移假设：设各点的位移，即位移场为其中为待定系数，由节点坐标和节点位移确定。将位移写成矩阵形式：1、位移假设：其中为待定系数，由节点坐标和节点位移确定。5在节点上满足：写成矩阵形式：在节点上满足：写成矩阵形式：6解线性代数方程组，得代入得解线性代数方程组，得代入得7节点位移与单元内位移的关系其中为形函数。用Machematica推导形函数：在编程中只关心节点位移的系数矩阵节点位移与单元内位移的关系其中为形函数。用Machemati82、应变（从位移求应变）只有轴向应变3、应力（从应变求应力，本构方程）2、应变（从位移求应变）3、应力（从应变求应力，本构方程）94、内力（从应力求内力）二力杆只有轴向力4、内力（从应力求内力）二力杆只有轴向力10用Machematica推导应力矩阵：(EA的两个字母连写作参数名)用内力表示的物理方程用Machematica推导应力矩阵：用内力表示的物理方程114、单元刚度矩阵（用虚功原理）给节点一组虚位移，单元内即产生相应的虚应变，实际节点力（轴向力、横向力、弯矩）在虚位移上所做的功等于实际应力（横向剪应力可忽略不计）在虚应变上产生的应变能。用Machematica推导单刚矩阵[Ke]：4、单元刚度矩阵（用虚功原理）用Machematica推导单12用Machematica推导，所得结果以矩阵形式显示用Machematica推导，所得结果以矩阵形式显示13三、总刚度矩阵、总载荷列阵、引入支撑1、坐标变换前面的推导是在局部坐标系中进行的，而总刚度矩阵、总载荷列阵、引入支撑等要在整体坐标中进行，因此，存在坐标变换问题。三、总刚度矩阵、总载荷列阵、引入支撑1、坐标变换14三、总刚度矩阵、总载荷列阵、引入支撑1、坐标变换三、总刚度矩阵、总载荷列阵、引入支撑1、坐标变换15三、总刚度矩阵、总载荷列阵、引入支撑2、在总结构和整体坐标中讨论问题，所需要的信息：（1）节点号；（2）单元长度；（2）单元倾角(单元起点到终点的有向线段与整体横轴之夹角）；（4）节点力所在的节点号和作用方向；（5）支撑所在的节点号和作用方向。3、总节点位移向量的排列规律按节点号和节点自由度的顺序（从小到大）排列。三、总刚度矩阵、总载荷列阵、引入支撑2、在总结构和整体坐标中16三、总刚度矩阵、总载荷列阵、引入支撑建立节点号数组（用Mathematica)：jm=Table[{0,0},{i,ne}];“jm数组名，ne总单元数;”jm={…，{i节点号，j节点号}，…};建立单元长度数组gc=Table[{0},{i,ne}];“gc数组名，ne总单元数;”gc={L1,L2,…Lne};建立单元方向角数组gj=Table[{0},{i,ne}];“gj数组名，ne总单元数;”gj={af1,af2,…afne};三、总刚度矩阵、总载荷列阵、引入支撑建立节点号数组（用Mat173、总节点位移向量的排列规律jm={{1,2}，{2,3}，…，{i节点号，j节点号}，…};单元节点位移排列规律：rr=2(k-1)+kkk=1,2(单元节点);kk=1,2,(节点自由度)k=1k=2kk=1kk=2总节点位移向量{d}的排列规律(从小到大排序）：ss=2(jm[[e,k]]-1)+kk;e=单元号对每个单元进行节点循环和自由度循环，从单元地址到结构地址搬家是单元到结构的映射(map)3、总节点位移向量的排列规律jm={{1,2}，{2,3183、总节点位移向量的排列规律单元节点位移号：rr=2(k-1)+kk；k=1,2(单元节点)；

kk=1,2(节点自由度)总节点位移号：ss=2(jm[[e,k]]-1)+kk；例如，已解得结构的总节点位移向量{d}求各单元的应力。用Mathematica：deg=Table[0,{i,4}];“整体坐标系中的单元节点位移向量”；del=Table[0,{i,4}];“局部坐标系中的单元节点位移向量”；del=se=Table[{0,0},{i,ne}];“单元内力数组ne行”;For[e=1,e<=ne,e++,For[k=1,k<=2,k++,For[kk=1,kk<=2,kk++,rr=2(k-1)+kk；ss=2(jm[[e,k]]-1)+kk；deg[[rr]]=d[[ss]]]];del=Transpose[T[gj[[e]]]].deg;resultant=De.BL.del;se[[e,1]]=e;se[[e,2]]=resultant[[j,1]]

];3、总节点位移向量的排列规律单元节点位移号：rr=2(k-194、总节点载荷向量的排列规律按节点号和节点自由度的顺序（从小到大）排列。设总节点力向量为{Rz}单元节点位移号：rr=2(k-1)+kk；k=1,2(单元节点)；kk=1,2(节点自由度)；总节点位移号：ss=2(jm[[e,k]]-1)+kk；4、总节点载荷向量的排列规律按节点号和节点自由度的顺序（从小204、总节点载荷向量的排列规律单元节点力号：rr=2(k-1)+kk；k=1,2(单元节点)；kk=1,2(节点自由度)；总节点力号：ss=2(jm[[e,k]]-1)+kk；用Mathematica编程：nj=ne+1总节点数Pk=Table[0,{i,4}];“单元节点力向量”；Rz=Table[0,{i,2nj}];“总节点力向量2nj行”;For[e=1,e<=ne,e++,Do[Pk[[i]]=….,{i,4}];“生成单元节点力”；For[k=1,k<=2,k++,For[kk=1,kk<=2,kk++,rr=2(k-1)+kk；ss=2(jm[[e,k]]-1)+kk；Rz[[ss]]=Rz[[ss]+Pk[[rr]];“叠加法”；]

]];4、总节点载荷向量的排列规律单元节点力号：用Mathemat215、总刚度矩阵的排列规律（1）总刚度矩阵中的行和列均按节点号和节点自由度的顺序（从小到大）排列；（2）在总刚度矩阵中，与同一节点相连的单元的刚度要叠加。对一个单元而言：ijijij5、总刚度矩阵的排列规律（1）总刚度矩阵中的行和列均按节点号225、总刚度矩阵的排列规律单刚行号：r=2(i-1)+ii；i=1,2(单元节点)；ii=1,2(节点自由度)；单刚列号：s=2(j-1)+jj；j=1,2(单元节点)；jj=1,2(节点自由度)；总刚行号：rr=2(jm[[e,i]]-1)+ii；总刚列号：ss=2(jm[[e,j]]-1)+jj；

jm=Table[{0,0},{i,ne}];“jm数组名，ne总单元数;”jm={…….，{i节点号，j节点号}，……};5、总刚度矩阵的排列规律单刚行号：r=2(i-1)+ii；23（2）在总刚度矩阵中，与同一节点相连的单元的刚度要叠加。用Mathematica编程：Kz=Table[0,{i,2nj},{j,2nj}];“开总刚度矩阵,nj总节点数”；For[e=1,e<=ne,e++,“生成单刚,变坐标系”；For[i=1,i<=2,i++,For[ii=1,ii<=2,ii++,r=2(i-1)+ii;rr=2(jm[[e,i+1]]-1)+ii；For[j=1,j<=2,j++,For[jj=1,jj<=2,jj++,s=2(j-1)+jj;ss=2(jm[[e,j+1]]-1)+jj；

“叠加法”；]]]]];单刚行号：r=2(i-1)+ii；i=1,2(单元节点)；ii=1,2(节点自由度)；单刚列号：s=2(j-1)+jj；j=1,2(单元节点)；jj=1,2(节点自由度)；总刚行号：rr=2(jm[[e,i]]-1)+ii；总刚列号：ss=32(jm[[e,j]]-1)+jj；（2）在总刚度矩阵中，与同一节点相连的单元的刚度要叠加。用M246、引入支撑(1)位移为零的约束前处理时输入信息：支撑数目nz、在整体坐标系中的自由度号zc={…,2(jm[[e,k]]-1)+kk,…};“k=1,2;kk=1,2;”如图nz=6;zc={3,4,5,6,7,8};6、引入支撑(1)位移为零的约束25For[i=1,i<=nz,i++,z=zc[[i]];Rz[[z]]=0;For[j=1,j<=2nj,j++,Kz[[z,j]]=0;For[k=1,k<=2nj,k++,Kz[[k,z]]=0]]；Kz[[z,z]]=1]6、引入支撑(1)位移为零的约束在Mathematica上实现如图:nz=6;zc={3,4,5,6,7,8};For[i=1,i<=nz,i++,z=zc[[i]]26(2)给定位移的约束如图，在节点1给一向右的水平位移u1。已知位移数目nf=1;总自由度号和位移值zf={1,u1};For[i=1,i<=nf,i++,z=zf[[i,1]];zz=zf[[i,2]];