显式有限元和隐式有限元

合集下载

13-2.静力隐式-显式和动力显式有限元列式

13-2.静力隐式-显式和动力显式有限元列式

弹塑性大变形有限元方法胡平§大变形弹塑性本构方程最简单的材料大变形弹塑性本构方程是不考虑加载历史的形变理论的全量本构方程。

这类本构方程基于加载历史是比例加载的基本假定。

对于许多金属冲压成形问题,比例加载的基本假定是近似正确的,因此,基于这类理论的所谓one step inverse algorithm,由于其高效高速的计算效率,近年来被广泛应用于汽车概念设计阶段的冲压工艺性快速校核。

关于全量理论的One-step inverse algorithm的理论知识将另行介绍(参见理论手册)。

然而,更精细和准确的成形问题数值分析应该是能够考虑加载历史的增量本构理论。

为了反映与加载历史的相关性,大变形弹塑性问题需要采用速率型(增量型)的本构方程。

首先,需要在大变形条件下,按照符合客观性的要求,建立增量型本构方程。

由(7.26)式可知,变形率张量[d]是客观张量,但由(7.36)式可知Cauchy应力张量的物质导数[]σ 不是客观张量,所以若用变形率张量和Cauchy应力张量来建立本构方程,则将不满足本构方程的客观性条件(7.42)式。

为此,需要另外定义Cauchy应力张量的导数。

1.Cauchy应力张量的Jaumann导数(率)而Cauchy应力的Jaumann导数定义为01lim tjdt dt(1)推导后得到ijij ik kj kj ki σσσωσω∇=--(2)其矩阵形式为2. 第一Piola 应力(名义应力)张量的本构导数(3)用Cauchy 应力Jaumann 导数表示的第一Piola 应力的本构导数为()ij ij ik kj kj ki ik jk ij kk t t d d l l σσσσσ∇=--++ (4)其矩阵形式为()][t t σ∇= (5)3. 第二Piola 应力张量的本构导数仿照第一Piola 应力本构导数的推导,有()ij ij ij kk kj ik ik jk T t l l l σσσσ=+--(6)把(7.179)式代入(7.182)式,便得用Cauchy 应力Jaumann 导数表示的第二Piola 应力的本构导数为()ij ij ik kj kj ki ij kk T t d d l σσσσ•∇=--+(7)其矩阵形式为(8)4. 大变形弹塑形本构方程三维大变形弹塑形问题属于大变形问题。

显式与隐式算法区别

显式与隐式算法区别

显式(explicit)和隐式(implicit)这两个词在有限元分析中大家可能经常看到,特别是涉及到动力学分析时。

但其实广义的说他们分别对应着两种不同的算法:显式算法(explicit method)和隐式算法(implicit method)。

所以不论在动力学或者静力学中都有涉及到。

显式算法:不直接求解切线刚度,不进行平衡迭代,计算速度快,时间步长只需要足够小,一般不存在收敛问题,需要的内存也小。

隐式算法:每一增量步都需要对静态方程进行平衡迭代,且每次迭代需要求解大量的线性方程组,这一特点使之占用大量的资源。

但该算法增量步可以很大,至少比显式算法大的多,实际计算中会受到迭代次数及非线性程度的影响我们都知道有限元分析FEA在计算微分方程(differential equations)时,由于计算本身的局限,比如计算机储存的位数有限,以及方程本身的复杂性,计算机运用的是数值算法(numerical algorithm)来逼近真实解的。

有限元分析中数值算法的基础是欧拉法(Euler method),欧拉法又分为forward Euler method 和backward Euler method,这两种方法被简称为显式法(explicit method)和隐式法(implicit method)。

中心差分法:(动力学分析)用有限差分代替位移对时间的求导,将运动方程中的速度与加速度用位移的某种组合来标示,这样就将常微分方程组的求解问题转化为代数方程组的求解问题,并假设在每个小的时间间隔内满足运动方程。

首先我们来看看这两种算法的区别。

显式算法(explicit method )(forward Euler method )考虑常微分方程:初始条件:设为每一步的时间步长, 在Tn 时刻,. (n=0,1,2,3...),在T(n+1)时刻有:所以在显式算法中,T(n+1)时刻的值由T(n)时刻决定,也就是说当前时刻的值由上一时刻的值决定。

有限元计算注意事项

有限元计算注意事项

1、计算原理任何有限元模拟的第一步都是利用一个有限单元的集合离散结构的实际几何形状,每个单元(element)代表这个实际结构的一个离散部分.这些单元通过公用的节点(node)来连接.节点和单元的集合成为网格(mesh).在一个特定网格中的单元数目称为网格密度(mesh density).在ansys计算过程中,程序以每个节点的每个自由度建立平衡方程,以节点的位移作为未知量,利用矩阵求解节点的位移.一旦节点位移求出,整个结构的应力和应变都很容易计算出来.这种计算的过程和方法,数学上称之为隐式方法.从上叙述来看,整个计算过程中就是求解n个n元一次方程组(n表示节点数量),当计算模型复杂而且庞大时,隐式求解方法的计算量还是相当大的.与之相对应的,显式求解方法.显式求解方法是通过动态方法从一个增量步前推到下一个增量步得到的.具体显式求解方法和隐式求解方法例子如下:(1)隐式求解(2)显式求解隐式求解中,计算的精度完全控制于计算步数,在一般的计算软件中(flac、abaqus),软件均是利用不平衡力来控制计算步数(当不平衡力<10-5时,停止计算).不平衡力=A+B.A表示施加在节点上的集中力;B表示:在n步数下,根据第n步计算出来的应力,求出节点的内力.Flac软件中 B,以上公式是根据虚功原理推倒而得到.具体推倒过程见《flac 原理》.2、Ansys计算注意事项:计算单位、参数、荷载、标准值、设计值,计算过程中系数的加入. (1)b eam单元对于beam单元.Ansys软件中我们常用的有两种梁单元:beam188和beam4.这两种单元均是三维的梁单元,每个节点都具有6个自由度(ux、uy、uz、mx、my、mz),并且单元坐标系x轴是i点指向j 点.Beam188单元是基于铁木辛哥理论的梁,beam4单元是我们常用的经典结构力学梁.(铁木辛哥理论考虑了梁的剪切变形,而我们常用的经典结构力学梁只考虑了弯矩对结构的变形影响)所以说,beam188可以更精确的计算梁单元,因此我们结构计算中,一般都采用beam188单元.当然还有beam189单元,189单元属于三维二次的梁单元(beam188属于三维一次梁单元),精度比beam188更加高.定义beam188单元,一般采用如下形式:!定义单元/prep7 !进入前处理et,1,beam188 !定义单元188号标号为1!定义材料属性mp,ex,1,2.55e7 !定义弹性模量(kn/m2)mp,nuxy,1,0.167 !定义泊松比mp,dens,1,2.5 !定义密度(KN/N*KG/M3)nummrg,all !合并重合节点numcmp,all !压缩编号!定义梁截面SECTYPE,1,BEAM,RECT,A1,0 ! 1表示梁编号 ; RECT表示是矩形梁(还有其他t型等等,具体见ansys帮助); A1 表示梁的名称 ; 0表示薄壁梁单元网格划分精细程度(0~5). SECDATA,1,3,4,12 !1表示梁b ; 3表示梁h ; 4和12定义对应宽长等分份数.SECOFFSET,CENT !cent质心 ; shrc剪切中心 ; origin原始中心 ; user用户定义;!注意:当梁单元和壳单元一起使用时,可以设置梁单元的偏心,使梁的一面和壳的一面共面.(secoffset,user,offset-y,offset-z),如下图:!划分网格LSEL,S,,,1 !选中编号为1的线.LATT,1,1,1,,,,1 !mp,r,et,,方向点,,SECTYPE截面号.LESIZE,ALL,0.2,,,,,,,1 !0.2是单元大小,1是确认细分规则.LMESH,ALL !用beam单元离散模型,形成网格.!对于划分网格,空间的beam单元,由于需要确定b、h的方向,ansys软件利用方向点来控制b、h的方向.方向点的编号最好定义的很大,如果定义太小,会影响后面的加载.具体方向点如何控制见上面的latt命令和ansys帮助.自己试两下就知道怎么用了. AllsFINISH!加载加约束/SOLUACEL,,,9.8 !重力加速度.注意方向,数值和整体坐标相反,比如重力指向z轴负向,则为正值.SFCUM,ALL,ADD !设置单元荷载是叠加还是替代,只对加在单元和节点上的荷载有效,对于加在面、线上的荷载,都只有替代作用(对同一个面,第二次加的荷载替代第一次加的荷载)!对于beam单元,只能根据sfbeam命令增加均布荷载①等大小的均布荷载.Lsel,s,,,1ESLL,S,1sfbeam,ALL,1,PRES,-161.5 !1表示作用在beam单元的①面上(如下图,③面表示beam单元的轴向,②面表示单元侧面,①面表示beam单元顶面),-161.5表示均布荷载大小,正负号可以控制作用力的方向.②梯形均布荷载Sfbeam命令是对每个单元进行加载.如果一根梁承受10~100的梯形均布荷载,而且这根梁被分成了10个beam单元,这样施加荷载就非常困难.因此我将这种加载过程写成命令流,让软件自动进行加载.命令流如下:LSEL,S,,,1 !选中要加载的那根梁(线)ESLL,S,1 !选中属于这根梁(线)的beam单元*GET,Nelem,ELEM,,COUNT, , , , !获得当前所选单元个数,赋予参数Nelem*GET,Ne,ELEM,,NUM,MIN, , , , !获得当前所选单元最小编号,赋予参数Ne*DO,I,1,Nelem !循环加载,循环次数=单元个数ESEL,S,,,NeNSLE,S,1*GET,Nnode,NODE,,COUNT, , , , !获得当前所选节点个数,赋予参数Nnode*GET,Nn,NODE,,NUM,MIN, , , , !获得当前所选节点最小编号,赋予参数NnNN1X=NX(NN) !将nn节点的x坐标赋予NN1X(NX表示x坐标,NY表示y坐标) NN=NDNEXT(NN) !NN=当前所选节点的下一个编号NN2X=NX(NN) !将nn节点的x坐标赋予NN2Xsfbeam,ALL,1,PRES,-1630.76/3.23*(3.23-NN1X),-1630.76/3.23*(3.23-NN2X)!以上荷载公式应根据实际情况进行调整LSEL,S,,,1ESLL,S,1NE=ELnext(NE) !NE=当前所选单元的下一个编号*ENDDO!对于此命令流,根据不同的实际情况,ABC部分需要修改,其他不需要修改.!后处理 (XY平面) (大拇指指向y,就是my)etable,ImY,smisc,2 !显示弯距etable,JmY,smisc,15PLLS,IMY,JMYETABLE,IFX,SMISC,1 !显示轴力ETABLE,JFX,SMISC,14PLLS,IFX,JFXETABLE,IFY,SMISC,5 !显示剪力ETABLE,JFY,SMISC,18PLLS,IFY,JFY!注意:beam单元的结果输出都是以单元坐标系输出的,且拉为正、压为负.前面我们已经知道,单元坐标系x轴就是i点指向j点,其他坐标可以根据整体坐标系推出.详细内容见ansys 帮助.(2)S hell单元对于shell单元应用的范围,ansys软件并没有强制规定,只是从字面上区分了薄壳和厚壳.我以前看过一本电子教案《仿真在线》,里面说一般规定壳体的主尺寸是厚度的10倍左右,都是可以用壳体来模拟的.一般高度与跨度之比(非与单元尺寸比较)<1/15,可以当作薄壳处理,>1/15 & <1/10,可以当作厚壳来处理.shell63是薄壳单元,他包含弯曲和薄膜效应,但是忽略横向剪切变形;shell43,shell143,shell181,shell91,shell93和shell99,都属于厚壳单元,不仅有弯曲、薄膜效应,他也包含了横向剪切效应.横向剪切被表示为整个厚度上的常剪切应变.这种一阶近似只适用于中等厚度壳体.线形分析时,如果不包含横向剪切应变,使用63,163单元;如果横向剪切变形重要,则遵守以下原则:均匀材料,使用43,93,143单元,复合材料使用91,99,181.我们土木工程中,一般利用shell43计算.!定义单元/prep7 !进入前处理et,1,shell43 !定义单元43号标号为1!定义材料属性mp,ex,1,2.55e7 !定义弹性模量(kn/m2)mp,nuxy,1,0.167 !定义泊松比mp,dens,1,2.5 !定义密度(KN/N*KG/M3)!定义墙体厚度!①等厚度板R,1,2 !1表示编号,2表示厚度(m)R,2,3Asel,s,,,1 !选中1号面Aatt,1,1,1 !mp,real,typeESIZE,0.2 !定义单元大小为0.2左右MSHAPE,0,2D !规定划分单元形状,0表示四边形(1表示三角形),2d表示划分面(3d表示划分体)MSHKEY,2 !指定是自由划分还是映射划分,2表示:尽量用映射划分,不符合要求就自动使用自由划分,具体参见ansys帮助的eshkey命令. AMESH,ALL !划分面单元.!!!注意:在网格剖面方面,最好全部用四边形,而且形状尽量规则、均匀!因为将来后处理内力提取的时候,提取出来的力和单元的大小有直接的关系。

显式有限元和隐式有限元

显式有限元和隐式有限元

显式算法和隐式算法,有时也称为显式解法和隐式解法,是计算力学中常见的两个概念,但是它们并没有普遍认可的定义,下面只是我的一些个人理解。

一、两种算法的比较1、显式算法基于动力学方程,因此无需迭代;而静态隐式算法基于虚功原理,一般需要迭代计算。

显式算法,最大优点是有较好的稳定性。

动态显式算法采用动力学方程的一些差分格式(如广泛使用的中心差分法、线性加速度法、Newmark法和wilson法等),不用直接求解切线刚度,不需要进行平衡迭代,计算速度快,时间步长只要取的足够小,一般不存在收敛性问题。

因此需要的内存也比隐式算法要少。

并且数值计算过程可以很容易地进行并行计算,程序编制也相对简单。

但显式算法要求质量矩阵为对角矩阵,而且只有在单元积分点计算尽可能少时速度优势才能发挥, 因而往往采用减缩积分方法,容易激发沙漏模式,影响应力和应变的计算精度。

静态显式法基于率形式的平衡方程组与Euler向前差分法,不需要迭代求解。

由于平衡方程式仅在率形式上得到满足,所以得出的结果会慢慢偏离正确值。

为了减少相关误差,必须每步使用很小的增量。

除了欧拉向前差分法外,其它的差分格式都是隐式的方法,需要求解线性方程组。

2、隐式算法隐式算法中,在每一增量步内都需要对静态平衡方程进行迭代求解,并且每次迭代都需要求解大型的线性方程组,这以过程需要占用相当数量的计算资源、磁盘空间和内存。

该算法中的增量步可以比较大,至少可以比显式算法大得多,但是实际运算中上要受到迭代次数及非线性程度的限制,需要取一个合理值。

二、求解时间使用显式方法,计算成本消耗与单元数量成正比,并且大致与最小单元的尺寸成反比,应用隐式方法,经验表明对于许多问题的计算成本大致与自由度数目的平方成正比,因此如果网格是相对均匀的,随着模型尺寸的增长,显式方法表明比隐式方法更加节省计算成本。

三、两种方法的应用范围:a)在求解动力学问题时,将方程在空间上采用有限元法(或其他方法)进行离散后,变为常微分方程组F=M(u)+C(u)+K(u)。

有限元各种时域计算方法

有限元各种时域计算方法

有限元各种时域计算方法有限元方法(FEM)是数值分析中一种常用的工程计算方法,用于求解连续介质的力学问题。

在时域情况下,FEM可以用于求解动力学问题,其中物体的响应随时间变化。

下面介绍几种常用的有限元时域计算方法:1. 爆炸分析方法(Explosion Analysis Method):用于模拟爆炸、冲击等快速载荷作用下的结构动力响应。

该方法将爆炸过程分解为多个离散时间步骤,并使用显式时间积分方法求解结构动力方程。

通过该方法可以得到结构的位移、速度、加速度等动态响应结果。

2. 频率域响应谱(Frequency Domain Response Spectrum):将时域问题转化为频域问题进行求解。

根据结构的固有频率和阻尼比,可以建立系统的频率响应函数,进而得到结构在特定载荷下的响应。

这种方法适用于大规模结构问题,可以有效地简化计算的复杂性。

3. 时间有限差分法(Time Finite Difference Method):该方法将时域问题转化为差分格式,用一系列离散时间步骤来近似连续时间。

通过在空间和时间上进行网格划分,可以利用差分格式求解结构动力方程。

这种方法对于线性和非线性问题都适用,并且可以实现高精度的模拟结果。

4. 显式时间积分法(Explicit Time Integration Method):该方法使用显式格式对结构动力方程进行时间积分,通过预测和修正的过程求解结构的动态响应。

显式时间积分法具有计算效率高的优点,适用于稳定性良好的问题,但在处理非线性和不稳定问题时可能出现数值耗散和不稳定现象。

5. 隐式时间积分法(Implicit Time Integration Method):与显式时间积分法相反,隐式时间积分法使用隐式格式进行时间积分,从而提高数值稳定性。

通过迭代求解非线性方程组,可以得到结构的准确动态响应。

隐式时间积分法对于非线性和不稳定问题的求解较为稳定,但计算效率较低。

以上是几种常用的有限元时域计算方法,每种方法都有各自的特点和适用范围。

有限元

有限元

有限元分析有限元分析是用较简单的问题代替复杂问题后再求解。

它将求解域看成是由许多称为有限元的小的互连子域组成,对每一单元假定一个合适的(较简单的)近似解,然后推导求解这个域总的满足条件(如结构的平衡条件),从而得到问题的解。

这个解不是准确解,而是近似解,因为实际问题被较简单的问题所代替。

由于大多数实际问题难以得到准确解,而有限元不仅计算精度高,而且能适应各种复杂形状,因而成为行之有效的工程分析手段。

有限元是那些集合在一起能够表示实际连续域的离散单元。

有限元的概念早在几个世纪前就已产生并得到了应用,例如用多边形(有限个直线单元)逼近圆来求得圆的周长,但作为一种方法而被提出,则是最近的事。

有限元法最初被称为矩阵近似方法,应用于航空器的结构强度计算,并由于其方便性、实用性和有效性而引起从事力学研究的科学家的浓厚兴趣。

经过短短数十年的努力,随着计算机技术的快速发展和普及,有限元方法迅速从结构工程强度分析计算扩展到几乎所有的科学技术领域,成为一种丰富多彩、应用广泛并且实用高效的数值分析方法。

在解偏微分方程的过程中, 主要的难点是如何构造一个方程来逼近原本研究的方程,并且该过程还需要保持数值稳定性.目前有许多处理的方法,他们各有利弊. 当区域改变时(就像一个边界可变的固体), 当需要的精确度在整个区域上变化, 或者当解缺少光滑性时, 有限元方法是在复杂区域(像汽车和输油管道)上解偏微分方程的一个很好的选择. 例如, 在正面碰撞仿真时, 有可能在"重要"区域(例如汽车的前部)增加预先设定的精确度并在车辆的末尾减少精度(如此可以减少仿真所需消耗); 另一个例子是模拟地球的气候模式, 预先设定陆地部分的精确度高于广阔海洋部分的精确度是非常重要的.2基本特点有限元方法与其他求解边值问题近似方法的根本区别在于它的近似性仅限于相对小的子域中。

20世纪60年代初首次提出结构力学计算有限元概念的克拉夫(Clough)教授形象地将其描绘为:“有限元法=Rayleigh Ritz法+分片函数”,即有限元法是Rayleigh Ritz法的一种局部化情况。

水下柔性结构流固耦合动力效应研究

水下柔性结构流固耦合动力效应研究

水下柔性结构流固耦合动力效应研究一、研究背景随着科技的不断发展,水下工程领域在船舶、海洋平台、海底隧道等诸多方面得到了广泛的应用。

然而由于水下环境的特殊性,如高压力、低温、盐度变化等,使得水下柔性结构在设计和施工过程中面临着诸多挑战。

为了提高水下柔性结构的可靠性和耐久性,研究其流固耦合动力效应显得尤为重要。

流固耦合是指物质在外力作用下发生的变形与流动现象,在水下柔性结构中,由于受到水流、波浪、潮汐等多种外部因素的影响,结构内部的应力分布和变形状态会发生动态变化。

因此研究水下柔性结构的流固耦合动力效应,有助于揭示其在不同工况下的响应特性,为优化设计提供理论依据。

近年来国内外学者对水下柔性结构的流固耦合动力效应进行了大量研究。

这些研究成果不仅为水下工程的设计提供了有力支持,还为实际工程应用提供了重要的参考价值。

然而现有研究成果主要集中在理论分析和数值模拟方面,对于实际工程中的具体问题解决能力有限。

因此进一步深入研究水下柔性结构的流固耦合动力效应具有重要的理论和实际意义。

1. 水下柔性结构的定义和分类梁式结构:梁式结构是最常见的一种水下柔性结构,主要包括横向梁和纵向梁。

横向梁主要用于承受横向水压力载荷,纵向梁则用于承受纵向拉力载荷。

这种结构形式简单、通用性强,适用于各种水下工程应用。

桁架结构:桁架结构是由许多相互支撑的杆件组成的空间框架结构。

在水下环境中,桁架结构可以通过调整杆件长度和间距来实现对受力状态的改变,从而适应不同的工况要求。

桁架结构具有较高的刚度和稳定性,但其制造工艺较为复杂。

索穹顶结构:索穹顶结构是一种以钢索为骨架,通过锚固在海底固定物上的穹顶状结构。

索穹顶结构具有良好的抗风蚀性能和抗冲击能力,同时能够承受较大的水压力载荷。

然而由于钢索的限制,索穹顶结构的刚度较低,且制造成本较高。

悬链网结构:悬链网结构是由一系列相互连接的链条组成的网状结构。

悬链网结构具有良好的柔韧性和抗拉强度,能够在受到外力作用时产生较大的形变,从而吸收部分能量,减小结构的应力集中。

有限元发展概况

有限元发展概况

有限元发展概况有限元发展概况⼀、有限元法介绍有限元法的基本思想是将结构离散化,⽤有限个容易分析的单元来表⽰复杂的对象,单元之间通过有限个节点相互连接,然后根据变形协调条件综合求解。

由于单元的数⽬是有限的,节点的数⽬也是有限的,所以称为有限元法(FEM,FiniteElementMethod)。

有限元法是最重要的⼯程分析技术之⼀。

它⼴泛应⽤于弹塑性⼒学、断裂⼒学、流体⼒学、热传导等领域。

有限元法是60年代以来发展起来的新的数值计算⽅法,是计算机时代的产物。

虽然有限元的概念早在40年代就有⼈提出,但由于当时计算机尚未出现,它并未受到⼈们的重视。

随着计算机技术的发展,有限元法在各个⼯程领域中不断得到深⼊应⽤,现已遍及宇航⼯业、核⼯业、机电、化⼯、建筑、海洋等⼯业,是机械产品动、静、热特性分析的重要⼿段。

早在70年代初期就有⼈给出结论:有限元法在产品结构设计中的应⽤,使机电产品设计产⽣⾰命性的变化,理论设计代替了经验类⽐设计。

⽬前,有限元法仍在不断发展,理论上不断完善,各种有限元分析程序包的功能越来越强⼤,使⽤越来越⽅便。

⼆、有限元法的孕育过程及诞⽣和发展⼤约在300年前,⽜顿和莱布尼茨发明了积分法,证明了该运算具有整体对局部的可加性。

虽然,积分运算与有限元技术对定义域的划分是不同的,前者进⾏⽆限划分⽽后者进⾏有限划分,但积分运算为实现有限元技术准备好了⼀个理论基础。

在⽜顿之后约⼀百年,著名数学家⾼斯提出了加权余值法及线性代数⽅程组的解法。

这两项成果的前者被⽤来将微分⽅程改写为积分表达式,后者被⽤来求解有限元法所得出的代数⽅程组。

在18世纪,另⼀位数学家拉格郎⽇提出泛函分析。

泛函分析是将偏微分⽅程改写为积分表达式的另⼀途经。

在19世纪末及20世纪初,数学家瑞雷和⾥兹⾸先提出可对全定义域运⽤展开函数来表达其上的未知函数。

1915年,数学家伽辽⾦提出了选择展开函数中形函数的伽辽⾦法,该⽅法被⼴泛地⽤于有限元。

瞬态接触_冲击问题的有限元分析

瞬态接触_冲击问题的有限元分析

李俊金咸定李维扬(上海交通大学船舶与海洋工程学院)摘要介绍瞬态接触2冲击问题的一种基于拉格朗日乘子法的有限元公式系统及其数值求解. 所建立的接触元能模拟接触的各个阶段, 如分离、粘附或有摩擦的滑移. 对于有摩擦的接触和释放给出了一种简洁的接触2释放条件, 为了更有效地求解耦合场运动方程, 引入了隐式2显式有限元算法. 该方法的有效性通过模拟核废料装运桶的冲击响应得到了验证.关键词固体; 接触2冲击问题; 非线性; 有限元分析中图法分类号O 313. 4固体之间或者固体与其他介质之间涉及接触2冲击的问题出现在许多领域中, 例如齿轮和滚筒的设计、核废料装运桶的冲击、飞机和汽车的防撞毁性、大型油轮液面的晃荡等. 由于接触边界条件的非线性, 实际感兴趣的问题几乎不能用解析法求解. 随着有限元法在50 年代末的出现, 数值求解接触问题已经被广泛地研究. Co n r y 和Se ireg 1 最早把接触问题处理成二次规划问题, 而C h a n 和T u b a2首先引入了一种普遍的滑移线方法. H u g h e s等3 建立的拉格朗日乘子法则对瞬态接触2冲击问题的求解作出了重大的贡献. 但文献3 中仅讨论了两维有摩擦无滑移的接触2冲击问题. 本文主要在某些方面修正和推广他们所做的工作.非线性瞬态问题的隐式2显式有限元分析1求解接触2冲击问题首先要建立支配变形体大运动的基本非线性方程. 本文的有限元离散法采用一种弱形式的支配偏微分方程, 它可由H am ilto n原理得到:∫t(∆E + ∆T - ∆W ) d t= 0 (1)t21式中: t1 和t2 表示两个时刻; ∆E , ∆T , ∆W分别为虚内功、虚动功和虚外功. 各符号的含义可参见文献4 .引入位移矢量并采用标准的有限元离散可得如下形式的离散运动方程:M ¨u + P (u,u )- R = 0 (2)近些年许多学者致力于研究瞬态耦合场数值求解的方法5 . 本文采用逐个单元的隐式2显式有限元算法求解上述非线性微分方程组. 假定M 为对角阵. 这一点与文献5 中的略有不同, 主要是因为对于冲击问题, 集总质量矩阵的性态优于相应的相容质量矩阵.接触2冲击问题的有限元模型的建立及求解2本节仅考虑两个物体的接触. 接触问题实际上是一个初边值问题. 可以把一个物体称作接触体, 另一个物体称作靶体.分析瞬态接触2冲击问题之所以极其困难, 主要是因为: 所考虑的物体的边界条件( 接触面积) 事先是未知的, 且接触面的大小和形状随时间而变化; 接触体之间的摩擦造成了随时间变化的粘附区域和滑移区域; 时间积分算法在接触面上失效. 前面两个问题的解决可通过适当地修正运动方程以反映非线性接触边界条件. 后一个问题的解决则通过修正时间积分算法以反映接触面上不连续的速度和加速度.收稿日期: 1997203215第一作者: 男, 1972 年生, 博士生. 邮编: 200030.考虑接触效应后, 式 (1) 变为2∑Α= 1t 2 ∫t{ (∆E Α + ∆T Α - ∆W Α) + ∆Q }d t = 0(3)1 式中: Α为相接触的物体; ∆Q 表示虚接触功. 采用拉格朗日乘子法强加接触协调性, 故取∆Q = ∫Σg d c(4)式中: Σ 为接触力矢量; g 为间隙矢量; c 为接触面积.对于普遍的两维有摩擦的接触2冲击问题, 采用双线性四边形单元来离散接触体和靶体. 接触面用接触体边界上的一组节点和靶体边界上的一组靶面来描述. 实际上任一接触体节点和任一靶面即构成了一个接触元. 离散接触泛函式 (4) 得N∫cΣg d c = ∑ΣΒg Β(5)Β= 1式中 N 为接触元的数目. 故基于拉格朗日乘子法的接触泛函的离散 式 ( 5) 将导得一混合变量的接触元, 这个接触元有 4 个节点, 每个节 点有两个自由度. 如图 1 所示.由上可知, 接触元的自由度是由如下矢量定义的广义位移: w = ( u 1 , Σ2 , u 3 , u 4 ) T , 相应的广义节点力为 F = (F 1 , g 2 , F 3 , F 4 ) T . 接触泛函 的离散及两个物体的标准离散可导得如下形式的运动方程:M ¨u + P (u , u α) + F (w ) = R ( t ) (6)·式中 P (u , u ) 的线性化可见文献4 . 下面仅考虑 F (w ) 的线性化. 不失一般性, 下面仅考虑一个接触元. 对 F (w ) 应用弗莱歇微分算子最后 导得图 1 接触元示意u1Σ2u 3 u4L0 0 L I - L- N 3L - N 4L0 N 3L 0 00 N 4L 0 - - c 0 0 1D F (w )w = K c w =(7), L =式中: K c 为接触切线刚度矩阵; I 为单位矩阵; N 3 , N 4 为节点 3 和 4 的插值函数. 需要指出的是: 对于无 摩擦的接触问题, F (w ) 的线性化可得到对称的接触刚度矩阵, 但对于涉及摩擦的滑移接触问题, F (w )的线性化得到的是非对称的接触刚度矩阵. 本文在适中滑移量的假设下, 引入小的完全滑移切线刚度矩阵, 从而得到了上面给出的对称的接触刚度矩阵.采用上节所述的隐式2显式算法即可求解系统运动方程 (6). 当运动方程的一个收敛解已经得到后,需要检查在当前时间步是否有节点进入或脱离接触. 如果有, 必须对它们应用接触释放条件, 以精确考虑接触和释放时产生的激波所造成的场变量的不连续性.3 瞬态接触和释放条件用隐式2显式算法求解两体接触问题时, 仅对位移场强加了接触和释放条件, 用此算法所得到的速度场和加速度场不满足不可穿透性约束. 文献6 中有一个简单的例子说明了时间积分算法在接触和释 放发生时的失效. 接触和释放条件计算的细节很复杂. 在文献3 中有详细的推导. 但文献3 中仅给出了有摩擦无滑移时的接触和释放条件. 下面将给出有摩擦的滑移时的接触和释放条件.不失一般性, 仅考虑一个接触元. 各符号见图 1, 接触和释放条件的速度的修正公式分别为d t ( 2 - 1 u αn ) N c γ (u αn ) 3 ( 3 - 11 3 1 - 1 Σn ) M + M - (u α3 ) + M 3+ M 1N 3N3M 1N 3N42 d t cγ c γ cγ c γ n = (8)M N c γN c γu αn ) M 1N 3 44 1 4 44 +( c γN c γM + (u αn ) M N c γ (u αn ) 4 4 - 1 1 4 1 - 1 (Σn ) 2 - 1+ - M 2接触和释放条件的加速度的修正公式分别为3+ M 1N 3N 3M 1N 3N 4(¨u 3 ) +M 3 (¨u 3 ) - + M 1N 3 (¨u 1 ) -n n c n cc c c(10) =M 1N 3 441 4 4( n )¨u 4 +4 4 -(n) - ¨u 1 4 1 -( n ) ¨u M + + c NcM N c N c M M N c(¨u 1 + ¨1 n ) (u n ) (11) = 接触和释放条件的接触力的修正公式分别为(Σ2 ) + = (Σ2 ) - - 0n M 1 [ (¨u 1 ) + - (¨u 1 ) - (12) (13)表示在当前时 ] n n n (Σ2 ) += n 上述诸式中: M 1 ,M 3 ,M 4 分别是节点 1, 3, 4 的集总质量; N Α= N Α Νc = N Α[ c Ν- 1 + Ν+c γ ( γ)( c ) /2 间步, 滑移的接触体节点的平均位置处的插值函数; 下标 n 为法线方向; ( ) - 1 为前一时间步的值; ( ) - 为当前时间步由隐式2显式算法计算所得到的值; ( ) + 为修正后的值. 与文献5 中的接触和释放条件的公 式相比, 本文的公式在形式上更简洁.实际上两维情况下的许多理论可直接推广到三维接触和冲击问题, 限于篇幅, 不作讨论.4 数值模拟本节介绍核废料装运桶与刚性地面的冲击的数值模拟, 并对结果进行讨论.为了安全输送装有核废料的桶, 需要研究具有有限初始速度的核装运桶 ( 有包壳的轴对称铅圆柱 体) 从 9 m 的高空自由落下与刚性地面发生冲击时的响应, 基本条件和数据列于表 1.表 1 核装运桶的基本条件和数据质量密度/( k g ·m - 3 ) 杨氏模量 /M P a屈服应力 /M P a应变硬化 率/M P a直径/c m厚度/c m长度/c m泊松比铅柱体不锈钢包壳30. 5- -0. 635 91. 491. 4 19. 5 1. 96×1040. 420. 33 29. 631018. 1 1. 91×1031. 13×10- 9 8. 00×10- 10模型的示意图如图 2 所示. 采用两维轴对称单元把铅柱体划分成 5×10 的网格 (显式单元) , 不锈钢包壳则被划分成 1×10 的网格 (隐式单元). 假定仅柱体的侧面被包壳覆盖, 滑移运动仅允许沿侧面的切 线方向发生. 还假定柱体的下表面保持与刚性地面接触 (没有反弹). 在界面上插入 10 个四节点接触元(隐式单元) , 其目的是模拟柱体2包壳之间可能发生的滑移.图 3 头部沉陷的时历图 2 核装运桶的冲击核装运桶在 0. 0 m s 时与地面发生冲击, 头部沉陷的时历如图 3 所示, 头部沉陷的最大值为6. 45c m , 发生于 5. 6 m s . 这与文献7 的计算结果基本吻合. 铅柱体的轴向压力分布示于图 4. 最大轴向应力 约为 18. 7 M P a . 包壳中的轴向应力状态则示于图 5. 结果显示了高频应力模式.在计算过程中可以发现柱体和包壳的分离, 这可能主要是因为柱体上下表面没有包壳, 柱体在径向图 4 10 m s 时铅柱体中轴向应力的分布图 5 10 m s 时包壳中轴向应力的分布参 考 文 献1 2 Co n ry T F , Se i reg A. A m a t h e m a t ica l p r o g ra mm ing m e tho d f o r de sign o f e l a st ic bo d ie s in c o n tac t . J A pp l M ech , 1971, 38: 387~ 392 C h an S K , T uba L S . A f in ite m e t ho d f o r c o n tac t p r o b le m s o f so lid bo d ies ——p a r t I : th eo ry and va lida t io n. I n t J M ech Sc i , 1971, 13: 615~ 625H ugh e s T J R , T ay lo r R L . A f in ite e le m en t m e tho d f o r la r ge d isp lace m en t c o n tac t and i m p ac t p r o b le m s . I n : B a th e K J , O den J T ed s .Fo r m u la t i o n s o f Com p u ta t io na l A lgo r ithm s , 1977. 469~ 495Zienk ie w icz O C , T ay l o r R L . T h e F in ite E l e m en t M eho d . L o n do n : M cG r a w 2H ill , 1991. B e ly t sch ko T , H ugh e s T J R ( ed s ). Com p u ta t io na l m e tho d s f o r t ran sien t ana ly sis in c om p u ta t i o na l m e tho d s in m ech n ice s . V o l . 1, N o r th 2H o lland, A m ste rda m , 1983.T ay lo r R L , P ap adopo u lo s P. A f in ite e le m en t m e tho d f o r dyna m ic c o n tac t p r o b le m s . I n: O na te E , P e r iaux J , Sa m ue lsso n A ed s . F in iteE le m en t s in th e 90’s, 1991. 212~ 222Ku lak R F . A dap t ive c o n tac t e le m en t s f o r th ree 2d i m en sio na l exp lic it t r an sien t ana ly sis . Com p u te r M e th in A pp l M ech and E ng , 1989,72: 125~ 1513 4 56 7 F in ite E l em e n t A na l ys i s fo r T ra ns ie n t C o n ta c t 2I m p a c t P ro b lem sL i J u n J in X ia n d in g L i W e i y a n g(Co llege o f N ava l A rch itec t u re an d O cean E n g i n e e r i n g , Sh a n g h a i J i ao to n g U n i ve r s ity )A b s t r a c t A f i n i te e lem en t fo rm u l a t i o n b a s ed o n lag r an g e m u lt i p li e r tech n i qu e s fo r t ran s i en t co n t ac t 2i m p ac t p ro b lem s is p re sen ted . T h e deve l op ed co n tac t e le m en t h a s th e ab ility to s i m u la te each o f th e co n t ac t s tage su ch a s sep a r a t i o n , s t i ck i n g an d f r i c t i o n a l s li d i n g . Fo r th e f r i c t i o n a l co n t ac t s an d re l ea s e s , a co n c i se fo rm u la t i o n o f co n t ac t 2re l ea s e co n d i t i o n s is g i ven . In o rde r to so lve co u p i ed 2f i e l d equ a t i o n s o f m o t i o n m o re eff i c i en t ly , th e i m p li c i t 2exp li c it f i n i te e l em en t a l go r ithm is a lso i n t ro d u c ed . T h e vad i lity o f th e tech n i qu e is dem o n s t r a t ed b y s i m u la t i n g i m p a c t re s po n s e o f ax isymm e t r i c sp e n t n u c l ea r fu e l 2sh i p 2p i n g ca s k . so li d s ; co n t ac t 2i m p a c t p ro b l em s ; no n 2li n e a r ; f i n i te e lem en t an a l y s isKe y w o rd s。

Abaqus隐式和显式求解的区别和应用

Abaqus隐式和显式求解的区别和应用

Abaqus隐式和显式求解的区别和应⽤Abaqus有限元计算要使⽤的求解器类型:隐式还是显式?求解器类型会影响求解的⽅程组、某些单元的可⽤性、运⾏时间,甚⾄是否获得收敛。

本⽂将解释Abaqus中可⽤的两个求解器之间的区别,它们的优缺点以及何时选择哪种求解器。

有什么区别?当谈论FEA中的隐式或显式时,谈论的是⽤于时间增量的算法,在两种情况下,模型的状态都是在多个时间点计算,新状态是根据旧状态计算的。

使⽤显式算法,可以直接从当前状态下的可⽤数据计算新状态,基本上是⼀种推断。

⽽使⽤隐式算法,不能直接从旧状态计算新状态,必须求解⼀个耦合的⽅程组,这需要⾮线性解算法,通常是Newton-Raphson⽅法。

计算成本和时间增量⼤⼩单个显式增量的计算成本很⼩,所需的所有信息都可⽤,计算简单快捷。

但是,时间增量不能太⼤,因为这样会变得不稳定,解上的误差会呈指数级增长。

如果我们将其视为外推,远远超出已知范围的推断往往会给出错误的结果,特别是如果错误有机会相加。

可以使⽤的最⼤时间增量由 Abaqus ⾃动计算,称为稳定时间增量,稳定的时间增量随着单元尺⼨的减⼩、密度的降低和刚度的提⾼⽽减⼩。

具有最⼩稳定时间增量的单元决定了整个分析的时间增量,因此,单个形状不良的单元会极⼤地影响仿真时间。

在整个分析过程中,稳定时间增量通常近似恒定,因此,⼀旦第⼀个增量完成,可以估计运⾏分析所需的时间。

隐式增量的计算成本要⼤得多,因为需要求解⼀个⽅程组,对于⾮线性分析,这甚⾄需要多次完成。

因此,解决单个隐式增量将需要更多的时间和内存。

另⼀⽅⾯,时间增量⼤⼩不受隐式算法稳定性问题的限制,通常允许⽐显式算法更⼤的时间增量。

Abaqus/Standard(隐式求解器)中的增量⼤⼩通常由 Abaqus 根据Newton-Raphson⽅案收敛的难易程度⾃动确定。

越是⾮线性,找到收敛解决⽅案的计算成本就越⾼。

图 1:显式和隐式求解器的时间增量。

lsdyna 动力计算原理

lsdyna 动力计算原理

lsdyna 动力计算原理
LS-DYNA是一种通用的显式和隐式有限元分析软件,主要用于
模拟动态加载、非线性变形、多物理场耦合等复杂问题。

在LS-DYNA中,动力计算原理涉及以下几个方面:
1. 有限元方法,LS-DYNA采用有限元方法来离散模拟对象,并
通过求解大型矩阵方程来计算结构的动力响应。

有限元方法将结构
划分为有限数量的单元,每个单元内的物理行为由一组局部方程描述,通过组装这些单元的方程,可以得到整个结构的全局方程。

2. 显式求解,LS-DYNA采用显式时间积分方法来求解结构的动
力响应。

这意味着在每个时间步内,软件会显式地计算结构的位移、速度和加速度,而不需要解耦结构的刚度矩阵和质量矩阵。

这种方
法适用于高速动态加载和非线性变形的情况。

3. 材料模型,LS-DYNA支持多种材料模型,包括线性弹性、塑性、损伤、断裂等模型。

在动力计算中,材料模型的选择对于准确
预测结构的动态响应至关重要。

4. 质量和边界条件,在动力计算中,需要准确描述结构的质量
分布和边界条件。

LS-DYNA允许用户定义不同类型的质量和边界条件,以便更准确地模拟真实的动力加载情况。

总的来说,LS-DYNA的动力计算原理涉及有限元方法、显式求解、材料模型和质量边界条件等方面,通过综合考虑这些因素,可以对结构在动态加载下的响应进行准确模拟。

有限元分析基础-隐式与显式算法的区别

有限元分析基础-隐式与显式算法的区别

有限元分析基础- 显式与隐式算法的区别所谓显式和隐式,是指求解方法的不同,即在数学上的计算方法不一样,是两种不同针对时间的积分方法。

显式求解是对时间进行差分,不存在迭代和收敛问题,最小时间步取决于最小单元的尺寸。

过多或过小的时间步都会导致求解时间非常漫长,但总能给出一个计算结果。

隐式求解和时间无关,采用的是牛顿迭代法(线性问题就直接求解线性代数方程组),因此存在一个迭代收敛问题,不收敛就得不到结果。

在某些情况下,显式算法的计算效率远高于隐式算法,尤其是在多处理器并行运算的场景下,对于自由度较大的三维结构,显式算法可能具有较高的计算效率。

然而,对于自由度较小的二维结构,隐式算法可能更适合,因为它在每个增量步内不需要进行迭代,从而减少了计算时间。

总结来说,隐式方法适用于需要高计算精度和稳定性的场景,而显式方法则适用于需要高计算效率的场景。

隐式与显式动力学的区别-弹性动力学有限元基本解法

隐式与显式动力学的区别-弹性动力学有限元基本解法

1.1. 弹性动力学有限元基本解法结构系统的通用运动学方程为:tR KU U C U M =++ (1) 求解该动力学振动响应主要有三类方法:(1)时域法(2)频域法(3)响应谱法时域法又可分为:(1)直接积分法,(2)模态叠加法。

直接积分法又可分为中心差分法(显式),Wilson θ(隐式)法以及Newmark (隐式)法等。

本文介绍中心差分法(显式)与Newmark (隐式)法。

1 中心差分法(显式)假定0,1t ,2t ,…,n t 时刻的节点位移,速度与加速度均为已知,现求解)(t t t n ∆+时刻的 结构响应。

中心差分法对加速度,速度的导数采用中心差分代替,即为:)2(12t t t t t t U U U t U ∆+∆-+-∆= )(21t t t t t U U tU ∆-∆+-∆= (2) 将(2)式代入(1)式后整理得到tt t R U M ˆˆ=∆+ (3) 式(3)中C tM t M ∆+∆=211ˆ2 t t t t t U C tM t U M t K R R ∆-∆-∆-∆--=)211()2(ˆ22 分别称为有效质量矩阵,有效载荷矢量。

R ,M ,C ,K 为结构载荷,质量,阻尼,刚度矩阵。

求解线性方程组(3),即可获得t t ∆+时刻的节点位移向量t t U ∆+,将t t U ∆+代回几何方程与物理方程,可得t t ∆+时刻的单元应力和应变。

中心差分法在求解t t ∆+瞬时的位移t t U ∆+时,只需t t ∆+时刻以前的状态变量t U 和t t U ∆-,然后计算出有效质量矩阵M ˆ,有效载荷矢量tR ˆ,即可求出t t U ∆+,故称此解法为显式算法。

中心差分法,在开始计算时,需要仔细处理。

t =0时,要计算t U ∆,需要知道t U ∆-的值。

因此应该有一个起始技术,因而该算法不是自动起步的。

由于0U ,0U ,0U 是已知的,由t =0时的(2)式可知: 02002U t U t U U t ∆+∆-=∆-中心差分法中时间步长t ∆的选择涉及两个方面的约束:数值算法的稳定性和计算时间。

catia有限元分析简述

catia有限元分析简述

前言运用固体力学理论(包括结构力学、弹性力学、塑性力学等)对结构进行强度和刚度分析,是工程设计的重要内容之一。

随着科学技术的进步和生产的发展,工程结构的几何形状和载荷情况日益复杂,新的材料不断出现,使得寻找结构分析的解析解十分困难,甚至不可能,因而人们转而寻求近似解。

1908年,W.Ritz提出一种近似解法,具有重要意义。

它利用带未知量的试探函数将势能泛函近似,对每一个未知量求势能泛函的极小值,得到求解未知量的方程组。

Ritz法大大促进了弹性力学在工程中的应用。

Ritz法的限制是试探函数必须满足边界条件。

对于几何形状比较复杂的结构来说,寻找满足整个边界条件的试探函数也非易事。

1943年,R.Couran对Ritz法做了极其重要的推广。

他在求解扭转问题时,将整个截面划分为若干个三角形区域,假设翘曲函数在各个三角形区域内做近似线性分布,从而克服了以前Ritz法要求整体近似函数满足全部边界条件的困难。

Couran这样应用Ritz法与有限元法的初期思想是一致的。

但是这种近似解法要进行大量数值计算,在当时还是个难题。

因此,未能得到发展。

有限单元法是采用计算机求解数学物理问题的一种数值计算近似方法。

它发源于固体力学,后迅速扩展到流体力学、传热学、电磁学、声学等其它物理领域。

固体力学有限元法的理论依据,从发展历史看,主要有三种途径,即结构矩阵法、变分法和加权余量法。

整个计算过程是泰国编制好的程序在电子计算机上自动进行。

它具有极大的通用性,在程序功能范围内,只要改变输入的数据,就可以求解不同的工程实际问题。

这种解法完全改变了解析法中针对一种实际问题寻找一种解法的局限性。

在1946年电子计算机诞生以后,首先采用它进行数值计算的是杆系结构力学。

它的理论依据是由结构力学位移法和力学演变成的矩阵位移法和矩阵力学,统称为结构矩阵法。

它采用矩阵代数运算,不仅能使算式书写简明,而且编制计算机程序非常方便。

结构矩阵法的力学概念清楚,全部理论公式按结构力学观点讲都是准确的,仅在数值计算过程中,由于计算机存储位数的限制,造成舍入误差。

接触问题的显式与隐式有限元方法研究

接触问题的显式与隐式有限元方法研究

接触问题的显式与隐式有限元方法研究接触问题是固体力学中的一个重要问题,它在机械、航空航天、汽车制造等领域都有着广泛的应用。

有限元方法是解决接触问题的一种有效手段,其中显式方法和隐式方法是两种常见的有限元方法。

显式方法指的是在求解过程中,通过显式地列出接触条件和摩擦力,从而求解接触问题。

显式方法的优点是计算过程简单,易于实现,且能够考虑接触面间的摩擦力。

但是,显式方法也存在一些缺点,例如难以处理复杂的接触形状和非线性材料性质,此外,显式方法还需要对摩擦系数进行假设,否则无法求解。

隐式方法指的是在求解过程中,通过隐式地考虑接触条件和摩擦力,从而求解接触问题。

隐式方法的优点是能够自适应地计算接触面间的摩擦力,且能够处理复杂的接触形状和非线性材料性质。

但是,隐式方法的计算过程比较复杂,需要对材料性质和接触面进行精确的建模,此外,隐式方法还需要对摩擦系数进行假设,否则无法求解。

为了验证显式方法和隐式方法的有效性,我们通过实例进行了计算。

我们选取了一个简单的接触问题,即两个平行的金属板之间的接触问题,其中一个金属板的长度为 100mm,宽度为 50mm,厚度为 1mm,另一个金属板的长度为 100mm,宽度为 50mm,厚度为 2mm。

两个金属板之间的接触面为 50mm×50mm。

我们分别采用了显式方法和隐式方法进行计算,并比较了计算结果。

显式方法的计算过程如下:我们首先建立有限元模型,然后列出接触条件和摩擦力,最后求解得到接触压力和摩擦力。

显式方法的计算结果如下:接触压力为 999.5N,摩擦力为 29.1N。

隐式方法的计算过程如下:我们首先建立有限元模型,然后采用隐式方法计算接触面间的摩擦力,最后求解得到接触压力和摩擦力。

隐式方法的计算结果如下:接触压力为 999.5N,摩擦力为 29.1N。

通过比较计算结果,我们可以发现,显式方法和隐式方法都能够准确地计算接触问题,但是,隐式方法的计算结果更加准确,能够考虑接触面间的摩擦力,而显式方法需要对摩擦系数进行假设,否则无法求解。

有限元显式时间步长的计算

有限元显式时间步长的计算

有限元显式时间步长的计算
有限元显式时间步长的计算指的是使用显式时间离散方法进行有限元计算。

通常而言,有限元分析的时间离散过程中会使用显式或隐式时间离散方法。

显式方法在计算过程中直接使用当前时刻的信息进行计算,因此需要满足时间步长的稳定性条件。

隐式方法则将未知量推迟到下一个时间步长进行计算,因此相对而言更为稳定。

在使用有限元方法进行求解时,有限元分析计算通常采用显式时间离散方法进行求解。

这种方法在计算复杂结构的动力学或热力学问题中具有很好的性能。

具体来说,有限元显式时间步长的计算过程通常包括以下几个步骤:
1. 构造有限元形式的时间离散格式,即将时刻$t^{n}$的解$u^{n}$表示为它的节点值$u_{i}^{n}$的集合。

2. 将时间区域$(0,T)$离散化为$N$个时间步长。

在每个时间步长上,采用有限元方法求解得到相应时刻下的节点值
$u_{i}^{n}$。

3. 使用显式时间离散格式更新下一个时间步长的节点值
$u_{i}^{n+1}$。

4. 重复步骤2和步骤3,直到计算结束。

需要注意的是,在显式时间离散格式中,由于每个时刻的解都由前一个时刻的解和该时刻的右端项组合而成,因此需要满足一定的时间步长限制条件。

具体来说,时间步长$\Delta t$需要满足以下稳定性条件:
$$\Delta t \leq \frac{2}{c_{max}}\frac{h}{u_{max}}$$
其中$c_{max}$是方程中的速度或波速极大值,$h$是网格尺寸,$u_{max}$是解的极值。

这个条件保证了显式时间离散格式的稳定性和收敛性。

abaqus显式方法和隐式方法

abaqus显式方法和隐式方法

abaqus是一个常用的有限元分析软件,广泛应用于工程领域。

在abaqus中,有两种常用的有限元分析方法,分别是显式方法和隐式方法。

这两种方法在不同的工程场景下有着不同的适用性和特点,本文将对abaqus显式方法和隐式方法进行介绍和比较。

一、abaqus显式方法1. 简介abaqus显式方法又称为动力显式有限元法,适用于求解瞬态非线性动力学问题。

其基本原理是通过将时间步长划分为相对较小的子步长,在每个子步长内计算结构的响应。

显式方法的特点是计算速度快,适用于求解高速碰撞、爆炸等暂态动力学问题。

2. 适用场景abaqus显式方法适用于以下工程场景:- 高速动态加载问题:如车辆碰撞、子弹撞击等。

- 爆炸冲击问题:如爆破工程、军事领域的冲击载荷等。

- 材料动态行为问题:如弹性-塑性材料在高速加载下的动态响应。

3. 优缺点优点:- 计算速度快,适用于大规模复杂结构的动力学仿真。

- 对动态非线性效应的处理能力强。

- 需要选择合适的时间步长和子步长,计算稳定性受到限制。

- 适用范围受到限制,不能用于稳态问题的分析。

二、abaqus隐式方法1. 简介abaqus隐式方法又称为静态隐式有限元法,适用于求解稳态和瞬态非线性静力学问题。

其基本原理是通过建立增量形式的非线性方程组,采用迭代算法求解非线性方程组的平衡解。

隐式方法的特点是计算稳定,适用于求解稳态和缓变动力学问题。

2. 适用场景abaqus隐式方法适用于以下工程领域:- 结构强度分析:如钢结构、混凝土结构的承载能力分析。

- 热-机-结构耦合问题:如热载荷下的构件稳定性分析、变温环境下的材料性能分析。

- 持续动态加载问题:如地震、风载等动态荷载下的结构响应分析。

3. 优缺点优点:- 计算稳定性好,适用于求解长时间的稳态和瞬态非线性问题。

- 对非线性效应的处理能力强。

- 计算速度相对较慢,在处理大规模结构的动力学仿真时耗时较长。

- 算法求解过程复杂,需要对非线性迭代收敛条件进行调整。

显式求解方法和隐式求解方法对比

显式求解方法和隐式求解方法对比

采用有限元方法开展结构的动力学分析最终归结为求解离散后的常微分方程组tR KU U C U M =++ 。

在时域内求解该方程最常用的方法是直接积分法,而又根据求解过程中是否需要迭代求解线性方程组,将直接积分法分为隐式积分方法和显式积分方法两类。

隐式积分法认为t+Δt时刻系统的状态不仅与t时刻状态有关,且与t+Δt时刻某些量有关。

因此隐式算法是根据t n 及t n-1...时刻体系的物理量值建立关于以t n+1时刻物理量为未知量的线性方程组,通过求解方程组确定t n+1时刻的物理量(常用的方法有线性加速度法、常平均加速度法、Newmark 方法、Wilson-θ法、Houbolt 方法等)。

而显式积分法认为t+Δt时刻系统的状态仅与t时刻状态有关可,因此可由t n 及t n-1...时刻体系的物理量值直接外推t n+1时刻物理量值(如中心差分法),不需要求解线性方程组,实现了时间离散的解耦。

两种算法的比较 :(1)隐式算法隐式算法基于虚功原理,要迭代计算。

隐式算法在每一增量步内都需要对静态平衡方程进行迭代求解,并且每次迭代都需要求解大型的线性方程组,这一过程需要占用相当数量的计算资源、磁盘空间和内存。

理论上在这个算法中的增量步可以很大,但是实际运算中上要受到接触以及摩擦等条件的限制。

随着单元数目的增加,计算时间几乎呈平方次增加。

由于需要矩阵求逆以及精确积分,对内存要求很高。

隐式算法的不利方面就是收敛问题不容易解决,且在开始起皱失稳时,在分叉点处刚度矩阵出现奇异。

(2)显式算法显示算法基于动力学方程,无需迭代,包括动态显式和静态显式算法。

动态显式算法采用动力学方程的中心差分格式,不用直接求解切线刚度,不需要进行平衡迭代,计算速度快,也不存在收敛控制问题。

该算法需要的内存也比隐式算法要少。

数值计算过程可以很容易地进行并行计算,程序编制也相对简单。

它也有一些不利方面。

显式算法要求质量矩阵为对角矩阵,而且只有在单元级计算尽可能少时速度优势才能发挥, 因而往往采用减缩积分方法,容易激发沙漏模式,影响应力和应变的计算精度。

显式有限元方法在时程分析和IDA分析中的应用

显式有限元方法在时程分析和IDA分析中的应用

Liu 梁单元,支撑采用杆单元(梁单元中的第3选项),在梁柱节点 处采用Mass 单元输人楼层质量。

模型中梁、柱划分为5个单元,支撑划分为1个单元。

计算中 考虑几何非线性与材料非线性。

地震记录采用El Centro 波,地震 波峰值加速度按照抗震设防烈度为8度的罕遇地震取值,即将地 震波的加速度峰值调整为〇. 4g(400 cm/S2),持时取地震记录的 前30 S 。

隐式计算的时间步长取为At = 0. 02 S ,显式计算积分步 长由软件根据中心差分法稳定条件自动选取。

计算结果对比:两种方法算出的时程曲线如图2所示,图2中时程曲线前10 s 用来消除重力荷载造成的振荡,地震波作用从10 s 开始,从中可以看出:1)两种方法的计算结果具有很好的一致性,时程曲线吻合非常好。

在最大顶点水平位移出现的t =22.14 s ,显式、隐式方法算得的位移值分别为50. 4 mm ,51. 1 mm ,仅仅相差1.3%,对于强烈 地震下的时程分析这是一个很小的误差,说明显式方法在时程分 析中能得到令人满意的计算精度。

0 10 20 30 40 0 10 20 30 40时间/s时间/sa )显式方法时程曲线b )隐式方法时程曲线图2显式、隐式方法时程分析结果对比2)显式方法有明显的速度优势。

本算例在同一计算平台上(CPU : Imer (R ) Core (TM )2 Duo P 8700 @ 2.53 GHz ,4G 内存,32 位winxp 系统)上使用显式、隐式方法所消耗的时间分别为33 s 和17 min ,显式方法的速度大致相当于隐式方法的30倍。

时程分 析非常耗时,显式方法在计算速度方面的巨大优势使其可以大显 身手。

2中心支撑钢框架IDA 倒塌分析本算例说明显式方法在IDA 倒塌分析中的特点。

采用隐式有限元进行IDA 分析,当地震波强度调整到倒塌临界点时,结构 的总刚矩阵行列式趋于零,水平位移趋于无限大,计算发散,于是认为结构倒塌[6]。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

按照计算每一时刻动力反应是否需要求解线性方程组,可将直接积分法分为隐式积分方法和显式积分方法两类。

隐式积分法是根据当前时刻及前几时刻体系的动力反应值建立以下一时刻动力反应值为未知量的线性方程组,通过求解方程组确定下一时刻动力反应。

隐式方法的研究和应用由来已久,常用的方法有线性加速度法、常平均加速度法、New mark方法、Wilson-θ法、Houbolt 方法等。

显式积分法可由当前时刻及前几时刻的体系动力反应值直接外推下一时刻的动力反应值,不需要求解线性方程组,实现了时间离散的解耦。

解方程组一般占整个有限元求解程序耗时的70%左右,因此,这一解耦技术对计算量的节省是可观的。

隐式方法大部分是无条件稳定的,显式方法为条件稳定。

显式方法的稳定性可以按满足精度要求的空间步距确定满足数值积分稳定性要求的时问步距来实现。

显式方法受条件稳定的限制,时间积分步长将取得较小,但计算经验表明,对于一些自由度数巨大且介质呈非线性的问题,显式法比隐式法所需的计算量要小得多。

因此,随着所考虑问题复杂性的增加,显式积分法得到重视。

对于显式与隐式有限元的理解
关键字: 有限元显式隐式
显式算法和隐式算法,有时也称为显式解法和隐式解法,是计算力学中常见的两个概念,但是它们并没有普遍认可的定义,下面只是我的一些个人理解。

一、两种算法的比较
1、显式算法
基于动力学方程,因此无需迭代;而静态隐式算法基于虚功原理,一般需要迭代计算。

显式算法,最大优点是有较好的稳定性。

动态显式算法采用动力学方程的一些差分格式(如广泛使用的中心差分法、线性加速度法、Newmark法和wilson法等),不用直接求解切线刚度,不需要进行平衡迭代,计算速度快,时间步长只要取的足够小,一般不存在收敛性问题。

因此需要的内存也比隐式算法要少。

并且数值计算过程可以很容易地进行并行计算,程序编制也相对简单。

但显式算法要求质量矩阵为对角矩阵,而且只有在单元积分点计算尽可能少时速度优势才能发挥, 因而往往采用减缩积分方法,容易激发沙漏模式,影响应力和应变的计算精度。

静态显式法基于率形式的平衡方程组与Euler向前差分法,不需要迭代求解。

由于平衡方程式仅在率形式上得到满足,所以得出的结果会慢慢偏离正确值。

为了减少相关误差,必须每步使用很小的增量。

除了欧拉向前差分法外,其它的差分格式都是隐式的方法,需要求解线性方程组。

2、隐式算法
隐式算法中,在每一增量步内都需要对静态平衡方程进行迭代求解,并且每次迭代都需要求解大型的线性方程组,这以过程需要占用相当数量的计算资源、磁盘空间和内存。

该算法中的增量步可以比较大,至少可以比显式算法大得多,但是实际运算中上要受到迭代次数及非线性程度的限制,需要取一个合理值。

二、求解时间
使用显式方法,计算成本消耗与单元数量成正比,并且大致与最小单元的尺寸成反比,应用隐式方法,经验表明对于许多问题的计算成本大致与自由度数目的平方成正比,因此如果网格是相对均匀的,随着模型尺寸的增长,显式方法表明比隐式方法更加节省计算成本。

有限元方法的显式及隐式
对时间积分的两种算法:
隐式方法:
大多数的有限元分析软件都是采用隐式方法,这种方法收敛速度较快。

cn+1=an+bn
优点是计算量比较小
缺点是有累计误差
n+1个时间步的量不可以由第n个时间步的量直接求得,称为隐式!
显式方法:
显示积分方法一般用在高度非线性有限元分析,如碰撞、爆炸、冲击等。

dyna 等软件一般采用显示有限元法。

这种方法的收敛较慢,为了保证收敛一般要取较短的时间步长。

关于显式积分与隐式积分的内容可以看一下《数值分析》中关于椭圆型、抛物线型或双曲型微分方程的差分方法等内容。

例如:
an+1+bn+1=cn
bn+1+cn+1=an
an+1+cn+1=bn
缺点是计算量比较大,需要通过方程组求解
优点是没有累计误差。

总之:用比较通俗的话说: 显式就是可以直接通过自变量求得因变量的解,自变量和因变量可以分离在等式的两侧;
隐式正好相反,因变量与自变量混和在一起,不能进行分离.
显式解法里,没有刚度矩阵的说法?!
显式解法基于牛顿第二定律,F=M*acce,
其中F由上一时步的外载,内力确定;
由acce --> velocity -->disp, 也就可相应求解应力,应变值了
根据我的理解,是这样的:
(1)显式算法包括动态显式和静态显式算法。

动态显式算法的最大优点是有较好的稳定性。

动态显式算法采用动力学方程的中心差分格式,不用直接求解切线刚度,不需要进行平衡
迭代,计算速度快,也不存在收敛控制问题。

该算法需要的内存也比隐式算法要少。

数值计算过程可以很容易地进行并行计算,程序编
制也相对简单。

它也有一些不利方面。

显式算法要求质量矩阵为对角矩阵,而且只有在单元积分点尽可能少时速度优势才能发挥, 因而往
往采用减缩积分方法,容易激发沙漏模式,影响应力和应变的计算精度。

动态隐式法还有一个重要特点是:对成形过程的仿真需要使用
者正确划分有限元网格和选择质量比例参数、速度和阻尼系数。

静态显式法基于微分形式的平衡方程组与Euler前插公式,不需要迭代求解。

由于平衡方程仅在微分形式上得到满足,所以得出的结果
会慢慢偏离正确值。

为了减少相关误差,必须每步使用很小的增量,通常一个仿真过程需要多达几千步。

由于不需要迭代,所以这种方
法稳定性好,但效率低。

(2)隐式算法
静态算法也是解决金属成形问题的一种方法。

在静态隐式算法中,在每一增量步内都需要对静态平衡方程而迭代求解。

理论上在这个算
法中的增量步可以很大,但是实际运算中上要受到接触以及摩擦等条件的限制。

随着单元数目的增加,计算时间几乎呈平方次增加。


于需要矩阵求逆以及精确积分,对内存要求很高。

隐式算法的不利方面还有收敛问题不容易得到解决以及当开始起皱失稳时,在分叉点
处刚度矩阵出现奇异。

另有一种静态隐式大增量步软件,也属于静态隐式算法,做出了某些改进,如在一些特殊接触条件处理上采用大增量时步,弯曲与拉伸
变形的非耦合求解算法,高精度的自适应网格划分等等。

这些专用于金属薄板成形的特征有时显得非常有效,但在某些方面不会那么准
确。

例如,它不能精确模拟接触和脱离接触的过程,无法有效预测起皱失稳。

相关文档
最新文档