任意拉格朗日—欧拉法.pdf

任意拉格朗日2欧拉描述法研究进展Ξ张 雄 陆明万 王建军(清华大学工程力学系,北京, 100084)摘 要 任意拉格朗日2欧拉(AL E)描述综合了纯拉格朗日和纯欧拉描述的优点,克服了各自的缺点,成为非线性连续介质力学中大变形分析的非常有效的方法本文论述了AL E法的研究进展及其在流体动力学、 流体2结构相互作用、 加工成型、 碰撞、 接触等大变形问题中的应用关键词 非线性连续介质力学;任意拉格朗日2欧拉描述;非线性有限元分类号 O33;O35;O 2421211 引 言非线性连续介质力学的有限元分析方法已取得了很多令人瞩目的进展,其中在大变形问题中绝大多数研究工作都采用拉格朗日(主要用在固体力学中)或欧拉(主要用在流体力学中)描述方法在拉格朗日描述(即物质描述或L描述)中,计算网格固定在物体上随物体一起运动,即网格点与物质点在物体的变形过程中始终保持重合,因此物质点与网格点之间不存在相对运动(即迁移运动,也称对流运动) 这大大地简化了控制方程的求解过程,而且能准确描述物体的移动界面,并可跟踪质点的运动轨迹 但在涉及到特大变形的问题中,物质的扭曲将导致计算网格的畸形而使得计算失败。

在欧拉描述(即空间描述或E描述)中,网格固定在空间中,即计算网格在物体的变形过程中保持不变,因此可很容易处理物质的扭曲 但对运动界面需要引入非常复杂的数学映射,将可能导致较大的误差 另外当使用普通的伽辽金离散时,由于迁移项的影响,有限元方程中的系数矩阵是非对称的,而且还可能得到振荡解纯拉格朗日和纯欧拉描述都存在严重的缺陷,但也具有各自的优势 如果能将二者有机地结合,充分吸收各自的优势,克服各自的缺点,则可解决一大批只用纯拉格朗日和纯欧拉描述所解决不了的问题 任意拉格朗日2欧拉(AL E)方法就是基于此目的最早由Noh(1964)[1]以耦合欧拉2拉格朗日的术语提出的,并用有限差分法求解带有移动边界的二维流体动力学问题在Noh的研究工作中,网格点可以随物质点一起运动,但也可以在空间中固定不动,甚至网格 点可以在一个方向上固定,而在另一个方向上随物体一起运动,因此AL E描述也被称为耦合欧拉2拉格朗日描述 例如在液体表面波的传播问题中,网格点在垂向随物质点一起运动,而在第14卷第1期计 算 力 学 学 报Vol . 14 No. 1 1997年2月CH I N ESE JOURNAL OF COM PU TA TLONAL M ECHAN ICSFebruary 1997Ξ核工业科学基金资助。

本文于1996年5月27日收到 © 1995-2004 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved.水平方向上固定不动这样可很容易描述液体表面的运动,而且网格不会发生扭曲 在AL E描述中计算网格可以在空间中以任意的形式运动,即可以独立于物质坐标系和 空间坐标系运动 这样通过规定合适的网格运动形式可以准确地描述物体的移动界面,并维持 单元的合理形状类似于E描述,在AL E描述下的控制方程中也将出现对流项,因此也可能 得到振荡解,需要进行相应的数值处理纯拉格朗日和纯欧拉描述实际上是AL E描述的两个 特例,即当网格的运动速度等于物体的运动速度时就退化为拉格朗日描述,而当网格固定于空 间不动时就退化为欧拉描述AL E最早是为了解决流体动力学问题而引入的,并且使用有限差分法[128]由于核反应堆结构安全分析的需要,Donea[12]、Belytschko[13, 14]等人分别将AL E法引入有限元法中,用以求 解流体与结构相互作用问题[12230]Hughes等人[9]建立了AL E描述的运动学理论,并使用有限 元法解决了粘性不可压缩流体流动和自由表面流动问题。

AL E描述中,参考构形是已知的,而初始构形和现时构形都是待求解的因此AL E法尤 其适合于在初始构形和现时构形都未知的问题中使用,如接触问题等目前,AL E法已被应用于固体力学领域中求解大变形问题,如碰撞[36][55]、 接触和弹性断 裂力学[31～34]、 路径相关材料(如弹塑性材料等)[35～39]、 加工成型[40～47]等另外对网格的运动以及控制方程的求解策略等问题也进行了多方面的研究[53～58]本文论述了AL E法的基本理论和研究进展及其在流体动力学、 流体2结构相互作用、 加工 成型、 碰撞、 接触等大变形问题中的应用2 AL E运动学描述将连续体在初始时刻t0的构形(称初始构形)记为8X,将t时刻的构形(称现时构形)记为8x 为了确定各物质点的位置,引入L agrange或物质坐标系OX1X2X3(简记为L坐标系),它是和物体固结在一起并随物体一起运动的 质点在初始时刻to的位置由其在L坐标系中的位置矢量X完全确定,因此X也可用来识别物体中的不同质点 为了描述现时构形,引入空间坐 标系或Euler坐标系ox1x2x3(简记为E坐标系),它是和空间固结在一起的 空间中各点的位置由其在E坐标系中的位置矢量x确定,因此可用x表示空间中的几何点。

拉格朗日描述以初始构形为参考构形研究物质点X在空间中的运动规律,即:x=x(X,t)(11a)而欧拉描述则是以现时构形为参考构形来研究空间点x上物质点的运动规律,即研究:X=X(x,t)(11b)式(11a)描述了同一质点X在不同时刻的空间位置,而式(11b)则描述了同一空间点x在 不同时刻被物质点占据的情况 不同于拉格朗日和欧拉描述,ALE描述另外引入了一个可以独立于初始构形和现时构形 运动的参考构形,记为 8Ν 在物体的变形过程中,观察者始终跟随参考构形运动,因而对观察者而言参考构形是固定不动的,而初始构形和现时构形则都相对于参考构形运动 为了确定参 考构形中各参考点的位置,引入参考坐标系oΝ1Ν2Ν3,参考构形中各点的位置由其在参考坐标系中的位置矢量 Ν确定,因此 Ν=Ν(X,t)(11c) 描述了物质点X在参考坐标系中的运动规律,而29计 算 力 学 学 报 14卷© 1995-2004 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved.x=x(Ν,t)(11d) 则描述参考点 Ν在空间中的运动规律。

在物质描述方法中的有限单元剖分是对物体进行剖分的,网格点就是物质点,即网格是随 物体一起运动的 在空间描述方法中的有限单元剖分是对空间进行的,网格点就是空间点,因 此网格是固定在空间中不随物体运动的 而在ALE描述方法中的有限单元剖分是对参考构形 进行的,网格点就是参考点,网格是独立于物体和空间运动的,可以根据需要自由选择 以上映射应该是一对一映射,雅可比行列式J(描述从初始构形到现时构形的映射关系) 和混合雅可比行列式Jδ(描述从参考构形到现时构形的映射关系)都应不等于0,即:J=ß5xi 5Xjß ≠0, Jδ=ß5xi 5Νjß ≠0(2)根据定义t时刻某质点X在空间中的运动速度v等于质点X在空间的位置矢量x(X,t)对 时间的导数,即:v=5x(X,t) 5tßX(3)参考构形中某点Ν在空间中的运动速度(也就是网格点的运动速度)vδ则等于参考点(即网格点)Ν在空间中的位置矢量x(Ν,t)对时间的导数,即vδ=5x(Ν,t) 5tßΝ(4)物质点X在参考坐标系中的位置矢量 Ν(X,t)对时间的导数w=5Ν(X,t) 5tßX(5)则为物质点X在参考坐标系中的速度 在AL E描述中,参考构形(也就是计算网格)的运动规律可以是任意给定的,指定特殊的 网格运动规律可以将AL E描述退化为L描述和E描述。

ı vδ= 0,即计算网格在空间中固定不动,退化为欧拉描述 ı vδ=v,即计算网格随同物体一起运动,退化为拉格朗日描述 ı vδ≠v≠0,即计算网格在空间中独立运动,对应于一般的AL E描述 在ALE描述中参考构形是已知的,各物理量用参考坐标Ν描述比较方便,即F=F(Ν, t) 因此各物理量的物质导数(即某固定质点X的物理量对时间的变化率)应通过该物理量的参 考导数(即某固定参考点 Ν处的物理量对时间的变化率,也称混合导数)来计算 为了表达方 便,本文采用指标求和法则 对x求物质导数,得:vi=vδ i+5xi(Ν,t) 5Νjıwj(6)将式(6)整理得:ci=5xi(Ν,t) 5Νjwj(7)式中:ci=vi-vδ i为物质点相对于网格点的运动速度,即迁移速度(也称对流速度) 从式(7)可 以看出,迁移速度c实质上是将质点在参考坐标系中的运动速度w变换到现时坐标系中而得 到的,也就是说c是在现时坐标系中物质点相对于参考构形的运动速度 式(7)也可化为:ci=5xi(Ν(X,t),t) 5tßX(8)391期 张 雄等:任意拉格朗日2欧拉描述法研究进展 © 1995-2004 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved.利用式(7)可将任意物理量F(Ν,t)的物质导数表示为:5F 5tßX=5F(Ν,t) 5tßΝ+ci5F 5xi(9)式(9)是AL E描述中很重要的关系式,利用它可将控制方程转化到参考坐标系下求解。

ALE描述的重要特征是可以根据需要给定合适的网格运动速度vδ i,以维持计算网格的合理形状并准确描述物体的移动界面 所以网格运动算法在ALE描述中占很重要的地位3 AL E描述的控制方程为了建立在AL E描述下的控制方程,利用非线性连续介质力学的方法以及高斯定理可以推出物理量G(t) =∫8Νg(Ν,t)dVΝ 和 H(t) =∫8xh(x,t)dVx(10)的物质导数为:5G 5tßX=∫ 8Ν(5g 5tßΝ+5wig 5Νi)dVΝ(111a)5H 5tßX=∫ 8x(5h5tßx+5vih 5xi)dVx(111b)311 质量守恒方程(连续性方程) 考察一连续体,用58X、 58x和58Ν分别表示物质域8X、 空间域8x和参考域8Ν的边界,而用Θ0、 Θ和 Θδ分别表示连续体各构形的密度 在不同构形中该连续体的质量M可以写为:M=∫ 8ΝΘδdV Ν=∫8xΘdVx=∫ 8XΘ0dVX(12)式中: Θδ(Ν,t) =JδΘ(x,t),Θ0(X,t) =JΘ(x,t)质量守恒定律表明,质量的整体变化率(即物质导数)应等于零 利用式(12)和式(111a)可以得到在参考坐标系下的质量守恒方程为:5Θδ5tßΝ+5Θδwi 5Νi= 0 在 8Ν域中(131a)类似于修正拉格朗日(UL)方法,有时在空间域8x中求解控制方程是比较方便的。

利用式(12)、 式(111b)和式(9)可以把质量守恒方程写成:5Θ 5tßΝ+ci5Θ 5xi+Θ5vi 5xi= 0 在 8x域中(131b)312 动量守恒方程(平衡方程) 动量守恒定律表明,在t时刻占据参考域 8Ν的物体的总动量的整体变化率等于施加在物体上的外力之和,即:5 5tßX∫ 8ΝΘδv idVΝ=∫58Νtδ idSΝ+∫8ΝΘδfidVΝ(14)式中tδ i是作用在参考域8Ν的边界58Ν的单位表面上的力,fi是作用于物体中单位质量的体力利用散度定理以及式(111a)、 式(131a)可将上式写成:Θδ5vi 5tßΝ+Θδwj5vi 5Νj=5Tδ ji 5Νj+Θδfi 在 8Ν域中(151a)49计 算 力 学 学 报 14卷© 1995-2004 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved.式中Tδ ji是定义在参考构形下的第一类皮奥拉 2克希荷夫应力张量(即拉格朗日应力张量) 它 与真实应力(即柯西应力)张量 Ρij之间的关系为:Tδ ij=Jδ5Νi 5xkΡkj(。