第二章平面桁架有限元分析及程序设计

有限单元法及程序设计第一章绪论第二章平面桁架有限元分析第三章平面刚架有限元分析第四章平面问题有限元分析有限单元法及程序设计第一章绪论§1.1有限单元法的概念§1.2有限单元法的基本步骤§1.3常用有限单元分析软件第一章绪论数值方法（模拟）有限元方法边界元方法有限差分方法有限体积法无网格方法§1.1有限单元法的概念基本思路：借助数学和力学知识,利用计算机技术解决工程技术问题有限元分析是利用数学近似分析方法对真实物理系统（几何、载荷工况）进行模拟，利用简单而又相互作用的元素，即单元，用有限数量的未知量去逼近无限未知量的真实系统。有限单元法(FEM)是20世纪50年代以来随着计算机的广泛应用而发展起来的一种现代数值解法。该方法首先应用在连续力学领域——飞机结构静、动态特性分析中。随后很快就广泛应用于求解传导、电磁场、流体力学等连续性问题。Rxyo真实半圆周长：近似半圆周长：数学：离散化的手段，将偏微分方程、变分方程变成代数方程；计算曲线长度：将曲线离散成N段，每段用直线长度代替曲线长度，直线长度总和即可曲线长度，当N为无穷大时，即为精确解。§1.1有限单元法的概念单元数计算值理论值误差(%)123.1415936.3422.828433.141599.9732.976073.141595.2743.035283.141593.3853.065963.141592.4163.084253.141591.8373.096193.141591.4583.10453.141591.1893.110543.141590.99103.115113.141590.841003.140763.141590.032003.14133.141590.013003.141433.141590.01§1.1有限单元法的概念圆周率的计算方法：单元数计算值单元数计算值单元数计算值121003.14076059100003.1415918220398022.8284271252003.141298567200003.1415923595933632.9760677433003.14143259300003.1415924935585943.035276184003.141488694400003.1415925496465253.0659554365003.141518269500003.1415925792140963.084252856003.141536068600003.1415925970102973.0961949547003.14154775700003.1415926086905583.1044960688003.141555901800003.1415926168403593.1105417379003.141561854900003.14159262279185103.11510595110003.1415663561000003.14159262729364小数点后面20位：3.14159265358979323846§1.1有限单元法的概念力学：离散化的手段，利用单元力学特性，将连续体离散成有限个单元的组合结构；齿轮的应力分布：将齿轮分割成很多单元，由结点相连并传递求解信息，用有限个连续体的组合代替原来的齿轮结构。§1.1有限单元法的概念电动铲运机§1.1有限单元法的概念液压挖掘机动臂§1.1有限单元法的概念在建的大连国际贸易中心大厦(78层，342米)63863个梁柱单元；34180个结点；§1.1有限单元法的概念在建的大连市体育中心§1.1有限单元法的概念§1.1有限单元法的概念溪洛渡拱坝(高297米)456个坝体单元，1040个地基单元，2163个结点；§1.1有限单元法的概念§1.1有限单元法的概念导管架码头结构§1.1有限单元法的概念高桩码头高桩码头下部结构高桩码头上部结构§1.1有限单元法的概念游艇码头结构§1.1有限单元法的概念重力式码头结构有限元模型是真实系统理想化的数学抽象。真实系统有限元模型有限元模型§1.1有限单元法的概念单元网格划分中每一个小的块体结点确定单元形状、单元之间相互联结的点结点力单元上结点处的结构内力载荷作用在单元结点上的外力（集中力、分布力）约束限制某些结点的某些自由度

单元单元载荷结点结点力约束有限元单元模型中几个重要概念§1.1有限单元法的概念§1.2有限单元法的基本步骤一、结构离散化将求解域离散成单元表示的组合体，单元以结点相连；二、选择插值函数（位移模式）选择插值函数来表达单元内场变量的变化规律三、形成单元性质的矩阵方程场变量：标量向量力、位移张量应力、应变单元刚度方程：单元结点力与单元结点位移之间的关系；§1.2有限单元法的基本步骤四、形成整体系统的矩阵方程综合所有单元性质的矩阵，形成整体系统的矩阵方程；五、约束处理，求解系统方程引入边界条件，求解结点上的未知场变量；六、其他参数的计算利用已经求解的场变量，计算其他场变量；§1.3常用的有限元分析软件(1)ANSYS功能强大、模块多、比较通用；土木工程：CivilFEM是商业化比较早的一个软件（收购了一些其他软件公司）(2)MSC产品系列多、通用软件、二次开发功能强；土木工程：MSC.MARC系列，Patran,MSCNastran及AdamsAdams------多体动力学；Actran------声学仿真；Easy5-----控制仿真工具；Marc------非线性；SimXpert------多学科仿真；MSCNastran------结构化与多学科FEA；Dytran------显式动力学与流固耦合；MSCFatigue------基于FE的耐久性仿真工具；Sinda------高级热分析解决方案；Digimat------非线性，多尺度的材料与结构建模平台；SimDesigner------CAD嵌入式多学科仿真；Patran------有限元分析解决方案；SimManager------仿真数据和流程管理§1.3常用的有限元分析软件(3)ADINA复杂非线性问题与复杂动力问题；具有直接求解和迭代求解两种流固耦合分析方法AutomaticDynamicIncrementalNonlinearAnalysis(4)ABAQUS（专注结构分析）主要应用于结构动力弹塑性分析；二次开发功能；HKS公司产品§1.3常用的有限元分析软件(5)SAP2000SystemsApplicationsandProductsinDataProcessing通用结构分析

第二章平面桁架有限元分析及程序设计§2.1平面桁架单元的离散§2.2平面桁架单元分析§2.3结点平衡与整体刚度矩阵的集成§2.4边界条件的处理§2.5单元内力与支座反力的计算§2.6平面桁架有限元程序设计有限单元法及程序设计解题方法ⅠⅠ方法1：节点法BaFP1ACDaaaFP方法2：截面法静定桁架回顾第二章平面桁架有限元分析及程序设计解题方法：力法和位移法超静定桁架P1234第二章平面桁架有限元分析及程序设计如图所示桁架，求各杆轴力。力的平衡条件：位移的协调方程：PN1N1N21杆和3杆位移：2杆位移：超静定桁架第二章平面桁架有限元分析及程序设计1杆轴力竖向分量：2杆轴力：式中：和为杆件的刚度系数；物理意义：4点产生单位位移，杆端产生的竖向杆端力；由杆件的物理性质和几何性质决定；V4为第4节点竖向位移P1234第二章平面桁架有限元分析及程序设计超静定桁架代入平衡方程：结构的整体刚度系数位移法求解超静定结构。离散原则：每个结点离散后还是一个结点，每个杆件离散后变成一个单元1结构的离散化：尽量将结构离散成数量最少的等截面直杆单元23456⑤①②③④⑥⑦⑧⑨§2.1平面桁架单元的离散9个单元，6个结点123456⑤①②③④⑥⑦⑧⑨78⑩15121413111616个单元，8个结点§2.2.1

局部坐标系下的单元刚度矩阵局部坐标系的建立iE，A，lje■

轴：沿单元的杆轴方向；

§2.2平面桁架的单元分析■

轴：从轴逆时针旋转90°。■

原点：以第一个结点为坐标原点；

杆端位移：ijeijeije杆端力：符号：与坐标系的方向一致为正，反之为负。单元右端杆端力：单元左端杆端力：单元应力：单元应变：右结点固定结点位移：左结点固定杆的受力分为两种情况：§2.2平面桁架的单元分析任意情况（左右结点均有变形）即为以上两种状态的叠加：杆端力为：式中为单元刚度矩阵(局部坐标系)杆单元轴力为：式中为单元应力(广义)矩阵；§2.2平面桁架的单元分析单元杆端力方程：杆端位移：ijeije杆端力：单元轴力：§2.2平面桁架的单元分析杆端位移和杆端力符号：与坐标系的方向一致为正，反之为负。杆端力：xy§2.2平面桁架的单元分析§2.2.2

整体坐标系下的单元刚度矩阵若局部坐标系与整体坐标系重合，则整体坐标系下的单元刚度矩阵与局部坐标下的单元刚度矩阵相同。若局部坐标系与整体坐标系不重合，如下图所示：杆端位移：杆端力：杆端位移：i

结点：j

结点：§2.2平面桁架的单元分析设杆件的长度为l，则：两边微分：由于杆件的变形产生位移：因此，杆件应变为：杆件轴力为：符号：杆件轴力以拉为正，压为负。杆件的结点力为：因此，杆件结点力向量为：式中是整体坐标系下的单元刚度矩阵；§2.2平面桁架的单元分析写成分块矩阵形式：式中：§2.2平面桁架的单元分析（1）单元刚度系数kij的意义j自由度(结点)产生的单位杆端位移引起的i自由度(结点)的杆端力（2）单元刚度矩阵是对称矩阵反力互等定理式中：杆件单元的应力矩阵为：单元刚度矩阵的性质§2.2平面桁架的单元分析（3）单元刚度矩阵一般是不可逆的§2.2.3单元坐标转换矩阵§2.2平面桁架的单元分析取任意杆件,建立如图所示的局部坐标系：杆端力：杆端位移：§2.2.3单元坐标转换矩阵xyxy§2.2平面桁架的单元分析杆端力：杆端位移：在上图中,建立如图所示的整体坐标系：以i结点为例：同理，对于j结点：§2.2平面桁架的单元分析写成矩阵形式：因此：其中，[T]为转换矩阵：转换矩阵的性质转换矩阵是正交矩阵；同理，位移也存在转换关系：代入局部坐标系下的刚度方程：§2.2平面桁架的单元分析与利用微分得到的单元在总体坐标下的刚度方程相同xyF1、对总体结点位移和单元进行编码；2、单元局部坐标系下的刚度矩阵;①②③123例：如图所示平面桁架，杆长为l，截面积为A，求三个单元在整体坐标系下的刚度矩阵。§2.2平面桁架的单元分析3、单元①整体坐标系分析:解：单元①整体坐标系下的刚度矩阵为:4、单元②整体坐标系分析:5、单元③整体坐标系分析:xyF①②③123例：如图所示平面桁架，杆长为l，截面积为A，求结构的刚度矩阵。§2.3结点平衡与整体刚度矩阵的集成1、单元整体坐标系下刚度矩阵分块解：§2.3.1结点的平衡方程§2.3结点平衡与整体刚度矩阵的集成xyF1①ij③ij2、结点1的平衡方程：结点1的受力状态为（如右图）:结点1的平衡条件为:由单元③的刚度方程：由单元①的刚度方程：§2.3结点平衡与整体刚度矩阵的集成代入结点1的平衡条件:3、结点2的平衡方程：同理，结点2的平衡条件为:由单元①的刚度方程：由单元②的刚度方程：§2.3结点平衡与整体刚度矩阵的集成代入结点2的平衡条件:4、结点3的平衡方程：同理，结点3的平衡条件为:由单元③的刚度方程：由单元②的刚度方程：§2.3结点平衡与整体刚度矩阵的集成代入结点3的平衡条件:5、系统的平衡方程：§2.3结点平衡与整体刚度矩阵的集成6、结构整体刚度方程写成矩阵形式，即可得到结构的整体刚度方程其中，[K]为结构的整体刚度矩阵；§2.3.2

整体刚度矩阵的集成步骤1、定位单元结点编号

2、累加整体结点编号

单元刚度系数

整体刚度系数§2.3结点平衡与整体刚度矩阵的集成单元定位向量xyF①②③123§2.3结点平衡与整体刚度矩阵的集成例：求上例平面桁架结构的整体刚度矩阵；1、定位单元①：1

2

3123(2)(1)(2)

(1)2、累加§2.3结点平衡与整体刚度矩阵的集成1、定位单元③：(2)(1)(2)

(1)2、累加1、定位单元②：2、累加1

2

3123

(1)

(2)(1)(2)1

2

3123§2.3结点平衡与整体刚度矩阵的集成结构的整体刚度矩阵为:将每个字块展开，结构的整体刚度矩阵为:下标表示自由度编号下标表示结点编号§2.3结点平衡与整体刚度矩阵的集成1、整体刚度矩阵的性质§2.3.2整体刚度矩阵的集成方法（1）刚度系数Kij的意义分块矩阵：（2）单元刚度矩阵是对称矩阵反力互等定理j结点产生的单位杆端位移引起的i结点的杆端力；j自由度产生的单位杆端位移引起的i自由度的杆端力；不分块矩阵：§2.3结点平衡与整体刚度矩阵的集成思考题：如何由刚度矩阵元素的意义确定整体刚度矩阵每个元素的组成？§2.3结点平衡与整体刚度矩阵的集成1234⑤①②③④⑥练习题：利用整体刚度矩阵的意义确定以下桁架刚度矩阵元素（分块）的组成。§2.3.2

整体刚度矩阵的集成步骤1、定位单元自由度编号

2、累加整体自由度编号

§2.3结点平衡与整体刚度矩阵的集成自由度定位向量单元刚度系数

整体刚度系数练习题：利用整体刚度矩阵的意义确定以下桁架刚度矩阵元素（自由度）的组成。§2.4.1

结点边界条件§2.4边界条件的处理结点可以自由变形，整体结点力等于对应的外荷载。1、自由变形的结点2、约束结点或给定了结点位移的数值约束结点：给定结点位移：§2.4.2

边界条件的处理方法§2.4边界条件的处理1）划行划列法处理方法：若第i个自由度位移为零，则将总刚第i行和第i列划掉，刚度矩阵相应降低一阶。000优点：简单易行，矩阵降阶，减小计算工作量；缺点：矩阵行列、位移、荷载向量需重新编号，程序实现比较复杂；只适用于约束结点情况；§2.4边界条件的处理2）0、1置换法（填0置1法）处理方法：将与约束自由度对应整体刚度矩阵对角线元素全部置换成1，相应行和列其他元素置换成0，将同一行荷载分量置换成0。000111000000000000000000000000000适用条件：只适用于约束结点，不适用给定位移边界条件；§2.4边界条件的处理3）乘大数法处理方法：将与约束自由度对应整体刚度矩阵对角线元素乘以一个大数N(1010-1015)，将同一行荷载分量置换成N与对角线元素的乘积与给定位移之积。NNNNK33b1NK44b2NK66b3优点：处理工作量小，适用于给定位移情况；000§2.5单元轴力及支座反力的计算单元轴力：(1)整体坐标系下的单元杆端位移(2)整体坐标系下的单元轴力§2.5.1

单元轴力的计算位移和荷载向量分解为自由结点和约束结点两部分，刚度矩阵相应分块，如下所示：§2.5.2

支座反力的计算其中：为自由位移;为约束位移;为外荷载;为约束自由度结点力;因此：为支座反力;为约束自由度结点荷载;其中：§2.5单元轴力及支座反力的计算1、计算分析题平面桁架结构网格如图所示，已知EA=1500kN，采用乘大数方法引入支撑条件，试求后处理法引入支撑条件后的满阵存贮的整体刚度矩阵[K]。结构坐标系下单元刚度矩阵计算公式如下：

例题式中：321①②③3myx4m2、计算分析题按照有限元法的计算步骤，求图示桁架结构各杆轴力。已知：EA=10kn。结构坐标系下的单元刚度矩阵计算公式如下：：

例题式中：§2.6程序设计程序设计原则完整性扩充性兼容性逻辑性可读性可维护性模块化§2.6程序设计程序设计流程图（程序框图）什么是流程图？

“程序流程图”常简称为“流程图”，是一种传统的算法表示法，程序流程图是人们对解决问题的方法、思路或算法的一种描述。它利用图形化的符号框来代表各种不同性质的操作，并用流程线来连接这些操作。2、如何画流程图编码和单元测试这个阶段的任务是程序员根据目标系统的性质和实际环境，选取一种适当的高级程序设计语言(必要时用汇编语言)，把详细设计的结果翻译成用选定的语言书写的程序，并且仔细测试编写出的每一个模块。

程序员在书写程序模块时，应使它的可读性、可理解性和可维护性良好。综合测试这个阶段的任务是通过各种类型的测试，使软件达到预定的要求。

最基本的测试是集成测试和验收测试。集成测试是根据设计的软件结构，把经单元测试的模块按某种选定的策略装配起来，在装配过程中对程序进行必要的测试。验收测试是按照需求规格说明书的规定，由用户对目标系统进行验收。

通过对软件测试结果的分析可以预测软件的可靠性；反之，根据对软件可靠性的要求也可以决定测试和调试过程什么时候可以结束。

在进行测试的过程中，应该用正式的文档把测试计划、详细测试方案以及实际测试结果保存下来，作为软件配置的一部分。

§2.6程序设计§2.6.1程序框图输入数据单元局部刚度坐标转换矩阵单元整体刚度集成整体刚度矩阵元素约束条件处理、解方程计算单元轴力、约束反力单元循环包括单元、结点、材料、荷载、约束数据§2.6

程序设计§2.6.2

程序说明1、总体刚度矩阵的半带宽存储总体刚度矩阵：对称稀疏矩阵；00半带宽：总体刚度矩阵集成：§2.6

程序设计行号：000列号：主对角线半带宽存储下三角：§2.6

程序设计2、先处理法处理边界条件单元定位向量：xyF①②③123(1,2)(3,4)(5,6)后处理：xyF①②③123(1,2)(0,0)(3,0)划行划列法先处理法：单元定位向量：

voidforce(){int

i,j,ie,m;floatdx,dy,dz,l,cx,cy,cz,ea,w[7];

for(ie=1;ie<=ne;ie++){

i=jm[ie][1];

j=jm[ie][2];

m=jm[ie][0];w[1]=f[2*i-2];

w[2]=f[2*i-1];

w[3]=f[2*j-2];

w[4]=f[2*j-1];

dx=xy[j][1]-xy[i][1];

dy=xy[j][2]-xy[i][2];

例题3、程序说明下面为一个平面桁架计算程序段，试在左端有编号的程序右面写出其注释。（10分）

l=sqrt(dx*dx+dy*dy);

cx=dx/l;cy=dy/l;ea=EA[m]/l;

dx=w[3]-w[1];

dy=w[4]-w[2];

l=ea*(cx*dx+cy*dy