第九届高级非线性培训一维和有限元演示文稿

第九届高级非线性培训一维和有限元演示文稿当前第1页\共有68页\编于星期六\16点优选第九届高级非线性培训一维和有限元当前第2页\共有68页\编于星期六\16点1引言

当前第3页\共有68页\编于星期六\16点1引言

非线性连续体一维模型（杆）的有限元方程在固体力学中，Lagrangian网格是最普遍应用的，其吸引力在于它们能够很容易地处理复杂的边界条件，并且能够跟踪材料点，因此能够精确地描述依赖于历史的材料。在Lagrangian有限元的发展中，一般采用两种方法：以Lagrangian度量的形式表述应力和应变的公式，导数和积分运算采用相应的Lagrangian（材料）坐标X，称为TotalLagrangian格式（TL）。2.以Eulerian度量的形式表述应力和应变的公式，导数和积分运算采用相应的Eulerian（空间）坐标x，称为UpdatedLagrangian格式（UL）。非线性与线性公式的主要区别是非线性需要定义积分赋值的坐标系（构形）和确定选择应力和应变的度量。当前第4页\共有68页\编于星期六\16点非线性两种格式的主要区别在于：TL，在初始构形上描述变量，UL，在当前构形上描述变量。不同的应力和变形度量分别应用在这两种格式中。TL，习惯于采用一个应变的完全度量，UL，常常采用应变的率度量。这些并不是格式的固有特点，在UL中采用应变的完全度量是可能的，并且在TL中可以采用应变的率度量。尽管TL和UL表面看来有很大区别，两种格式的力学本质是相同的；因此，TL可以转换为UL，反之亦然。1引言

当前第5页\共有68页\编于星期六\16点

对于每一种公式，将建立动量方程的弱形式，即虚功原理（或虚功）。这种弱形式是通过对变分项与动量方程的乘积进行积分来建立。在TL格式中，积分在X坐标上进行；在Eulerian和UL格式中，积分在x上进行。我们将说明如何处理力边界条件，因此近似（试）解不需要满足力边界条件。这个过程与在线性有限元分析中的过程是一致的，在非线性公式中的主要区别是需要定义积分赋值的坐标系和确定选择应力和应变的度量。1引言

推导有限元近似计算的离散方程。对于考虑加速度（动力学）或那些包含率相关材料的问题，推导离散有限元方程为普通微分方程（ODEs）。因为有限元仅将空间微分运算转化为离散形式，而没有对时间导数进行离散，这个空间的离散过程称为半离散化。对于静力学与率无关材料问题，离散方程独立于时间，有限元离散将导致一组非线性代数方程。当前第6页\共有68页\编于星期六\16点2完全的Lagrangian格式当前第7页\共有68页\编于星期六\16点2.2TL的控制方程初始构形参考构形当前构形变形构形2完全的Lagrangian格式当前第8页\共有68页\编于星期六\16点物体的运动由Lagrangian坐标和时间的函数描述是在初始域与当前域之间的映射

当材料坐标在初始位置

2完全的Lagrangian格式位移差

或者

变形梯度

偏微分的意义？

当前第9页\共有68页\编于星期六\16点2完全的Lagrangian格式定义Jacobian：作为变形物体的无限小体积相对于变形前物体微段体积的比值

应变的度量

在变形前构形中上式为零，它等效于工程应变

应力的度量

Cauchy应力

名义应力

在多维上没有工程应力的定义。

工程应力

物理应力

初始值，J0＝1当前第10页\共有68页\编于星期六\16点2完全的Lagrangian格式推导方程应用下面方程推导非线性杆：

1.质量守恒

2.动量守恒

3.能量守恒

4.变形度量：应变-位移方程

5.本构方程：应力-变形度量的关系

另外，要求变形保持连续性，称为协调性要求。当前第11页\共有68页\编于星期六\16点质量守恒

2完全的Lagrangian格式对于Lagrangian格式，质量守恒方程为对于一维杆动量守恒

由名义应力P和Lagrangian坐标给出(单位长度的力)

如果初始横截面面积A0在空间保持常数，则动量方程成为应力在坐标方向的分量

b－单位质量的力－体力

当前第12页\共有68页\编于星期六\16点平衡方程

2完全的Lagrangian格式平衡意味着物体处于静止或者以匀速运动

能量守恒

内部功率由变形率的梯度和名义应力P的乘积给出

本构方程

不计惯性力，则动量方程成为平衡方程etc.

表示影响应力的其他变量，如温度，夹杂等。

是变形历史的函数。

完全形式

率形式

当前第13页\共有68页\编于星期六\16点2完全的Lagrangian格式本构方程的例子

1)线弹性材料

完全形式

率形式

2)线性粘弹性材料

etc.

表示影响应力的其他变量，如温度，夹杂等。

是变形历史的函数。

完全形式

率形式

当前第14页\共有68页\编于星期六\16点2完全的Lagrangian格式边界条件

位移边界

力边界n0单位法线（＋,－）

一端固定一端自由杆

边界条件满足

初始条件

动量方程是关于X二阶的(偏微分方程)。因此在每一端，必须描述u或者作为边界条件。当前第15页\共有68页\编于星期六\16点内部连续条件

跳跃条件函数的连续性

如果函数的第n阶导数是连续函数，该函数为连续函数是连续可导的(它的一阶导数存在并且处处连续)在函数中，导数只是分段可导，一维函数不连续发生在点上，二维函数不连续发生在线段上，三维函数不连续发生在面上。

函数本身不连续，xi是不连续点。动量平衡要求式中[[f]]表示在f(x)中的跳跃，即如在材料界面处，应变不连续，需要补充该条件，跳跃函数给出内部连续条件当前第16页\共有68页\编于星期六\16点关于泛函和变分的概念变

量函数函数泛函

泛函－函数的函数（functional,functionoffunction)当虚位移是真实位移的增量时，虚位移原理We＝V中的V就是泛函V的变分。微分是函数的增量，变分是泛函的增量。w(x)是x函数V(w(x))是w(x)的泛函当前第17页\共有68页\编于星期六\16点

自然变分原理是对物理问题的微分方程和边界条件建立对应的泛函，使泛函取驻值得到问题的解答，但是其未知场函数需要满足一定的附加条件。

广义变分原理（或称约束变分方程）不需要事先满足附加条件，采用Lagrange乘子法和罚函数法将附加条件引入泛函，重新构造一个修正泛函，将问题转化为求修正泛函的驻值。称为无附加条件的变分原理。(如第10章的接触弱形式）对于罚函数方法，将罚参数取正值，对修正泛函得到的近似解只是近似地满足附加条件，罚参数值越大，附加条件的满足程度就越好。而在实际计算中，罚函数只能取有限值，所以利用罚函数求解只能得到近似解。2完全的Lagrangian格式当前第18页\共有68页\编于星期六\16点

有限元方法不能直接离散动量方程。为了离散这个方程，需要一种弱形式(weakform)，称为变分形式，即虚功原理或者虚功率，通过对变分项与动量方程的乘积进行积分来建立的。虚功原理或者弱形式是等价于动量方程和力边界条件的。后者称为经典强形式。有限元的理论基础是泛函、变分和虚功原理。2完全的Lagrangian格式2.3TL的弱形式强形式到弱形式弱形式到强形式当前第19页\共有68页\编于星期六\16点

对于动量方程和力边界条件,现在建立弱形式，要求：满足所有位移边界条件并足够平滑，因此确切定义了动量方程中的所有导数。也假设足够光滑，这样确切定义了所有的后续步骤，并在指定的位移边界条件上为零。这是标准和经典的建立弱形式的方法。尽管它所导致的连续性要求比在有限元近似中更加严格，在看到以较少的强制连续性要求所得到的结论之前，仍采用这种方法。2完全的Lagrangian格式试函数变分项强形式到弱形式当前第20页\共有68页\编于星期六\16点取动量方程与变分项的乘积并在全域内积分得到弱形式，给出

2完全的Lagrangian格式强形式到弱形式名义应力P是一个试位移函数。展开第一项乘积的导数，整理得到分布积分在指定位移边界处变分项消失，第二行服从边界互补条件和力边界条件。给出完全的Lagrangian格式的动量方程和力边界条件的弱形式

当前第21页\共有68页\编于星期六\16点问：为什么消除关于应力P的导数？答：在这种弱形式中出现了应力的空间导数，由本构关系，应力是连续函数，也应该是位移的导数的平滑函数，如果应力是C1连续，位移和速度就不得不是C2函数；在高于一维的情况下C2函数是不容易构造的。而且，不得不随之构造C2试函数以便于满足面力边界条件，这也是困难的。通过分部积分消去应力的导数，结果是对应力函数降低了平滑性，可以是C0函数。在线性化方程中也导致了某些对称性，这将在第6章中见到。另外，消除了关于应力P的导数，补充了力边界条件项(第4项)，否则，力边界条件不得不施加在试位移函数P上(第1项)。建立弱形式中的关键步骤是分部积分：d(uv)/dx=udv/dx+vdu/dx当前第22页\共有68页\编于星期六\16点弱形式到强形式2完全的Lagrangian格式弱形式给出

由虚位移的任意性，试证明得到强形式

(参考节)：动量方程力边界条件内部连续条件

如果允许较低平滑的变分项和试函数，在强形式中将附加一个方程－内部连续条件。如果选取的变分项和试函数满足较高的平滑条件，在强形式中则没有内部连续条件。当前第23页\共有68页\编于星期六\16点

对于平滑的变分项和试函数，弱形式仅采用动量方程和力边界条件。

较低平滑性要求的变分项和试函数仅是连续，需要处理在横截面上和材料参数中的不连续点。

在材料界面处，应变，即位移场的导数是不连续的，因此，经典强形式是不适用的，因为它假设任何点的二阶导数是唯一定义的。采用粗糙的变分项和试函数，在这些界面上自然出现附加条件－内部连续条件。在TL弱形式中，所有的积分都是在材料域上进行，比如参考构形。由于求导是对材料坐标X进行，所以在材料域上应用分部积分是最方便的。2完全的Lagrangian格式当前第24页\共有68页\编于星期六\16点2完全的Lagrangian格式虚功项的物理名称

外力虚功内力虚功惯性虚功虚功原理

虚功方程是动量方程、力边界条件和应力跳跃条件的弱形式。当前第25页\共有68页\编于星期六\16点

弱形式中包含强形式，强形式中包含弱形式，所以弱和强形式是等价的。对于动量方程，强和弱形式的这种等价称为虚功原理。2完全的Lagrangian格式当前第26页\共有68页\编于星期六\16点

以弱形式作为虚功表达式的观点提供了统一性，对于在不同坐标系上和不同类型问题中建立弱形式是很有用途的：为了获得弱形式，只需要写出虚能量方程。因此，可以避免前面所做的由变分项与方程相乘并进行各种处理的过程。

从数学观点来看，没有必要考虑变分函数作为虚位移：它们是简单的变分函数，满足连续条件和在位移边界上为零。对于有限元方程的离散，变分函数与方程的乘积没有物理意义。但是，对于物理问题满足自然变分原理还是有其科学意义的。

建立弱形式中的关键步骤是分部积分，从而消除了关于应力P的导数。如果没有这一步，力边界条件就不得不强加在试函数上，或构造C1的应力场函数；作为弱形式，由分部积分和降低对应力和试位移函数平滑性的要求（如C0）是更方便的。2完全的Lagrangian格式当前第27页\共有68页\编于星期六\16点3有限元离散，单元和总体矩阵当前第28页\共有68页\编于星期六\16点3.1TL的有限元离散有限元近似

通过对变分项和试函数应用有限元插值，由虚功原理得到有限元模型的离散方程。

有限元试函数

是连续插值函数，称为形函数。形函数满足条件：是Kroneckerdelta或单位矩阵：当I=J时；当IJ时；运动学条件，试函数u要满足连续性和位移边界条件。方程表示变量分离：解的空间相关性由形函数表示，而时间相关性归属于节点变量。3有限元离散，单元和总体矩阵当前第29页\共有68页\编于星期六\16点节点力

3有限元离散，单元和总体矩阵为了建立有限元方程，要为每个虚功项定义节点力

外力虚功内力虚功惯性虚功这些名称给节点力赋予了物理意义：内部节点力对应于在材料内部的应力，外部节点力对应于外部施加的荷载，动态或惯性节点力对应于惯性力。节点力与节点位移是功共轭的，一个节点位移增量与节点力的标量积给出功的增量，一旦违背，质量和刚度矩阵的对称性将被破坏。当前第30页\共有68页\编于星期六\16点节点力

内部节点力是由固体对变形的阻力而引起的节点力；

外部节点力

惯性节点力

3有限元离散，单元和总体矩阵每一个虚功项的节点力表达式代入虚功原理给出由于的任意性，在所有节点除了位移边界外，即节点1，它服从

当前第31页\共有68页\编于星期六\16点运动方程－半离散方程

3有限元离散，单元和总体矩阵

在模型中，节点1的加速度并不是未知的，它是一个给定位移的节点。可以通过给定节点位移对时间求二次导数，得到给定位移节点的加速度。这个给定的位移必须足够光滑，可求导二次；这要求它是时间的C1函数。MIJ－J处位移对I处惯性力贡献的质量。

在矩阵形式中，不能简单地表示给定位移边界条件，所以必须考虑指标形式(上式)以补充。

如果质量矩阵不是对角阵时，给定边界位移对不在边界的节点也作贡献；对于对角质量阵，不出现下式右端项。

当前第32页\共有68页\编于星期六\16点运动方程－半离散方程

3有限元离散，单元和总体矩阵

从弱形式的一致性推导出的质量矩阵称为一致质量矩阵。在许多应用中，采用对角质量矩阵(集中质量矩阵)更有优势。质量矩阵对角化的过程是相当特殊的，这些过程没有理论。最常用的一种过程是对矩阵的行求和式中用到了这样的事实，形函数对行求和必须等于1。这种对角化的过程使物体的总动量守恒，例如对于任意的节点速度，对角质量的系统动量应该等价于一致质量的系统动量：对角质量矩阵也可以由下式赋值当前第33页\共有68页\编于星期六\16点

运动方程在空间是离散的，在时间上是连续的，有时简称离散方程。在有限元离散中，质量矩阵常常为非对角阵(一致质量矩阵)，此时运动方程区别于牛顿第二定律，当MIJ≠0时，节点J处的力可以在节点I处产生加速度。而集中质量矩阵的运动方程等价于牛顿第二定律。3有限元离散，单元和总体矩阵上式为在质点I上的静力。由牛顿第三定律，作用在节点上的力大小相等，而方向相反，因此内部节点力需要一个负号。

运动方程－半离散方程(矩阵形式)

当前第34页\共有68页\编于星期六\16点单元和总体矩阵

在有限元程序中，通常以一个单元水平计算节点力和质量矩阵，将单元节点力结合入总体矩阵，称为集合或矢量组合。

组合单元的质量矩阵和其它方阵到总体矩阵，称为矩阵装配。

通过计算可以从总体矩阵中提取单元节点位移，称为离散。3有限元离散，单元和总体矩阵当前第35页\共有68页\编于星期六\16点2节点单元一维网格的集合和离散运算的描述，两组单元，计算节点位移：位移根据单元节点编号离散；计算节点力：节点力根据节点编号返回总体力矩阵。3有限元离散，单元和总体矩阵当前第36页\共有68页\编于星期六\16点2节点单元一维网格的单元形函数Ne(X)和总体形函数N(X)3有限元离散，单元和总体矩阵单元节点位移与总体节点位移的关系为

Le为连接矩阵。类似的获得单元节点力。

应用连接矩阵还可以建立单元形函数和总体形函数之间的关系，总体位移场可以由所有单元的位移求和得到：对单元形函数求和得到总体形函数

当前第37页\共有68页\编于星期六\16点3有限元离散，单元和总体矩阵当前第38页\共有68页\编于星期六\16点例题2.12节点线性位移单元

3有限元离散，单元和总体矩阵单元的初始长度为l0，横截面面积为常数A0，随时间变形后，长度为

l(t)，面积为A(t)，位移场，应变和B0矩阵：由线性插值和材料坐标式中，当前第39页\共有68页\编于星期六\16点3有限元离散，单元和总体矩阵以节点位移的形式为应变度量赋值由上式定义B0矩阵节点内力：假定横截面面积和名义应力P

为常数，被积函数是常数，则有单元节点力大小相等、方向相反，满足单元平衡。单元节点力的平衡性质应用于所有发生移动但没有变形的单元，但不能应用于轴对称单元。因为，节点力等于单元承担的载荷T。当前第40页\共有68页\编于星期六\16点3有限元离散，单元和总体矩阵节点外力：如果用线性Lagrange插值近似体积力，，则有代入上面公式，取A0为常数，得到节点外力值当前第41页\共有68页\编于星期六\16点3有限元离散，单元和总体矩阵单元质量矩阵：上式表明质量矩阵与时间无关，它仅取决于初始密度，初始横截面面积和初始长度。由对行求和技术的公式得到对角质量矩阵为上式表明在2节点杆单元的对角质量矩阵中，每个节点分配单元的一半质量。称其为集中质量矩阵，每个节点集中一半的质量。当前第42页\共有68页\编于星期六\16点4更新的Lagrangian格式当前第43页\共有68页\编于星期六\16点4.1UL的控制方程初始构形参考构形当前构形变形构形4更新的Lagrangian格式的控制方程,弱形式，单元方程当前第44页\共有68页\编于星期六\16点4更新的Lagrangian格式的控制方程,弱形式，单元方程UL格式是TL格式的一个简单转换。在数值上，离散方程是相同的，而实际在同一程序中，对某些节点力我们可以应用TL格式，而对其它的节点力应用UL格式。为什么采用两种方法，而它们基本上是一致的。4.1UL的控制方程

主要原因是它们都在被广泛地应用，因此，为了理解程序和文献，有必要熟悉两种格式。

当前第45页\共有68页\编于星期六\16点应变的度量由变形率给出

应力的度量

Cauchy应力

4UL格式的控制方程，弱形式，单元方程4.1UL的控制方程以Eulerian坐标表述相关变量，空间坐标速度应变

UL格式的两个相关变量－速度和Cauchy应力

当前第46页\共有68页\编于星期六\16点质量守恒

对于杆

动量守恒

4UL格式的控制方程，弱形式，单元方程能量守恒

－热流量－热源本构方程

变形度量

4.1UL的控制方程当前第47页\共有68页\编于星期六\16点边界条件

速度边界等价位移边界

力边界n

单位法线（＋,－）

一端固定一端自由杆

边界条件满足

初始条件

4UL格式的控制方程，弱形式，单元方程4.1UL的控制方程当前第48页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程4.2UL的弱形式由动量方程乘以变分函数

弱形式－虚功率原理

强形式－虚功率原理的逆过程：动量方程， 力边界条件 内部连续条件

积分在当前域上完成

当前第49页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程4.2UL的弱形式内部虚功率

外力虚功率

惯性力虚功率

弱形式

当前第50页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程4.3UL的单元方程

在一个单元的水平上建立方程，通过装配获得总体方程。相关变量为速度和应力。建立本构方程、质量守恒方程，动量方程。由于质量守恒是一个代数方程，可以容易地计算任意一点的密度。建立半离散方程。单元的速度场为单元的加速度场为

将形函数表示成为材料坐标的函数是非常关键的，它与时间无关，可分离变量求解。如果将形函数由Eulerian坐标表示为形函数的材料时间导数不为零（注意与TL区别），并且不能将加速度表示为同样形函数与节点加速度乘积的形式。

当前第51页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程4.3UL的单元方程Eulerian坐标与单元坐标之间的映射为

当前第52页\共有68页\编于星期六\16点位移可以由相同的形函数进行插值

4UL格式的控制方程，弱形式，单元方程4.3UL的单元方程

形函数与时间无关，通过位移的导数得到速度和加速度，它们与变分函数都可以由同一形函数给出。变形率可以表示为形函数的形式为当前第53页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程4.3UL的单元方程通过B矩阵，将变形率表示为节点速度的形式变形率可以表示为形函数的形式为形函数的空间导数由链规则得到

当前第54页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程4.3UL的单元方程

与TL格式相同，在UL格式中，质量矩阵不随时间变化，在程序中仅计算一次即可。当前第55页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程例2.53节点二次位移单元以单元坐标的形式写出位移和速度场

单元坐标当前第56页\共有68页\编于星期六\16点以单元坐标的形式写出位移和速度场

4UL格式的控制方程，弱形式，单元方程例2.53节点二次位移单元其中:

B矩阵为

变形率给出

当前第57页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程例2.53节点二次位移单元变形率

如果是常数，单元中的变形率是线性地变化，这是节点2位于其它两节点中间时的一种情况(

)。然而，当由于单元的畸变，节点2偏离中间位置时，成为的线性函数，而变形率成为一个有理函数。而当节点2从中间移开时，有可能成为负数，或为零，在这种情况下，空间坐标和单元坐标的映射将不再一一对应。当前第58页\共有68页\编于星期六\16点4UL格式的控制方程，弱形式，单元方程例2.53节点二次位移单元内部节点力

其中

这个表达式与TL格式的内力表达式是相同的。当前第59页\共有68页\编于星期六\16点例2.53节点二次位移单元-检查网格畸变当单元的节点2位于离节点1的1/4单元长度时

在

有

Jacobian

在该点处的当前密度为无穷大。若节点2移动并接近节点1，在部分单元上J<0，这意味着是负密度值,违背了质量守恒。这些情况经常隐藏在数值积分中，因为在高