有限元与有限差分法基础超详细版本精选PPT.ppt

《有限元与有限差分法基础超详细版本精选PPT.ppt》由会员分享，可在线阅读，更多相关《有限元与有限差分法基础超详细版本精选PPT.ppt（162页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、2023/4/31/23关于有限元与有限差分法基础超详细版本第1页，讲稿共162张，创作于星期二2023/4/32/162有限元法基础n有限元发展过程n有限元应用n有限元发展方向第2页，讲稿共162张，创作于星期二2023/4/33/162有限元法的基本思想n基本思想 1）将连续的求解系统离散为一组由节点相互联在一起的单元组合体 2）在每个单元内假设近似函数来分片表示系统的求解场函数 第3页，讲稿共162张，创作于星期二2023/4/34/162有限元法的基本思想第4页，讲稿共162张，创作于星期二2023/4/35/162有限元法的基本思想第5页，讲稿共162张，创作于星期二2023/4/3

2、6/162有限元法的基本思想第6页，讲稿共162张，创作于星期二2023/4/37/162有限元法的基本思想第7页，讲稿共162张，创作于星期二2023/4/38/162有限元法的基本思想n离散为单元网格的件仍然要保证是一个连续体，单元与单元之间没有裂缝、不能重叠，所有单元通过单元节点相互关联着n变形体无论产生多大的塑性变形，单元与单元之间依然不会产生裂缝、交叉和重叠，关联单元的节点也不能脱开第8页，讲稿共162张，创作于星期二2023/4/39/162有限元法的基本思想n不合格单元单元裂缝单元裂缝单元重叠单元重叠第9页，讲稿共162张，创作于星期二2023/4/310/162有限元法的基本思

3、想n变形前后单元之间都是连续的变形前的网格变形前的网格变形后的网格变形后的网格第10页，讲稿共162张，创作于星期二2023/4/311/162有限元法的基本思想n基本思想n通过在单元内假设不同的插值函数，建立不同的单元模型，适应各种各样的变形模式和受力模式XFXF第11页，讲稿共162张，创作于星期二2023/4/312/162有限元法的基本思想n有限元法分类 1）位移法：基于最小势能原理或虚功原理 2）力法：基于最小余能原理 3）混合法：基于修正余能原理 第12页，讲稿共162张，创作于星期二2023/4/313/162有限元法的基本思想n基本过程离散化过程 约束处理过程 单元平衡方程组装

4、过程 应变、应力回代过程 方程组求解过程 第13页，讲稿共162张，创作于星期二2023/4/314/162离散化过程n最小势能原理 弹性体的势能为弹性体变形后所具有的内能 为弹性体所受的外力功 第14页，讲稿共162张，创作于星期二2023/4/315/162离散化过程 为弹性体的应变 为弹性体的应力 u为弹性体的可容位移 弹性体处于平衡状态时，其势能应为最小 0第15页，讲稿共162张，创作于星期二2023/4/316/162离散化过程n单元插值关系 n单元几何关系 n单元本构关系 N为单元形函数矩阵 L为单元几何微分算子 为单元弹性矩阵 单元节点自由度向量第16页，讲稿共162张，创作于

5、星期二2023/4/317/162离散化过程B 称为应变矩阵 单元平衡方程或单元刚度方程 k 称为单元刚度矩阵 f 称为单元载荷向量 第17页，讲稿共162张，创作于星期二2023/4/318/162单元刚度矩阵的特性 n对称性 n奇异性 n主元恒正且对角占优 离散化过程第18页，讲稿共162张，创作于星期二2023/4/319/162线弹性问题几何方程三维问题 三维问题第19页，讲稿共162张，创作于星期二2023/4/320/162线弹性问题几何方程二维问题 二维问题平面应力和平面应变状态 第20页，讲稿共162张，创作于星期二2023/4/321/162线弹性问题几何方程二维问题 二维问

6、题轴对称状态 第21页，讲稿共162张，创作于星期二2023/4/322/162线弹性问题几何方程一维问题 一维问题第22页，讲稿共162张，创作于星期二2023/4/323/162线弹性问题本构方程三维问题 三维问题E为弹性模量；为泊松比 第23页，讲稿共162张，创作于星期二2023/4/324/162线弹性问题本构方程平面应力 二维问题平面应力状态 第24页，讲稿共162张，创作于星期二2023/4/325/162线弹性问题本构方程平面应力 平面应力状态 第25页，讲稿共162张，创作于星期二2023/4/326/162线弹性问题本构方程平面应变 二维问题平面应变状态 第26页，讲稿共1

7、62张，创作于星期二2023/4/327/162线弹性问题本构方程平面应变 平面应变状态 第27页，讲稿共162张，创作于星期二2023/4/328/162线弹性问题本构方程轴对称 二维问题轴对称状态 第28页，讲稿共162张，创作于星期二2023/4/329/162线弹性问题本构方程轴对称 二维问题轴对称状态 第29页，讲稿共162张，创作于星期二2023/4/330/162线弹性问题本构方程轴对称 轴对称状态 第30页，讲稿共162张，创作于星期二2023/4/331/162线弹性问题本构方程一维问题 一维问题第31页，讲稿共162张，创作于星期二2023/4/332/162常用单元模型

8、n单元模型插值关系一一对应n单元类型一维单元、二维单元、三维单元等参单元、超参单元、次参单元第32页，讲稿共162张，创作于星期二2023/4/333/162常用单元模型n一维单元 2节点线单元3节点线单元梁单元第33页，讲稿共162张，创作于星期二2023/4/334/162常用单元模型n二维单元3节点三角形线性单元6节点三角形二次单元第34页，讲稿共162张，创作于星期二2023/4/335/162常用单元模型n二维单元10节点三角形三次单元4节点四边形双线性单元第35页，讲稿共162张，创作于星期二2023/4/336/162常用单元模型n二维单元8节点四边形二次单元12节点四边形三次单

9、元第36页，讲稿共162张，创作于星期二2023/4/337/162常用单元模型n三维单元4节点四面体线性单元10节点四面体二次单元第37页，讲稿共162张，创作于星期二2023/4/338/162常用单元模型n三维单元8节点六面体线性单元20节点六面体二次单元第38页，讲稿共162张，创作于星期二2023/4/339/162常用单元模型n准三维空间单元桁架单元一维2节点线单元+单元局部随体坐标系 为什么要建立单元局部随体坐标系？1.简化分析问题的复杂程度。2.在局部坐标系中，空间桁架的每根杆每变成了一维2节点线单元第39页，讲稿共162张，创作于星期二2023/4/340/162常用单元模型

10、n准三维空间单元框架单元三维梁单元+一维2节点线单元+单元局部随体坐标系 两端都是刚性联结 可以要承受拉压、弯曲、扭转3种变形模式 框架单元的特点第40页，讲稿共162张，创作于星期二2023/4/341/162常用单元模型n准三维空间单元板单元薄板单元中厚板单元弯曲和横向剪切2种变形模式抵抗板的变形如果板很薄，忽略横向剪切抗力，认为抵抗载荷的主要因素是弯矩第41页，讲稿共162张，创作于星期二2023/4/342/162常用单元模型n准三维空间单元壳单元 抵抗拉压变形的二维单元+板单元+单元局部随体坐标系。适合于薄壳单元和中厚壳单元从几何上分为薄壳单元和中厚壳单元 组合单元第42页，讲稿共1

11、62张，创作于星期二2023/4/343/162常用单元模型n准三维空间单元 壳理论单元 由空间壳理论严格构造的壳单元。适合于薄壳单元和中厚壳单元 退化单元 由三维实体单元退化成的壳单元。只适合于中厚壳单元 第43页，讲稿共162张，创作于星期二2023/4/344/162单元模型构造 n有限元法的基本思想 通过单元分片近似，在每个单元内假设近似函数来分片表示系统的场函数 n选择近似函数简单、实用的原则在有限元法中，近似函数称为插值函数 第44页，讲稿共162张，创作于星期二2023/4/345/162单元模型构造n插值函数 一般都采用多项式函数，主要原因是:采用多项式插值函数比较容易推导单元

12、平衡方程，特别是易于进行微分和积分运算。随着多项式函数阶次的增加，可以提高有限元法的计算精度。从理论上说，无限提高多项式的阶数，可以求得系统的精确解。第45页，讲稿共162张，创作于星期二2023/4/346/162单元模型构造方法 n整体坐标系法n局部坐标系法 nLagrange插值方法nHermite插值方法第46页，讲稿共162张，创作于星期二2023/4/347/162单元模型构造方法n2节点线单元12 oxu1u2x1x2ux1.假设插值多项式2.利用节点值求 a0 和 a1 第47页，讲稿共162张，创作于星期二2023/4/348/162单元模型构造方法3.代入a0 和 a1，得

13、插值多项式 u(x)4.按u1 和 u2合并同类项，设 l=x2-x1第48页，讲稿共162张，创作于星期二2023/4/349/162单元模型构造方法n关键 如何构造插值多项式 u?二维问题三维问题，如何构造插值多项式?第49页，讲稿共162张，创作于星期二2023/4/350/162n收敛性条件 在单元内，场函数必须是连续的；完备性：插值多项式的阶次必须由低到高依次增加，不能出现跳跃现象；协调性：各单元边界必须连续，单元边界不能出现开裂现象。插值多项式收敛性条件 n收敛：当单元逐渐缩小时，如果插值多项式满足收敛性条件，则数值解将收敛于精确解 第50页，讲稿共162张，创作于星期二2023/

14、4/351/162插值多项式收敛性条件n协调单元 满足插值多项式收敛性条件和的单元 n完备单元 满足插值多项式收敛性条件的单元ncr 阶连续性 插值多项式的第r阶导数是连续的 第51页，讲稿共162张，创作于星期二2023/4/352/162插值多项式收敛性条件n非协调单元与部分协调单元 对于一般固体力学问题来说，协调性要求单元在变形时，相邻单元之间不应引起开裂、重叠或其它不连续现象。例如，梁、板、壳等单元，在单元边界不但要求位移是连续的，而且其一阶导数也必须是连续的。板、壳单元位移函数沿单元边界的法向导数（转角）的连续性一般比较难实现，因此出现了许多不完全满足协调性要求的“非协调单元”或“部

15、分协调单元”，有时它们的精度也很好。第52页，讲稿共162张，创作于星期二2023/4/353/162插值多项式选择条件 n插值多项式应该尽可能满足其收敛性条件（收敛性）n由插值多项式所确定的场函数变化应该与局部坐标系的选择无关（各向同性）n假设的插值多项式系数的数量应该等于单元的节点数（解的唯一性）选择条件第53页，讲稿共162张，创作于星期二2023/4/354/162插值多项式选择条件n深入分析由收敛性条件可知，插值多项式中必须含有常数项（刚体位移项），高阶项的次数必须依次增加，不允许有跳跃第54页，讲稿共162张，创作于星期二2023/4/355/162插值多项式选择条件由选择条件可知

16、，插值多项式函数在所有自由度方向上要满足各向同性性，这样就不会随局部坐标系变化而改变了 n深入分析第55页，讲稿共162张，创作于星期二2023/4/356/162插值多项式选择条件n深入分析选择条件是为了能由单元节点值唯一确定插值多项式 4节点四边形的插值多项式应该是 插值多项式系数i(i=0,1,2,3)也是4个 第56页，讲稿共162张，创作于星期二2023/4/357/162单元模型构造整体坐标系法n基本思想 针对弹性体有限元网格建立一个统一的坐标系，每个单元的插值多项式都在这个坐标系上建立 第57页，讲稿共162张，创作于星期二2023/4/358/162单元模型构造整体坐标系法n2

17、节点线单元12 oxu1u2x1x2ux1.假设插值多项式2.利用节点值求 a0 和 a1 第58页，讲稿共162张，创作于星期二2023/4/359/162单元模型构造整体坐标系法3.代入a0 和 a1，得插值多项式 u(x)4.按u1 和 u2合并同类项，设 l=x2-x1第59页，讲稿共162张，创作于星期二2023/4/360/162单元模型构造整体坐标系法N1 和 N2 称为单元的形函数；N 称为单元的形函数矩阵；ue 称为单元节点位移向量。n2节点线的单元形函数第60页，讲稿共162张，创作于星期二2023/4/361/162单元模型构造整体坐标系法n二维3节点三角形单元 建立整体

18、坐标系oxy 第61页，讲稿共162张，创作于星期二2023/4/362/162单元模型构造整体坐标系法1.假设插值多项式2.首先，利用节点值求 0、1 和 2 n二维3节点三角形单元 第62页，讲稿共162张，创作于星期二2023/4/363/162单元模型构造整体坐标系法A为单元面积第63页，讲稿共162张，创作于星期二2023/4/364/162单元模型构造整体坐标系法3.将 0、1 和 2 代入插值多项式，按u1、u2、u3合并同类项第64页，讲稿共162张，创作于星期二2023/4/365/162单元模型构造整体坐标系法4.同理可得第65页，讲稿共162张，创作于星期二2023/4/

19、366/162单元模型构造整体坐标系法5.单元插值多项式为第66页，讲稿共162张，创作于星期二2023/4/367/162单元模型构造整体坐标系法6.单元插值多项式写成矩阵形式（常用）第67页，讲稿共162张，创作于星期二2023/4/368/162单元模型构造整体坐标系法7.单元插值多项式的另一种矩阵形式（不常用）第68页，讲稿共162张，创作于星期二2023/4/369/162单元模型构造整体坐标系法n4节点四面体单元第69页，讲稿共162张，创作于星期二2023/4/370/162单元模型构造整体坐标系法1.假设插值多项式2.插值多项式为第70页，讲稿共162张，创作于星期二2023/

20、4/371/162单元模型构造整体坐标系法（i=1,2,3,4）循环轮换脚标1、2、3、4，相应可以得到a2，b2，c2，d2、a3，b3，c3，d3、a4，b4，c4，d4 第71页，讲稿共162张，创作于星期二2023/4/372/162单元模型构造整体坐标系法3.单元插值多项式写成矩阵形式（常用）第72页，讲稿共162张，创作于星期二2023/4/373/162单元模型构造整体坐标系法4.单元插值多项式另一种矩阵形式（不常用）第73页，讲稿共162张，创作于星期二2023/4/374/162单元模型构造整体坐标系法n从理论上讲，整体坐标系法可以求任意单元的形函数，但计算过程太复杂n只能求

21、一维2节点线单元、二维3节点三角形单元和三维4节点四面体单元3种简单单元的形函数n复杂的或二次以上的单元必须采用局部坐标系法求n位移场 u 是形函数 Ni 的线性组合，因此形函数Ni同样具有插值多项式的特性第74页，讲稿共162张，创作于星期二2023/4/375/162单元刚度矩阵2节点线单元n一维2节点线单元n单元插值关系 n单元几何关系 n单元本构关系 N=N1 N2 De=E第75页，讲稿共162张，创作于星期二2023/4/376/162单元刚度矩阵2节点线单元n单元刚度矩阵A为单元截面积；l为单元长度n矩阵B第76页，讲稿共162张，创作于星期二2023/4/377/162单元刚度

22、矩阵三角形单元n二维3角形单元n单元插值关系 第77页，讲稿共162张，创作于星期二2023/4/378/162单元刚度矩阵三角形单元n单元几何关系 第78页，讲稿共162张，创作于星期二2023/4/379/162单元刚度矩阵三角形单元n单元本构关系 平面应力问题第79页，讲稿共162张，创作于星期二2023/4/380/162单元刚度矩阵三角形单元n矩阵B第80页，讲稿共162张，创作于星期二2023/4/381/162单元刚度矩阵三角形单元n单元刚度矩阵h为单元厚度k为对称的6*6常数矩阵A为单元面积第81页，讲稿共162张，创作于星期二2023/4/382/162作业n求4节点四面体单

23、元的单元刚度矩阵第82页，讲稿共162张，创作于星期二2023/4/383/162单元模型构造整体坐标系法n单元形函数的特性正规性：单元形函数之和等于1。正交性：形函数在本节点的值等于1，在其它节点的值等于0。例如:2节点线单元形函数第83页，讲稿共162张，创作于星期二2023/4/384/162单元模型等参单元n等参单元 单元内任意一点的位移u与单元节点位移ue之间的关系为 一般单元坐标的插值关系也采用与位移插值关系相同的变换关系即单元内任意一点的坐标x与单元节点坐标xe之间的关系为 第84页，讲稿共162张，创作于星期二2023/4/385/162单元模型等参单元n等参单元凡是几何形状和

24、位移场采用同阶同参数插值关系来描述的单元，称为等参单元 前面介绍的所有单元都属于等参单元 在描述单元的几何形状和位移场时，并不一定非采用同阶插值关系 第85页，讲稿共162张，创作于星期二2023/4/386/162单元模型等参单元n等参单元3节点三角形等参单元 第86页，讲稿共162张，创作于星期二2023/4/387/162单元模型等参单元n超参单元如果几何形状插值函数的阶数高于位移场插值函数的阶数，称为超参单元 n次参单元如果几何形状插值函数的阶数低于位移场插值函数的阶数，称为次参单元 第87页，讲稿共162张，创作于星期二2023/4/388/162单元平衡方程组装过程 n为什么要组装

25、?消除内力n组装的原则是什么?单元自由度与结构自由度对应第88页，讲稿共162张，创作于星期二2023/4/389/162单元平衡方程组装过程 2 F 1 3 U3U4U2U1U5U6结构自由度向量U第89页，讲稿共162张，创作于星期二2023/4/390/162单元平衡方程组装过程 2 1 U3U4U2U11u1u2u3u4 3 U611 U2U12u1u2U5u3u42第90页，讲稿共162张，创作于星期二2023/4/391/162单元平衡方程组装过程2 1 U3U4U2U11u1u2u3u42第91页，讲稿共162张，创作于星期二2023/4/392/162单元平衡方程组装过程组装单

26、元第92页，讲稿共162张，创作于星期二2023/4/393/162单元平衡方程组装过程 3 U611 U2U12u1u2U5u3u4第93页，讲稿共162张，创作于星期二2023/4/394/162单元平衡方程组装过程再组装单元总体刚度方程 K 称为总体刚度矩阵 U 称为位移向量 F 称为载荷向量 第94页，讲稿共162张，创作于星期二2023/4/395/162总体刚度矩阵K的特性 n对称性 n奇异性 n稀疏性 n非零元素带状分布 第95页，讲稿共162张，创作于星期二2023/4/396/162约束处理过程 n为什么要约束处理?总体平衡方程组是奇异的消除无限制的刚体运动 使总体平衡方程组

27、存在唯一一组解第96页，讲稿共162张，创作于星期二2023/4/397/162约束处理过程边界条件n 边界条件分类 力(载荷)边界条件位移边界条件 集中载荷力 表面分布力 自重力热交换引起的温度载荷 固定位移约束 强制位移约束 关联位移约束 第97页，讲稿共162张，创作于星期二2023/4/398/162约束处理过程模型简化xy第98页，讲稿共162张，创作于星期二2023/4/399/162约束处理过程模型简化yx xy第99页，讲稿共162张，创作于星期二2023/4/3100/162约束处理过程约束方程123456789101112yx第100页，讲稿共162张，创作于星期二2023

28、/4/3101/162约束处理过程约束处理方法n位移约束处理方法 赋0赋1法 乘大数法 第101页，讲稿共162张，创作于星期二2023/4/3102/162约束处理过程赋0赋1法n强制位移约束条件处理 U4=C第102页，讲稿共162张，创作于星期二2023/4/3103/162约束处理过程赋0赋1法n强制位移约束条件处理 U4=C第103页，讲稿共162张，创作于星期二2023/4/3104/162约束处理过程赋0赋1法n有6个方程，5个未知数，如果约束方程可以消除有限元平衡方程组的奇异性，则取任意5个方程联立求解，都会得到方程组的唯一一组解。n系数矩阵由原来的对称的变成了非对称的，这对于

29、大规模有限元方程组求解是十分不利的，采用相同的求解方法，在求解时间和矩阵存贮容量方面都增加了一倍。第104页，讲稿共162张，创作于星期二2023/4/3105/162约束处理过程赋0赋1法为了保证系数矩阵的对称性，去掉方程组第4行第105页，讲稿共162张，创作于星期二2023/4/3106/162约束处理过程赋0赋1法引入强制位移约束方程 U4=C，使方程组求解时直接将自由度U4求出第106页，讲稿共162张，创作于星期二2023/4/3107/162约束处理过程赋0赋1法n固定位移约束条件处理 U4=0第107页，讲稿共162张，创作于星期二2023/4/3108/162约束处理过程赋0

30、赋1法n基本原理 利用初等变换对求解方程组进行相同的行列变换，既保证方程组解不会改变，又可以保持方程组系数矩阵的对称性。在进行初等变换时，只要保证对方程组系数矩阵做相同的行列变换，就可以保持方程组系数矩阵的对称性。第108页，讲稿共162张，创作于星期二2023/4/3109/162约束处理过程乘大数法n乘大数法n基本原理 利用矩阵的初等变换不改变方程组解的思想。第109页，讲稿共162张，创作于星期二2023/4/3110/162约束处理过程乘大数法n强制位移边界条件 第110页，讲稿共162张，创作于星期二2023/4/3111/162约束处理过程乘大数法强制约束方程 A是一个大数，是系数

31、矩阵中对角线元素K44的1010倍量级以上 为什么要乘以大数A？放大位移约束方程的优势第111页，讲稿共162张，创作于星期二2023/4/3112/162约束处理过程乘大数法n强制位移边界条件 第112页，讲稿共162张，创作于星期二2023/4/3113/162约束处理过程乘大数法n固定位移边界条件 C=0 约束后的方程组简化为 第113页，讲稿共162张，创作于星期二2023/4/3114/162约束处理过程两种方法比较n赋0赋1法在约束处理过程中是严格精确的，而乘大数法是一种近似约束处理方法，它的精度取决于所乘大数A值 两种方法都可以消除有限元平衡方程的奇异性，得到符合实际边界条件的唯

32、一一组解。但两种方法还是有很大的区别 第114页，讲稿共162张，创作于星期二2023/4/3115/162约束处理过程两种方法比较n采用乘大数法约束处理后的有限元平衡方程在求解时可能造成解的失真，大数A值越大可能解的偏差会越大，而赋0赋1法就不会出现类似的问题，它在约束过程和求解过程都是精确的n乘大数法相对于赋0赋1法在约束处理过程上简单一些 第115页，讲稿共162张，创作于星期二2023/4/3116/162约束处理过程两种方法比较n赋0赋1法实际上是将关联位移约束方程代入到有限元平衡方程中的，是代入法。而乘大数是将占绝对优势的关联位移约束方程合并到有限元平衡方程中的，是罚方法，计算误差

33、来自于合并过程，计算精度取决于关联位移约束方程的优势大小n商业软件中，位移边界条件的约束处理都采用赋0赋1法，乘大数很少被采用主要原因是它是一种近似方法，而且大数的大小也不好确定，有时还会造成求解失败 第116页，讲稿共162张，创作于星期二2023/4/3117/162约束处理过程弹簧单元假设柔性弹簧 kOXYU4f f=kU4k第117页，讲稿共162张，创作于星期二2023/4/3118/162弹簧约束方程 f=kU4约束处理过程弹簧单元第118页，讲稿共162张，创作于星期二2023/4/3119/162方程组求解过程特点n方程组求解有限元计算过程中很重要的一部分，在有限元法的发展过程

34、中，有限元方程的求解效率一直是其应用的最大瓶颈之一 n有限元方程组的特点：有限元方程组的系数矩阵具有对称、稀疏、带状分布以及正定、主元占优。有效地利用这些特点，以减少系数矩阵的存贮量，提高方程组求解效率 第119页，讲稿共162张，创作于星期二2023/4/3120/162方程组求解过程分类比较n线性方程组的解法主要分两大类:直接解法：以高斯消去法基础，以等带宽或变带宽方式存贮系数矩阵内元素，对于求解规模比较大的问题，要存贮的元素非常巨大。迭代解法：只需要存贮系数矩阵中非零元素，存贮量很小，一般是变带宽存贮量的20%或更少，有些算法的求解效率也非常高，适合求解大规模线性方程组。但是这种解法对接

35、近病态的方程组很难保证收敛性。第120页，讲稿共162张，创作于星期二2023/4/3121/162方程组求解过程带宽定义n有限元方程组系数矩阵是稀疏的、非零元素呈带状分布，带宽就是它的宽度，带宽的大小是由系统有限元网格的节点号排序决定的，具体求法是 带宽=(单元最大节点号之差+1)*节点自由度数 n带宽是网格节点标注方法直接决定的，不同标注方法带宽可能相关很大 第121页，讲稿共162张，创作于星期二2023/4/3122/162方程组求解过程带宽n带宽是网格节点标注方法直接决定的，不同标注方法带宽可能相关很大 第122页，讲稿共162张，创作于星期二2023/4/3123/162方程组求解

36、过程带宽n所示四边形网格的三种节点号标注方法，每个节点是2个自由度n结构的带宽分别是12，18，56，相差很大，其中12和56之间相差近5倍，这就意味着系数矩阵的存贮量也是相差5倍，因此，对于大规模复杂系统的节点号优化是十分必要的 第123页，讲稿共162张，创作于星期二2023/4/3124/162方程组求解过程系数矩阵存贮 n系数矩阵存贮 如果节点号排序优化的比较好，系数矩阵的存贮量就会减少很多。根据系数矩阵的对称性，一般都是按半带宽存贮。n系数矩阵存贮的方法 二维等带宽存贮 一维变带宽存贮第124页，讲稿共162张，创作于星期二2023/4/3125/162方程组求解过程二维等带宽存贮

37、n二维等带宽存贮 第125页，讲稿共162张，创作于星期二2023/4/3126/162方程组求解过程二维等带宽存贮n二维等带宽存贮消除了最大带宽以外的全部零元素，节省了系数矩阵元素的存贮量。但是由于取最大带宽为存贮范围，因此不能排除在带宽内的大量零元素。当系数矩阵的各行带宽变化不大时，适合采用二维等带宽存贮，方程组求解过程中系数矩阵元素的寻址也比较方便，求解效率较高。n当出现局部带宽特别大的情况时，采用二维等带宽存贮时，将由于局部带宽过大而使整体系数矩阵的存贮大大增加。第126页，讲稿共162张，创作于星期二2023/4/3127/162方程组求解过程一维变带宽存贮 n一维变带宽存贮 一维变

38、带宽存贮方法就是把变化的带宽内的元素按一定的顺序存贮在一个一维数组中。由于它不按最大带宽存贮，因此比二维等带宽存贮更节省内存。按照解法可分为按行一维变带宽存贮和按列一维变带宽存贮。第127页，讲稿共162张，创作于星期二2023/4/3128/162n按行一维变带宽存贮 方程组求解过程一维变带宽存贮 辅助的寻址数组M 第128页，讲稿共162张，创作于星期二2023/4/3129/162n一维变带宽存贮是最节省内存的一种方法，但是由于要借助于寻址数组寻找系数矩阵元素的位置，相对二维等带宽存贮方法来说要复杂一些，而且在程序实现时也要复杂得多，方程组求解过程中也要消耗一些数组寻址时间。因此，在选用

39、存贮方法时要权衡二者的利弊，统盘考虑。一般当带宽变化不大，计算机内存允许时，采用二维等带宽存贮方法是比较合适的。方程组求解过程一维变带宽存贮 第129页，讲稿共162张，创作于星期二2023/4/3130/162方程组求解过程求解方法n方程组求解方法 高斯消去法 三角分解法 雅可比（Jacobi）迭代法 高斯-赛德尔（Gauss-Seidel）迭代法 第130页，讲稿共162张，创作于星期二2023/4/3131/162应变、应力回代过程 n单元应变和应力回代求解 通过求解有限元平衡方程得到有限元节点位移后，就可以进行系统的刚度校核。如果所分析问题要进行强度校核，就要回代求解单元的应变和应力。

40、n由插值关系和几何关系可得单元应变，再通过本构关系得到单元应力第131页，讲稿共162张，创作于星期二2023/4/3132/162有限差分法有限差分法第132页，讲稿共162张，创作于星期二2023/4/3133/162 从弹性力学的基本方程建立以来，这些方程在各种问题的边界条件下如何求解，一直是很多数学工作者和力学工作者研究的内容。即弹性力学的经典解法存在一定的局限性，当弹性体的边界条件和受载情况复杂一点，往往无法求得偏微分方程的边值问题的解析解，许多工程重要问题，不能够得出函数式的解答。因此，弹性力学问题的各种数值解法便具有重要的实际意义。第133页，讲稿共162张，创作于星期二2023

41、/4/3134/162 工程中常用的数值解法有有限单元法和差分法。有限单元法 是以有限个单元的集合体来代替连续体，属于物理上的近似。差分法 是把弹性力学的基本方程和边界条件（一般均为微分方程）近似地改用差分方程（代数方程）来表示，把求解微分方程的问题改换成为求解代数方程的问题，属于数学上的近似。第134页，讲稿共162张，创作于星期二2023/4/3135/162 第一节 差分方程 第二节 应力函数的差分解 第三节 深梁应力函数的差分解第135页，讲稿共162张，创作于星期二2023/4/3136/162第一节 差分方程 差分法是沿用已久的一种数值解法。随着计算机的普及和相应的软件发展，此法成

42、为解弹性力学问题的一种有效的方法。第136页，讲稿共162张，创作于星期二2023/4/3137/162 我们在弹性体上，用相隔等间距h而平行于坐标轴的两组平行线织成正方形网格，x=y=h，如图。设f=f(x,y)为弹性体内的某一个连续函数。该函数在平行于x轴的一根网线上，如在-上,它只随x坐标的改变而变化。在邻近结点处，函数f可展为泰勒级数如下：第137页，讲稿共162张，创作于星期二2023/4/3138/162 我们将只考虑离开结点充分近的那些结点，即（x-x0）充分小。于是可不计（x-x0）的三次及更高次幂的各项，则上式简写为：在结点，x=x0-h,在结点1,x=x0+h，代入(b)得

43、：第138页，讲稿共162张，创作于星期二联立(c),(d),解得差分公式：同理，在网线-上可得到差分公式第139页，讲稿共162张，创作于星期二 差分公式(-)及(-)是以相隔2h的两结点处的函数值来表示中间结点处的一阶导数值，可称为中点导数公式。以相邻三结点处的函数值来表示一个端点处的一阶导数值，可称为端点导数公式。应当指出：中点导数公式与端点导数公式相比，精度较高。因为前者反映了结点两边的函数变化，而后者却只反映了结点一边的函数变化。因此，我们总是尽可能应用前者，而只有在无法应用前者时才不得不应用后者。第140页，讲稿共162张，创作于星期二 以上(-)(-)是基本差分公式,从而可导出其

44、它的差分公式如下：第141页，讲稿共162张，创作于星期二2023/4/3142/162第二节 应力函数的差分解 当不计体力时，我们已把弹性力学平面问题归结为在给定边界条件下求解双调和方程的问题。用差分法解平面问题，就应先将双调和方程变换为差分方程，而后求解之。第142页，讲稿共162张，创作于星期二 一旦求得弹性体全部节点的值后，就可按应力分量差分公式（对节点0）算得弹性体各节点的应力。第143页，讲稿共162张，创作于星期二 可见，用差分法解平面问题，共有两大任务：一、建立差分方程 将（1-68）代入双调和方程 对于弹性体边界以内的每一结点，都可以建立这样一个差分方程。整理即得第144页，

45、讲稿共162张，创作于星期二2023/4/3145/162 二、联立求解这些线性代数方程，就能求得各内结点处的值。为了求得边界上各结点处的值，须要应用应力边界条件，即：一般建立和求解差分方程，在数学上不会遇到很大困难。但是，当对于边界内一行的（距边界为h的）结点，建立的差分方程还将涉及边界上各结点处的值，并包含边界外一行的虚结点处的值。第145页，讲稿共162张，创作于星期二代入上式，即得：l1=cos(N,x)=cos=dy/ds,l2=cos(N,y)=sin=-dx/ds,于是，式（a）可改写为：由右图可见，第146页，讲稿共162张，创作于星期二关于边界上任一点处由此得：的值，可将（b

46、）式从A点到B点对s积分得到：第147页，讲稿共162张，创作于星期二 将此式亦从A点到B点沿s进行积分，就得到边界上任一点B处的值。为此利用分部积分法，得：由高等数学可知，第148页，讲稿共162张，创作于星期二将式(b),(c)代入，整理得：由前知，把应力函数加上一个线性函数，并不影响应力。因此，可设想把应力函数加上a+bx+cy，然后调整a，b，c三个数值，使得由式(d)及式(c)可见，设即可根据面力分量及导数求得为已知，第149页，讲稿共162张，创作于星期二 从图易看出，式（23）右边的积分式表示A与B之间的x方向的面力之和；式（24）右边的积分式表示A与B之间的y方向的面力之和；式

47、（25）右边的积分式表示A与B之间的面力对于B点的矩。于是式(d)，式(c)即简化为：第150页，讲稿共162张，创作于星期二 至此，我们解决了怎样计算边界上各结点的值的问题。至于边界外一行虚结点处 的值，则可用边界上结点处的 或 值和边界内一行相应结点处 的值来表示。例如，对于图1中的虚结点14，因为有 第151页，讲稿共162张，创作于星期二2023/4/3152/162所以有 当求出全部结点上的值以后，我们就可按应力分量的差分公式（21）计算应力分量。用差分法解弹性平面问题时，可按下列步骤进行：（1）在边界上任意选定一个结点作为基点A，取第152页，讲稿共162张，创作于星期二 然后由面

48、力的矩及面力之和算出边界上所有各结点处 的值，以及所必需的一些 及 值，即垂直于边界方向的导数值。（2）应用公式（26），将边界外一行虚结点处的值用边界内的相应结点处的值来表示。（3）对边界内的各结点建立差分方程（22），联立求解这些结点处的值。（4）按照公式（26），算出边界外一行的各虚结点处的值。第153页，讲稿共162张，创作于星期二2023/4/3154/162（5）按照公式（21）计算应力的分量。说明：1以上是针对单连体导出的结果。对于多连体，情况就不象这样简单。2.如果一部分边界是曲线的，或是不与坐标轴正交，则边界附近将出现不规则的内结点。对于这样的结点，差分方程（22）必须加以修

49、正。第154页，讲稿共162张，创作于星期二2023/4/3155/162第三节 例 深梁的应力函数差分解 现以如图所示的混凝土深梁为例，应用应力函数的差分解求出应力分量。已知混凝土深梁上边受有均布向下的铅直荷载q，并由下角点处的反力维持平衡。第155页，讲稿共162张，创作于星期二(1）计算边界上各结点的、和 值。取A为基点，且由上面公式所得的计算结果见下表。第156页，讲稿共162张，创作于星期二（2）计算边界外以行各虚结点处的值。由式（2-6）及前表可得第157页，讲稿共162张，创作于星期二(3）边界内各结点的差分方程,由式（2-2）可知第158页，讲稿共162张，创作于星期二联立求解

50、上式，可得（以qh2为单位)（4）计算结点外一行各结点处的值。由(a)、(b)、(c)可得（5）计算应力。对于结点M，由式（1）可知 第159页，讲稿共162张，创作于星期二2023/4/3160/162同理可得 沿着梁的中线MA，的变化如下图和右图所示。第160页，讲稿共162张，创作于星期二差分方法的优缺点n优点：1.差分法是解微分方程边值问题和弹性力学问题的有效方法；2.差分法简便易行；3.对于某些结构，为了更精确地分析局部的应力状态，可以用差分法进行分析；n缺点：1.对于曲线边界和斜边界等产生的不等间距网格的处理，比较麻烦和易于出错；2.比较比较适用于求解二维问题或平面问题；3.比较适