塑性加工过程CAE-02FEM基础.ppt

塑性加工过程CAE 板料冲压成形CAE分析,鲍 益 东 南京航空航天大学 机电学院 航空宇航制造工程系,第三章 板料成形的有限元基础,3.1 序言 3.2 杆和梁单元 3.2 二维问题 3.3 有限元求解技术 3.4 板壳单元,3.1 序言,3.1.1 基本概念 3.1.2 矩阵代数的回顾 3.1.3 弹簧元,3.1 序言,有限元方法或有限元分析的基本思想就是将一个复杂的物体分解成许多易于处理的小片来处理。这种思路在日常生活和工程实际中经常被使用。,3.1.1 基本概念,搭积木游戏,建筑,3.1.1 基本概念,圆面积的近似求解：,注意：一个复杂的连续的物体可以被许多小片（单元）来近似代替。,工程中有限元的应用 机械、航空、土木、汽车工程 结构分析（静态/动态、线性/非线性） 热/液体流动 电磁 地质力学 生物力学 ,齿轮式弹性轴接的模型,3.1.1 基本概念,有限元的简要发展历史 1943 Courant 1956 Tuner, Clough, Martin, Topp (Stiffness) 1960 Clough (“Finite Element” 平面问题) 1970s 在大型机上得到应用 1980s 在微机上得到应用，前后处理软件 1990s 大型机构系统分析,3.1.1 基本概念,3.1.1 基本概念,易拉罐的跌落测试,结构分析中的有限元（过程） 将结构体划分成小片（单元，节点） 形成描绘物理特性的单元刚度 将单元组装成一个整体结构的近似方程组 求解已经引入位置物理量（位移）的方程组 计算单元中用户所关心的物理量（应变，应力）,3.1.1 基本概念,计算机 前处理（建立有限元模型，载荷和约束） 有限元求解器（组装和求解系统方程） 后处理（显示计算结果） 商业化有限元软件包 ANSYS，NASTRAN，ALGOR(通用目的) ABAQUS(非线性动力分析) PATRAN，HyperMesh(前后处理) LS-DYNA(碰撞动力分析) DynaForm(前后处理) ,3.1.1 基本概念,本章有限元基础课的目的 理解有限元分析的基本原理和思路 掌握本章中所涉及到的单元模型的推导和适用范围 对一个给定的问题能够建立适当的有限元模型 能够解释并评价有限元分析结果的优劣（知道问题的物理意思） 明白有限元的局限性（不要错误地应用有限元数值工具）,3.1.1 基本概念,线性代数方程组系统： 矩阵形式： A称为nn矩阵，x和b分别是n维的列向量,3.1.1 矩阵代数的回顾,矩阵的加法减法： 矩阵的乘法： 矩阵的转置： 对称矩阵： 单位矩阵：,3.1.1 矩阵代数的回顾,矩阵行列式的值： 奇异矩阵： 如果 det A =0，那么系统存在问题（非唯一解，发散等） 矩阵的逆：,3.1.1 矩阵代数的回顾,线性方程组系统的求解技术： 高斯消去法 迭代法 正定矩阵： 对于所有非零向量X，有XTAX0，A为正定矩阵 正定矩阵为非奇异矩阵 矩阵的导数和积分：,3.1.1 矩阵代数的回顾,弹簧元：,3.1.3 弹簧单元,考虑弹簧元的力的平衡条件： 节点i： 节点j： 矩阵形式： 注意：K为对称矩阵，K是非奇异的还是奇异的？,3.1.3 弹簧单元,弹簧元系统： 单元1： 单元2：,3.1.3 弹簧单元,对整个系统进行单元刚度矩阵的组装： 考虑节点力的平衡条件： 节点1： 节点2： 节点3： 矩阵形式： K为该弹簧系统的刚度矩阵（结构矩阵）,3.1.3 弹簧单元,单元刚度矩阵另一种组装方法： 分别扩大单元1和2的刚度矩阵 这与根据节点力平衡得出的矩阵是一样的。,3.1.3 弹簧单元,引入边界条件和力的条件： 假设u1=0，F2=F3=P，那么 未知量为： u2，u3和F1 求解方程可得：,3.1.3 弹簧单元,检查计算结果： 结构变形后的形状 外力平衡 有关弹簧元的注意事项： 适于刚度分析的计算 不适合用于弹簧本身的应力分析计算 在弹簧元的横向是否具有刚度，弹簧元是否具有扭转刚度,3.1.3 弹簧单元,例子1： 已知： 求： (a) 整体刚度矩阵 (b) 节点2和3的位移 (c) 节点1和4的支反力 (d) 弹簧2的力,3.1.3 弹簧单元,（a）问：分别求出单元刚度矩阵： 单元1： 单元2： 单元3： 组装后：,3.1.3 弹簧单元,组装后： 注意到整体刚度矩阵是对称并带状分布的。 该系统的平衡方程为：,3.1.3 弹簧单元,（b）问： 将边界条件（u1=u4=0）应用到平衡方程中，去掉1行1列，4行4列后： 求解得： （c）问： 从平衡方程组中的1和4可得：,3.1.3 弹簧单元,（d）问： 弹簧元2的平衡方程为： 式中i=2，j=3，可以计算出弹簧力为：,3.1.3 弹簧单元,例子2： 问题描述：对上述一个具有任意弹簧元节点和单元的系统，求其整体刚度矩阵。,3.1.3 弹簧单元,单元拓扑关系： 上表中为每个弹簧元的局部节点号和整体节点号的对应关系。 然后依次求出每个弹簧元的刚度矩阵：,3.1.3 弹簧单元,组装后： 整体刚度矩阵是对称并带状分布的。,3.1.3 弹簧单元,3.2 杆和梁单元,3.2.1 线性静力分析 3.2.2 杆单元 3.2.3 梁单元,大部分结构分析问题都可以看作是线性静力分析问题，它们都基于以下假设： 1. 小变形（加载方式不会因为变形而改变） 2. 线弹性材料（不存在塑性和破裂） 3. 静力载荷（在结构上的载荷是慢速平稳地施加上去的） 线性静力分析可以解决结构分析问题中大部分的问题，对于大多数的结构分析问题，静力分析可以得到一个近似的结果。 线性静力分析是非线性分析的基础。,3.2.1 线性静力分析,杆单元：,3.2.2 杆单元,刚度矩阵直接法： 假设位移u沿着杆的轴向线性分布： K即为杆单元的刚度系数，杆单元和弹簧元类似。,3.2.2 杆单元,单元刚度矩阵为： 单元的平衡方程组为： 节点自由度： 对一维的杆单元，每个节点就只有一个自由度 刚度矩阵K中系数的物理意义： K中第j列的系数表示在节点j施加单位位移而其他节点固定不动的时候，杆上所承受的力。,3.2.2 杆单元,刚度矩阵正规推导方法： 定义两个形函数： 位移u沿着杆的轴向线性分布： 可得应变为： B为单元的应变-位移矩阵：,3.2.2 杆单元,单元应力为： 杆单元上的应变能为： 杆单元2个节点所做的功为： 对于保守系统，U = W ，故有：,3.2.2 杆单元,上式等价于： 其中k即为单元刚度矩阵 上述方法就是正规的推导过程，该方法也可以用来推导其他类型的单元刚度矩阵。 单元刚度矩阵也可以通过其他严格的方法获得，比如最小势能原理，伽辽金方法等。 于是我们可以得到杆单元的刚度矩阵： 上式结果和直接法的结果一样。,3.2.2 杆单元,例子： 问题描述：在节点2处施加F，求杆1和杆2上的应力。 求解方法：用1-D杆单元。 单元1和单元2的刚度矩阵分别为： 假设在节点2处是用一个无摩擦的铰链将杆1和杆2所连接。,3.2.2 杆单元,组装后： 载荷和边界条件为： 删除第1行第1列和第3行第3列后得：,3.2.2 杆单元,最后可以得到单元1内的应力： 同理，可以得到单元2内的应力： 负号表示单元2所受的是压应力。,3.2.2 杆单元,注意： 1 ) 在这个例子中，根据一维线性理论所计算出来的单元1和2中的应力是精确解，所以如果我们将单元细化不会提高精度。 2 ) 如果对于阶梯杆结构，其中的A要采用横截面的平均面积。 3 ) 为了得到单元1和2中的应力，我们首先要得到节点的位移，因为采用的基于位移场的有限元法。,3.2.2 杆单元,二维空间上的杆单元： 注意：在线弹性理论的前提下，横向位移 对杆单元的拉伸没有贡献。,3.2.2 杆单元,坐标变换： 矩阵形式：,3.2.2 杆单元,对于2个节点的杆单元有： 节点力也采用相同的变换方法： 二维空间下的局部坐标系下的刚度矩阵：,3.2.2 杆单元,扩大后： 矩阵形式： K为整体坐标系下的刚度矩阵,3.2.2 杆单元,K的显式表达式为： 单元应力：,3.2.2 杆单元,三维空间上的杆单元： 先在局部坐标系下求出单元刚度矩阵，然后再转换并组装到整体坐标系下进行计算。,3.2.2 杆单元,简单的平面梁单元： 基本的梁理论：,3.2.3 梁单元,3.3 二维问题,3.3.1 基本理论回顾 3.3.2 二维问题的有限元,应力和应变： 在特定的条件下，应力应变状态可以简化，因此，一般的三维问题可以简化成二维问题来分析,3.3.1 基本理论回顾,平面应力： 平面应变：,3.3.1 基本理论回顾,本构关系（平面应力）： 本构关系（平面应变）：,3.3.1 基本理论回顾,应变与位移的关系： 边界条件：,3.3.1 基本理论回顾,常应变单元（CST or T3） 位移函数：,3.3.2 二维问题的有限元,位移插值函数：,3.3.2 二维问题的有限元,应变： CST单元刚度矩阵为： t为单元厚度，k为一个66的对称矩阵。,3.3.2 二维问题的有限元,双线性四边形单元（Q4） 位移函数：,3.3.2 二维问题的有限元,例子：,3.3.2 二维问题的有限元,网格划分：,3.3.2 二维问题的有限元,计算结果：,3.3.2 二维问题的有限元,载荷的变换： 集中力(点载荷)，面力(压力载荷)，体力(重力)为三种主要的外载荷。 面力和体力都需要经过变换以后才能加到有限元中。 变换的基本思想就是：等效功。,3.3.2 二维问题的有限元,其中，t为单元厚度，L为单元边长，un为垂直于边界AB的位移分量。 对于Q4单元来说，有：,3.3.2 二维问题的有限元,应力计算： 单元内的应力计算公式为： B为应变与节点位移之间的关系矩阵，d为计算后的节点位移向量。 我们可以计算出单元内部任意一点的应力，包括单元中心和节点位置。 应力分布的等值线可以通过后处理软件显示出来。,3.3.2 二维问题的有限元,讨论： 1）需要了解每个单元类型的特性： T3和Q4：线性位移，常量应变和应力 T6和Q8：二次位移，线性的应变和应力 2）对于一个特定的问题选择一种合适的单元 毫无疑问，应该尽可能地采用高阶单元和细网格。 3）要避免单元具有大形状比和尖角：,3.3.2 二维问题的有限元,讨论： 4）单元之间要互相连通： 在有限元模型中不能存在间隙或自由单元。,3.3.2 二维问题的有限元,3.3 有限元求解技术,3.3.1 方程组求解 3.3.2 有限元方法的实质 3.3.3 数值误差 3.3.4 有限元求解的收敛性,直接法（高斯消去法）： 求解时间是和NB2成正比（N为矩阵的维数，B为矩阵的带宽） 适合小到中等规模的问题，小带宽的问题。 对于多重载荷的情况，容易处理。 迭代法： 求解时间事先不可预知 降低对存储空间的要求 适合大型规模问题的求解（带宽可以很大，收敛快） 对于不同载荷的工况需要重新求解。,3.3.1 方程组求解,高斯消去法例子： 回代求解：,3.3.1 方程组求解,迭代法（高斯-赛德尔）例子： 从一个初始解X(0)开始按照以下公式计算： 一直到X满足以下条件，迭代终止。 其中为收敛控制的容忍误差。,3.3.1 方程组求解,有限元模型是基于很多近似以后对实际结构的一个数学模型。 实际结构具有无限个节点，所以具有无限个自由度。 有限元模型具有有限个节点，所以具有有限个自由度。 位移场是被有限个节点位移的值所控制的：,3.3.2 有限元方法的实质,刚度效果： 有限元模型比实际结构要更刚一些。 一般来说，位移解要比实际情况偏小一些。 所以，有限元计算的位移是精确解的一个下限：,3.3.2 有限元方法的实质,误差错误（模型和计算）。 误差的类型： 模型误差（梁，板，理论） 离散误差（有限，分段， ） 数值误差（在求解有限元方程组的时候） 例子（数值误差）：,3.3.3 数值误差,有限元方程组为： 如果K2 K1的话，方程组是奇异的，病态的：,3.3.3 数值误差,如果K2 K1的话，方程组是非病态的： 在有限元模型中不同部分的刚度相差很大的话有可能导致有限元方程组的病态，产生很大的误差。 病态方程组的情况下，很小的输入变化量（右端项向量）求解时会引起很大的变化量。,3.3.3 数值误差,随着有限元模型中的网格不断地被加密和细化，有限元的计算结果将收敛于被求解的问题的数学模型的精确解。 有限元网格自适应类型： h-refinement: 减小单元尺寸（h指单元的尺寸）。 p-refinement: 增加单元多项式插值函数的阶次（p指更高阶次的多项式插值函数）。 r-refinement: 重新划分有限元网格。 hp-refinement: h方法和p方法互相结合的方法。,3.3.4 有限元求解的收敛性,网格自适应类型例子：,3.3.4 有限元求解的收敛性,网格自适应类型例子：,3.3.4 有限元求解的收敛性,网格自适应类型例子：,3.3.4 有限元求解的收敛性,3.4 板壳单元,3.4.1 板理论 3.4.2 板单元 3.4.3 壳和壳单元,板理论： 平板 横向加载 弯曲效应为主 注意： 应用范围： 剪力墙 地板 架子等,3.4.1 板理论,板单元的力和弯矩： 应力：,3.4.1 板理论,板单元的力和应力之间的关系： 单位长度的弯矩： 单位长度的扭矩： 单位长度的剪力： 最大弯曲应力： 最大弯曲应力总是在 在中性层没有弯曲应力（与梁单元类似）,3.4.1 板理论,薄板理论（Kichhoff 板理论）： 直法线假设（与梁单元类似）： 与中性层垂直的直法线在变形后仍然与中性层保持垂直，也就是说没有剪切变形： 位移：,3.4.1 板理论,应变： 注意：在中性层没有拉伸变形。 应力（平面应力）： 剪力和弯矩：,3.4.1 板理论,厚板理论（Mindlin 板理论）： 如果板厚不够薄，即 t/L 1/10，应该采用厚板理论。 直线假设：该理论认为在变形过程中截面角度发生改变，即： 与中性层垂直的直线在变形后与中性层不再保持垂直，但仍然是直线。,3.4.1 板理论,新的独立变量： x和y ：分别为直线沿x轴和y轴的旋转角度，变形前该直线与中性层垂直。 新的关系： 注意：如果我们引入以下条件，厚板理论又可以退化成薄板理论。,3.4.1 板理论,Kichhoff 板单元： 四节点四边形单元 每个节点的自由度： 每个单元的z方向的位移为：,3.4.2 板单元,Kichhoff 板单元是一种不协调单元，但是其刚度矩阵仍旧可以写成： B为