《清华大学计算固体力学第二次课件一维Lagrangian和Eulerian有限元.ppt》由会员分享，可在线阅读，更多相关《清华大学计算固体力学第二次课件一维Lagrangian和Eulerian有限元.ppt（70页珍藏版）》请在三一办公上搜索。
1、计算固体力学第2章 Lagrangian和Eulerian有限元,第2章 一维Lagrangian和Eulerian有限元,引言完全的Lagrangian格式的控制方程,弱形式有限元离散，单元和总体矩阵更新的Lagrangian格式的控制方程,弱形式，单元方程求解方法Eulerian格式的控制方程，弱形式,有限元方程,1 引言,1 引言,非线性连续体一维模型（杆）的有限元方程,在固体力学中，Lagrangian网格是最普遍应用的，其吸引力在于它们能够很容易地处理复杂的边界条件，并且能够跟踪材料点，因此能够精确地描述依赖于历史的材料。,在Lagrangian有限元的发展中，一般采用两种方法：以L
2、agrangian度量的形式表述应力和应变的公式，导数和积分运算采用相应的Lagrangian（材料）坐标X，称为完全的Lagrangian格式（TL）。2.以Eulerian度量的形式表述应力和应变的公式，导数和积分运算采用相应的Eulerian（空间）坐标x，称为更新的Lagrangian格式（UL）。,非线性与线性公式的主要区别是前者需要定义积分赋值的坐标系和确定选择应力和应变的度量。,两种格式的主要区别在于：TL，在初始构形上描述变量，UL，在当前构形上描述变量。不同的应力和变形度量分别应用在这两种格式中。TL，习惯于采用一个应变的完全度量，UL，常常采用应变的率度量。这些并不是格式的
3、固有特点，在UL中采用应变的完全度量是可能的，并且在TL中可以采用应变的率度量。,尽管TL和UL表面看来有很大区别，两种格式的力学本质是相同的；因此，TL可以转换为UL，反之亦然。,1 引言,对于每一种公式，将建立动量方程的弱形式，已知为虚功原理（或虚功）。这种弱形式是通过对变分项与动量方程的乘积进行积分来建立。在TL格式中，积分在所有材料坐标上进行；在Eulerian和UL格式中，积分在空间坐标上进行。也将说明如何处理力边界条件，因此近似（试）解不需要满足力边界条件。这个过程与在线性有限元分析中的过程是一致的，在非线性公式中的主要区别是需要定义积分赋值的坐标系和确定选择应力和应变的度量。,1
4、 引言,推导有限元近似计算的离散方程。对于考虑加速度（动力学）或那些包含率相关材料的问题，推导离散有限元方程为普通微分方程（ODEs）。这个空间的离散过程称为半离散化，因为有限元仅将空间微分运算转化为离散形式，而没有对时间导数进行离散。对于静力学与率无关材料问题，离散方程独立于时间，有限元离散将导致一组非线性代数方程。,2 完全的Lagrangian格式,2.2 TL的控制方程,初始构形参考构形,当前构形变形构形,2 完全的Lagrangian格式,物体的运动由Lagrangian坐标和时间的函数描述,是在初始域与当前域之间的映射,当材料坐标在初始位置,2 完全的Lagrangian格式,位移
5、差,或者,变形梯度,偏微分的意义？,2 完全的Lagrangian格式,定义Jacobian：作为变形物体的无限小体积,相对于变形前物体微段体积,的比值,应变的度量,在变形前构形中上式为零，它等效于工程应变,应力的度量,Cauchy 应力,名义应力,在多维上没有工程应力的定义。,工程应力,物理应力,初始值，J01,2 完全的Lagrangian格式,推导方程 应用下面方程推导非线性杆：1.质量守恒 2.动量守恒 3.能量守恒 4.变形度量，也常称为应变-位移方程 5.本构方程，描述材料应力与变形度量的关系 另外，要求变形保持连续性，称为协调性要求。,质量守恒,2 完全的Lagrangian格式
6、,对于Lagrangian格式，质量守恒方程为,对于一维杆,动量守恒,由名义应力P 和Lagrangian坐标给出(单位长度的力),如果初始横截面面积在空间保持常数，则动量方程成为,应力在坐标方向的分量,b单位质量的力体力,平衡方程,2 完全的Lagrangian格式,平衡意味着物体处于静止或者以匀速运动,能量守恒,内部功率由变形率的梯度,和名义应力P 的乘积给出,本构方程,不计惯性力，则动量方程成为平衡方程,etc.表示影响应力的其他变量，如温度，夹杂等。,是变形历史的函数。,完全形式,率形式,2 完全的Lagrangian格式,本构方程的例子,1)线弹性材料,完全形式,率形式,2)线性粘弹
7、性材料,etc.表示影响应力的其他变量，如温度，夹杂等。,是变形历史的函数。,完全形式,率形式,2 完全的Lagrangian格式,边界条件,位移边界,力边界 n0 单位法线（,）,一端固定一端自由杆,边界条件满足,初始条件,动量方程是关于X 二阶的（偏微分方程）。因此在每一端，必须描述u 或者,作为边界条件。,2 完全的Lagrangian格式,内部连续条件,跳跃条件,函数的连续性,如果函数的第n 阶导数是连续函数，该函数为,连续,函数是连续可导的（它的一阶导数存在并且处处连续）,在函数中，导数只是分段可导，一维函数不连续发生在点上，二维函数不连续发生在线段上，三维函数不连续发生在面上。,函
8、数本身不连续，xi是不连续点。,动量平衡要求,关于泛函和变分的概念,变 量,函 数,函 数,泛 函,泛 函函数的函数（functional,function of function),当虚位移是真实位移的增量时，虚位移原理We V 中的 V 就是泛函V 的变分。微分是函数的增量，变分是泛函的增量。,w(x)是 x 函数,V(w(x)是 w(x)的泛函,自然变分原理是对物理问题的微分方程和边界条件建立对应的泛函，使泛函取驻值得到问题的解答，但是其未知场函数需要满足一定的附加条件。广义变分原理（或称约束变分方程）不需要事先满足附加条件，采用Lagrange乘子法和罚函数法将附加条件引入泛函，重新构
9、造一个修正泛函，将问题转化为求修正泛函的驻值。称为无附加条件的变分原理。对于罚函数方法，将罚参数取正值，对修正泛函得到的近似解只是近似地满足附加条件，罚参数值越大，附加条件的满足程度就越好。而在实际计算中，罚函数只能取有限值，所以利用罚函数求解只能得到近似解。,2 完全的Lagrangian格式,有限元方法不能直接离散动量方程。为了离散这个方程，需要一种弱形式，称为变分形式，即虚功原理或者虚功率，通过对变分项与动量方程的乘积进行积分来建立的。虚功原理或者弱形式是等价于动量方程和力边界条件的。后者称为经典强形式。,2 完全的Lagrangian格式,2.3 TL的弱形式,强形式到弱形式,弱形式到
10、强形式,对于动量方程和力边界条件,现在建立弱形式，要求：,满足所有位移边界条件并足够平滑，因此确切定义了动量方程中的所有导数。,也假设足够光滑，这样确切定义了所有的后续步骤，并在指定的位移边界条件上为零。这是标准和经典的建立弱形式的方法。尽管它所导致的连续性要求比在有限元近似中更加严格，在我们看到以较少的强制连续性要求所得到的结论之前，我们仍继续采用这种方法。,2 完全的Lagrangian格式,试函数,变分项,强形式到弱形式,取动量方程与变分项的乘积并在全域内积分得到弱形式，给出,2 完全的Lagrangian格式,强形式到弱形式,名义应力P 是一个试位移函数。展开第一项乘积的导数，整理得到
11、,分布积分,在指定位移边界处变分项,消失，第二行服从边界互补条件和力边界条件。,给出完全的Lagrangian格式的动量方程和力边界条件的弱形式,弱形式到强形式,2 完全的Lagrangian格式,弱形式给出,由虚位移的任意性，试证明得到强形式(参考4.3.2节)：,动量方程,力边界条件,内部连续条件,可以看出，如果允许较低平滑的变分项和试函数，在强形式中将附加一个方程内部连续条件。如果选取的变分项和试函数满足经典的平滑条件，在强形式中则没有内部连续条件。对于平滑的变分项和试函数，弱形式仅采用动量方程和力边界条件。较低平滑性要求的变分项和试函数仅是 连续，需要处,理在横截面上和材料参数中的不连
12、续点。在材料界面，经典强形式是不适用的，因为它假设任何点的二阶导数是唯一定义的。然而，在材料界面处，应变，即位移场的导数是不连续的。采用粗糙的变分项和试函数，在这些界面上自然出现附加条件内部连续条件。在TL弱形式中，所有的积分都是在材料域上进行，比如参考构形。由于求导是对材料坐标X进行，所以在材料域上应用分部积分是最方便的。,2 完全的Lagrangian格式,2 完全的Lagrangian格式,虚功项的物理名称,外力虚功,内力虚功,惯性虚功,虚功原理,方程是动量方程、力边界条件和应力跳跃条件的弱形式。,弱形式中包含强形式，并且强形式中包含弱形式，所以弱和强形式是等价的。对于动量方程，强和弱形
13、式的这种等价称为虚功原理。,2 完全的Lagrangian格式,以弱形式作为虚功表达式的观点提供了统一性，对于在不同坐标系上和不同类型问题中建立弱形式是很有用途的：为了获得弱形式，只需要写出虚能量方程。因此，可以避免前面所做的由变分项与方程相乘并进行各种处理的过程。从数学观点来看，没有必要考虑变分函数作为虚位移：它们是简单的变分函数，满足连续条件和在位移边界上为零。对于有限元方程的离散，方程与变分函数的乘积没有物理意义。建立弱形式中的关键步骤是分部积分，从而消除了关于应力P 的导数。如果没有这一步，力边界条件就不得不强加在试函数上。作为弱形式，由分部积分和降低对应力和试位移平滑性的要求是更方便
14、的。,2 完全的Lagrangian格式,3 有限元离散，单元和总体矩阵,3.1 TL的有限元离散,有限元近似,通过对变分项和试函数应用有限元插值，由虚功原理得到有限元模型的离散方程。,有限元试函数,是,连续插值函数，称为形函数。,形函数满足条件：,是Kronecker delta或单位矩阵：,当I=J时；,当I J时；,运动学条件，试函数u要满足连续性和基本边界条件。方程表示变量分离：解的空间相关性由形函数表示，而时间相关性归属于节点变量。,3 有限元离散，单元和总体矩阵,节点力,3 有限元离散，单元和总体矩阵,为了建立有限元方程，要为每一个虚功项定义节点力,外力虚功,内力虚功,惯性虚功,这
15、些名称给节点力赋予了物理意义：内部节点力对应于在材料内部的应力，外部节点力对应于外部施加的荷载，而动态或惯性节点力对应于惯性。节点力与节点位移是功共轭的，一个节点位移的增量与节点力的标量积给出功的增量，一旦违背，质量和刚度矩阵的对称性将被破坏。,节点力,内部节点力是由固体对变形的阻力而引起的节点力；,外部节点力,惯性节点力,3 有限元离散，单元和总体矩阵,每一个虚功项节点力表达式,代入虚功原理,给出,由于 的任意性，在所有节点除了位移边界外，即节点1，它服从,运动方程半离散方程,3 有限元离散，单元和总体矩阵,在模型中，节点1的加速度并不是未知的，它是一个给定位移的节点。可以通过给定节点位移对
16、时间求二次导数，得到给定位移节点的加速度。这个给定的位移必须足够光滑，可求导二次；这要求它是时间的C1函数(细长梁模型）。,当质量矩阵不是对角阵时，给定位移对没有在边界的节点也作出贡献。对于对角质量阵的情况，不出现下式右端项。,MIJJ处位移对I处惯性力贡献的质量。,在矩阵形式中，不能简单地表示给定位移边界条件，所以必须考虑指标形式(上式)以补充。,运动方程半离散方程(矩阵形式),运动方程在空间是离散的，在时间上是连续的，有时简称离散方程。在有限元离散中，质量矩阵常常为非对角阵(一致质量矩阵)，运动方程区别于牛顿第二定律，当MIJ0时，节点I处的力可以在节点J处产生加速度。而集中质量矩阵的运动
17、方程等价于牛顿第二定律。,3 有限元离散，单元和总体矩阵,为在质点I上的净力。由牛顿第三定律，作用在节点上的力大小相等，而方向相反，因此内部节点力需要一个负号。,单元和总体矩阵,在有限元程序中，通常以一个单元水平计算节点力和质量矩阵，将单元节点力结合入总体矩阵，称为离散或矢量组合。组合单元的质量矩阵和其它方阵到总体矩阵，称为矩阵装配。通过计算可以从总体矩阵中提取单元节点位移，称为集合。,3 有限元离散，单元和总体矩阵,2节点单元一维网格的集合和离散运算的描述，两组单元节点位移的集合：位移根据单元节点编号集合；计算节点力的离散：节点力根据节点编号返回总体力矩阵。,3 有限元离散，单元和总体矩阵,
18、2节点单元一维网格的单元形函数Ne(X)和总体形函数N(X),3 有限元离散，单元和总体矩阵,单元节点位移与总体节点位移的关系为,Le为连接矩阵。类似的获得单元节点力。,应用连接矩阵还可以建立单元形函数和总体形函数之间的关系，总体位移场可以由所有单元的位移求和得到：,对单元形函数求和得到总体形函数,3 有限元离散，单元和总体矩阵,例题,3 有限元离散，单元和总体矩阵,3 有限元离散，单元和总体矩阵,3 有限元离散，单元和总体矩阵,3 有限元离散，单元和总体矩阵,4 更新的Lagrangian格式,4.1 UL的控制方程,初始构形参考构形,当前构形变形构形,4 更新的Lagrangian格式的控
19、制方程,弱形式，单元方程,4 更新的Lagrangian格式的控制方程,弱形式，单元方程,UL格式是TL格式的一个简单转换。在数值上，离散方程是相同的，而实际在同一程序中，对某些节点力我们可以应用TL格式，而对其它的节点力应用UL格式。为什么采用两种方法，而它们基本上是一致的。,4.1 UL的控制方程,主要原因是它们都在被广泛地应用，因此，为了理解程序和文献，有必要熟悉两种格式。,应变的度量由变形率给出,应力的度量,Cauchy 应力,4 UL格式的控制方程，弱形式，单元方程,4.1 UL的控制方程,以Eulerian坐标表述相关变量，空间坐标,速度应变,UL格式的两个相关变量速度和Cauch
20、y应力,质量守恒,对于杆,动量守恒,4 UL格式的控制方程，弱形式，单元方程,能量守恒,热流量,热源,本构方程,变形度量,4.1 UL的控制方程,边界条件,速度边界等价位移边界,力边界 n 单位法线（,）,一端固定一端自由杆,边界条件满足,初始条件,4 UL格式的控制方程，弱形式，单元方程,4.1 UL的控制方程,4 UL格式的控制方程，弱形式，单元方程,4.2 UL的弱形式,由动量方程乘以变分函数,弱形式虚功率原理,强形式虚功率原理的逆过程：动量方程，力边界条件 内部连续条件,积分在当前域上完成,4 UL格式的控制方程，弱形式，单元方程,4.2 UL的弱形式,内部虚功率,外力虚功率,惯性力虚
21、功率,弱形式,4 UL格式的控制方程，弱形式，单元方程,4.3 UL的单元方程,在一个单元的水平上建立方程，通过装配获得总体方程。相关变量为速度和应力。建立本构方程、质量守恒方程，动量方程。由于质量守恒是一个代数方程，可以容易地计算任意一点的密度。建立半离散方程。,单元的速度场为,单元的加速度场为,将形函数表示成为材料坐标的函数是非常关键的，它与时间无关。如果将形函数由Eulerian坐标表示为,形函数的材料时间导数不为零（注意与TL区别），并且不能将加速度表示为同样形函数与节点加速度乘积的形式。,4 UL格式的控制方程，弱形式，单元方程,4.3 UL的单元方程,Eulerian坐标与单元坐标
22、 之间的映射为,位移可以由相同的形函数进行插值,4 UL格式的控制方程，弱形式，单元方程,4.3 UL的单元方程,形函数与时间无关，通过位移的导数得到速度和加速度，变分函数由同一形函数给出,变形率可以表示为形函数的形式为,4 UL格式的控制方程，弱形式，单元方程,4.3 UL的单元方程,通过一个B矩阵，将变形率表示为节点速度的形式,变形率可以表示为形函数的形式为,形函数的空间导数由链规则得到,4 UL格式的控制方程，弱形式，单元方程,4.3 UL的单元方程,与TL格式相同，在UL格式中，质量矩阵不随时间变化，在程序中仅计算一次即可。,4 UL格式的控制方程，弱形式，单元方程,例2.5 3节点二
23、次位移单元,以单元坐标的形式写出位移和速度场,以单元坐标的形式写出位移和速度场,4 UL格式的控制方程，弱形式，单元方程,例2.5 3节点二次位移单元,其中:,B矩阵给出为,变形率给出为,4 UL格式的控制方程，弱形式，单元方程,例2.5 3节点二次位移单元,如果 是常数，单元中的变形率是线性地变化，这是节点2位于其它两节点中间时的一种情况。然而，当由于单元的畸变，节点2偏离中间位置时，变成为 的线性函数，而变形率变成为一个有理函数。而当节点2从中间移开时，有可能成为负数，或为零，在这种情况下，当前空间坐标和单元坐标的映射将不再一一对应。,变形率,4 UL格式的控制方程，弱形式，单元方程,例2
24、.5 3节点二次位移单元,内部节点力,其中,这个表达式与TL格式的内力表达式是相同的。,4 UL格式的控制方程，弱形式，单元方程,例2.5 3节点二次位移单元-检查网格畸变,当单元的节点2是位于离节点1的1/4单元长度时,在,有,Jacobian,在该点处的当前密度为无穷大。若节点2移动并接近节点1，在部分单元上Jacobian成为负数，这意味着是负密度值,违背了质量守恒。这些情况经常隐藏在数值积分中，因为在高斯积分点，当Jacobian成为负数时,畸变是非常严重的。,由,不能满足一一对应条件也可能导致变形率 出现奇异。当分母 为零或成为负数，我们难以得到势能。,例2.5 3节点二次位移单元-
25、检查网格畸变,4 UL格式的控制方程，弱形式，单元方程,在节点1处变形率为无穷大。在断裂力学中，利用这种二次位移单元的性质建立包含裂纹尖端奇异应力的单元，称为四分之一点单元。但是在大位移分析中，这种行为会出现问题。在一维单元中，网格畸变的影响不像在多维问题中那么严重。事实上，应用变形梯度F 作为这种单元的变形度量多少可以减轻网格畸变的影响。在3节点单元中，如果X2 的初始位置位于中点，那么变形梯度F 绝不会成为奇异。,例2.5 3节点二次位移单元-检查网格畸变,4 UL格式的控制方程，弱形式，单元方程,例2.5 3节点二次位移单元-检查网格畸变,4 UL格式的控制方程，弱形式，单元方程,3节点
26、1/4点二次位移单元,例2.5 3节点二次位移单元-检查网格畸变,5 求解方法,5 求解方法,为了求解非线性问题，最简单的方法，即时间显式积分。最广泛应用的显式方法是中心差分方法，采用对角或集中质量矩阵。,速度,加速度,在时间间隔中点的导数值由在间隔端点处函数值的差得到，顾名思义为中心差分公式。,从t0出发，取时间步长t,5 求解方法,5 求解方法,对于位移的更新不需要代数方程的任何解答，因此，在某种意义上，显式积分比静态线性应力分析更加简单,不需要矩阵求逆解刚度方程。如在流程图中看到，对于控制方程和时间积分公式，大多数的显式程序是直接向前赋值。程序从施加初始条件开始,第一个时间步多少与其它时
27、间步的不同在于它仅取半步，这使程序能正确地解释关于应力和速度的初始条件。大部分程序运算时间是在计算单元节点力，尤其是内部节点力。节点力是逐个单元进行计算的。在开始计算前，从总体的列矩阵中集合出单元节点速度和位移。如流程图所示，内部节点力的计算包括应变方程和本构方程的应用。通过应力为内部节点力赋值。当完成了单元节点力的计算，根据它们的节点编号将其离散到总体列矩阵。,5 求解方法,稳定性准则 显式积分的缺陷在于时间步长必须低于一个临界值，否则由于数值不稳定将使解答隆起。对于采用对角质量的2节点单元的临界时间步长为,是单元的初始长度,是波速,稳定性准则 能量守恒(断裂力学中为能量平衡),增加时间步长
28、：放大质量，调整单元尺寸。,6 Eulerian格式的控制方程,6 Eulerian格式的控制方程，弱形式，有限元方程,在Eulerian格式中，节点在空间固定，相关变量为Eulerian空间坐标 x 和时间 t 的函数，应力度量为Cauchy（物理的）应力，变形度量为变形率，运动由速度描述。在Eulerian格式中，因为不能建立未变形、初始的构形，所以不能将运动表示为参考坐标的函数。,TL格式比UL格式需要更多的存储空间，以存储形函数及其导数值；而UL则需要在每一个时间步重复搜索和计算形函数，也会影响计算效率。因此，在实际问题中应有所选择，例如对于大变形的瞬态问题，或者路径无关材料，可采用TL，而与变形历史有关的路径相关材料，如弹塑性和粘弹塑性材料，则可采用UL。,
