1、2018/12/16,1,内容提要,引言 非线性的来源 非线性问题的求解 求解非线性有限元问题应当考虑的因素,2018/12/16,2,第一章 引言,2018/12/16,3,1.1 线性有限元回顾,2018/12/16,4,前提:材料处于弹性状态,应力-应变关系是线性的。位移和应变是微小的。因此,2018/12/16,5,有限元法实施的步骤,一、建立有限元模型离散化,前处理,二、单元分析位移假设,几何方程,物理方程,等效节点力,单元刚度矩阵,三、总刚度矩阵、总载荷列阵、引入支撑最小势能原理、 集成技术(叠加法、对号入座),四、求解、结果分析,后处理,2018/12/16,6,一、建立有限元模
2、型,注意:在有限元模型中,单元与单元、单元与支座只在节点相连,载荷作用在节点上。,2018/12/16,7,1.几何划分,二维 线 (直线 折线 曲线)三维 面 (平面 折面 曲面) 由边界和切割线(面)形成有限元网格,使得连续域成为离散域每一小块:单元(element)结点:场变量在该点的值为未知量。结点位移:位移元的基本未知量。,2018/12/16,8,一般来说,商业软件提供一系列有限单元模型,如:,一维,二维,三维,2018/12/16,9,以及,一系列结构单元模型,如:,杆单元,梁单元,管单元,膜单元,平面单元,2018/12/16,10,通常把三维实体划分成4面体或6面体单元的网格
3、,平面问题划分成三角形或四边形单元的网格。,2018/12/16,11,平面问题的四边形单元和三角形单元划分,2018/12/16,12,空间问题的四面体形单元划分,2018/12/16,13,空间问题的六面体形单元划分,2018/12/16,14,实际物体厚度较小时,选择板壳单元为好,2018/12/16,15,例如,对一个实际物体按照3维问题进行离散处理,2018/12/16,16,当实际物体太大,可以考虑先离散成子结构(Substructures),用超级单元处理,2018/12/16,17,3、形成网格1) 合理的疏密,变化剧烈的地方可密,变化不剧烈的地方 可疏。2) 合理的过渡(协调
4、) *保证界面连续,*单元形态的合理性。,2、节点和单元 1)结点坐标 2)结点两种编号 3)结点位移为位移元的基本未知量,2018/12/16,18,4、单元分析与单元刚度矩阵以平面三角形常应变单元为例,注意:一个单元的节点数目( 3);一个节点的自由度数目(可能位移,2)。,19,形成单元刚度矩阵从虚功原理出发:给节点一组虚位移,单元内即产生相应的虚应变,那么,实际节点力在虚位移上所做的功就等于实际应力在虚应变上产生的应变能。,Ke =,2018/12/16,20,5、等效节点载荷等效原则:等效节点载荷在节点虚位移上所做的功等于实际载荷在相应虚位移上所做的功。,2018/12/16,21,
5、6、形成总刚度矩阵、总载荷列阵、引入边界条件,按照总刚度矩阵形成方式以及代数方程德求解方案不同,可以分为: (1)直接求解(稀疏带状求解和波阵求解); (2)迭代求解。,采用稀疏带状求解器,要求节点号排列规律应当优化,以使任一单元内的节点号的差值最小。,采用波阵求解器,要求单元号排列规律应当优化,以使任一节点所连接的单元号的差值最小。,2018/12/16,22,7、后处理,对求解的数据用等值线和云图表示,2018/12/16,23,1.2 线性与非线性分析,2018/12/16,24,线性分析施加的载荷和系统响应间为线性关系 在边界条件相同的情况下,满足(载荷)叠加原理。(载荷)叠加原理对线
6、性动力学也适用),Ku=P,定系数,2018/12/16,25,非线性分析 结构问题的非线性是指结构的刚度随其变形而改变 实际上所有的物理结构均是非线性的。(线性分析只是一种方便的近似,这对设计来说通常是足够精确的。)必须做非线性分析的例子:冲压等加工过程、碰撞分析以及轮胎和发动机垫圈一类的橡胶部件的分析问题等。,2018/12/16,26,非线性分析,2018/12/16,27,非线性分析,每种载荷都必须作为独立的分析进行定义及求解。,K(u, , , )u=P,不定系数,.,u=K(u, , , )-1 P,.,2018/12/16,28,1.3 非线性的来源,2018/12/16,29,
7、在结构力学模拟中有三种非线性的来源: 材料非线性边界非线性几何非线性,2018/12/16,30,材料非线性,弹塑性材料轴向拉伸应力应变曲线,2018/12/16,31,材料非线性,橡胶类材料应力应变曲线,2018/12/16,32,材料非线性材料的非线性也可能与应变以外的其它因素有关。应变率相关材料的材料参数和材料失效都是材料非线性的表现形式。材料性质也可以是温度和其它预先设定的场变量的函数。,2018/12/16,33,边界非线性,将碰到障碍物的悬臂梁,2018/12/16,34,边界非线性,边界非线性是极度不连续的,在模拟分析中发生接触时,结构的响应特性会在瞬时发生很大的变化。,2018
8、/12/16,35,边界非线性,另一个边界非线性的例子是将板材材料冲压入模具的过程。在与模具接触前,板材在压力下的伸展变形是相对容易产生的,在与模具接触后,由于边界条件的改变,必须增加压力才能使板材继续成型。,2018/12/16,36,几何非线性,与分析过程中模型的几何形状的改变相联系。几何非线性发生在位移的大小影响到结构响应的情形。这可能由于:大挠度或转动“突然翻转”初应力或载荷强化,2018/12/16,37,几何非线性,悬臂梁的大挠度,2018/12/16,38,几何非线性,平板的突然翻转,2018/12/16,39,1.4 非线性问题的求解,2018/12/16,40,非线性载荷位移
9、曲线,2018/12/16,41,一般在商业软件中均使用Newton-Raphson法来求解非线性问题。,在非线性分析中的求解不能像线性问题中那样,只求解一组方程即可;而是逐步施加给定的载荷,以增量(increment)形式趋于最终解。因此非线性有限元方法将计算过程分为许多载荷增量步(load increments),并在每个载荷增量步结束时寻求近似的平衡构形。通常要经过若干次迭代才能找到某一载荷增量步的可接受的解。所有这些增量响应累加的和就是非线性分析的近似解。,2018/12/16,42,(a) 模拟计算中的外部载荷 (b) 作用于节点上的内部作用力,物体上的外部载荷和内部作用力,为了使物
10、体处于平衡状态,每个节点上施加的净作用力必须为零。因此平衡的基本判据为内部作用力I 和外部作用力P必须互相平衡: PI 0,2018/12/16,43,分析步(steps),增量步(increments)和迭代步(iterations),分析步(steps)模拟计算的加载过程包含单个或多个步骤,所以要定义分析步。它一般包含分析过程选择、载荷选择和输出要求选择。每个分析步都可以采用不同的载荷、边界条件、分析过程和输出要求。例如: 步骤一:将板材夹于刚性夹具上。 步骤二:加载使板材变形。 步骤三:确定变形板材的自然频率。,2018/12/16,44,分析步(steps),增量步(increment
11、s)和迭代步(iterations),增量步(increments)增量步是分析步的一部分。在非线性分析中,一个分析步中施加的总载荷被分解为许多小的增量,这样就可以按照非线性求解步骤来进行计算。当给定初始增量步的大小后,软件一般会自动选择后继的增量步的大小。每个增量步结束时,结构处于(近似)平衡状态,结果可以写入输出数据库文件(.odb)、重启动文件(.res)、数据文件(.dat)或结果文件(.fil)中。如果选择在某一增量步将计算结果写入输出数据库文件,我们就给这个增量步起个名字叫做帧(frames)。,2018/12/16,45,分析步(steps),增量步(increments)和迭代
12、步(iterations),迭代步(iterations)迭代步是在每一个增量步中寻找平衡解的一种尝试。如果模型在迭代结束时不是处于平衡状态,将进行新一轮迭代。随着每一次迭代,得到的解将更接近平衡状态;有时需要进行许多次迭代才能得到平衡解。当平衡解得到以后一个增量步才完成,结果只能在一个增量步的末尾才能获得。,2018/12/16,46,平衡迭代和收敛性,结构对于一个小的载荷增量P的非线性响应示于图中。非线性有限元方法利用基于构形u0时的结构初始刚度K0,和P来计算结构的位移修正值ca。利用ca将结构的构形更新为ua。,2018/12/16,47,平衡迭代和收敛性,收敛性(convergenc
13、e)基于结构新的构形ua,形成新的刚度Ka。利用Ka来计算更新后的构形中结构的内部作用力Ia。所施加的总载荷P和Ia的差值可如下计算: RaPIa 其中Ra是迭代的作用力残差值。,2018/12/16,48,平衡迭代和收敛性,收敛性(convergence)如果Ra在模型的每一自由度上均为零,图中的a点将位于载荷位移曲线上,结构将处于平衡状态。在非线性问题中,几乎不可能使Ra等于零,因此将Ra与容许残差进行比较。如果Ra比作用力容许残差小,就接受结构的更新构形作为平衡结果。默认的容许残差设置为结构中对时间进行平均的作用力的0.5%。在整个模拟过程中自动在空间分布的这个对时间的平均值。,2018
14、/12/16,49,平衡迭代和收敛性,收敛性(convergence)若Ra比目前的容许残差小,就认为P和Ia处于平衡状态,ua就是结构在当前载荷下合理的平衡构形。而在接受此解前,还要检查位移修正值ca与总的增量位移uauau0相比是否是一小量。若ca大于增量位移的1%,将重新进行迭代。只有这两个收敛性检查都得到满足,才认为此载荷增量下的解是收敛的。上述收敛判断规则有一个例外,即所谓线性(linear)增量情况。线性增量的定义是指增量步内最大的力残差小于时间平均力乘以10-8的增量步,凡严格满足这个定义的增量步被认为是线性的并无需再进行迭代,无需进行任何检查即可认为其解是可接受的。,2018/
15、12/16,50,平衡迭代和收敛性,收敛性(convergence)若迭代结果不收敛,将进行下一次迭代以使内部和外部作用力达到平衡。第二次迭代采用前面迭代结束时计算得到的刚度Ka和Ra一起来确定新的位移修正值cb,使得系统更加接近平衡状态(见图中的点b)。,2018/12/16,51,平衡迭代和收敛性,收敛性(convergence)利用结构新构形ub中的内部作用力计算新的作用力残值Rb,再次将任意自由度上的最大作用力残值与作用力残值容许值进行比较,将第二次迭代的位移修正值cb与增量位移ubub u0进行比较。如果需要的话将进行进一步的迭代。,2018/12/16,52,平衡迭代和收敛性,收敛
16、性(convergence)对于非线性分析中的每次迭代,都要重新形成模型的刚度矩阵并求解方程组。从计算费用的角度来说,这意味着每次迭代等价于进行一次完整的线性分析。现在你就可以清楚地看到非线性分析的计算费用可能要比线性问题大许多倍了。可以在每一收敛的增量步都保存结果,所以非线性模拟计算中得到的输出数据量将是线性分析中得到的数据量的很多倍。因此在规划计算机资源时,必须考虑这些因素及所想进行的非线性模拟计算的类型。,2018/12/16,53,自动增量控制,自动调整载荷增量步的大小,因此它能便捷而有效地求解非线性问题。用户只需在每个分析步计算中给出第一个增量步的大小,就会自动调整后续增量步的大小。
17、若用户未提供初始增量步大小,会试图将该分析步的全部载荷都作为第一增量步载荷来施加,这样在高度非线性的问题中不得不反复减小增量步大小,从而导致CPU时间的浪费。一般来说,提供一个合理的初始增量步大小(见例题中的“修改模型”)将是很有利的;只有在很平缓的非线性问题中才可能将分析步中的所有载荷施加于一个增量步中。,2018/12/16,54,自动增量控制,在一个载荷增量里得到收敛解所需的迭代步数会随系统的非线性程度而变化。默认情况下,如果在16次迭代中仍不收敛或出现发散,会放弃当前增量步,并将增量步大小置为先前值的25%,重新开始计算,即利用比较小的载荷增量来尝试找到收敛的解。若此增量仍不收敛,将再
18、次减小增量步大小。允许一增量步最多有五次被减小,否则就会中止分析。如果增量步的解在少于五次迭代时就收敛,这表明找到解答相对很容易。因此如果连续两个增量步只需少于五次的迭代就可以得到收敛解,自动将增量步大小提高50%。自动载荷增量方案的详细内容在信息文件(.msg)中给出。,2018/12/16,55,1.5 用商业软件进行非线性分析,2018/12/16,56,从二十世纪60年代中期以来,大量的理论研究不但拓展了有限元法的应用领域,还开发了许多通用或专用的有限元分析软件。目前应用较多的通用有限元软件如下表所列:,通用有限元软件,2018/12/16,57,几何非线性,只需对模型做些小的修改就可
19、以将几何非线性效应包含于分析中。首先要在定义分析步时考虑几何非线性效应,要给出分析步中允许包含增量步的最大数目。如果有限元需要比此数目更多的增量步来完成分析,它将中止分析并给出出错信息。分析步中默认的增量步最大数目是100,但如果题目有显著的非线性,可能会需要更多的增量步。用户给出的增量步数目是有限元可以采用的增量步数的上限,而不是它所必须使用的增量步数。,2018/12/16,58,几何非线性,在非线性分析中,一个分析步是发生于一段有限的“时间”内的,除非惯性效应或率相关效应作为重要因素进入分析,否则这里的“时间”并没有物理含义。用户是在这个理解背景下指定初始时间增量和此分析步的总时间的。这
20、些数据也指定了第一个增量步中所施加的载荷的相对于总体载荷的大小。初始载荷增量如下确定:,2018/12/16,59,几何非线性,初始时间增量的选择对于某些非线性模拟计算可能会很关键,但对大多数分析来说,初始时间增量的大小介于总分析步时间的5%到10%之间通常是足够的。除非有类似模型中包含率相关材料效应或阻尼器等情况出现,在静态模拟计算时,为了方便,总分析步时间通常均置为1.0。当总分析步时间为1.0时,所施加载荷的比例总是等于当前的时间步大小,也就是说,时间步为0.5时施加的载荷也就是总载荷的50%。,2018/12/16,60,几何非线性,尽管初始增量大小必须指定,但后面的增量大小却由有限元
21、自动控制。尽管也可以对增量大小进行进一步的人工控制,这种自动控制对于大多数非线性模拟计算来说是适合的。如果收敛性问题造成过多的增量减小,使得增量值降到了最小容许时间增量以下,有限元就会中止分析。默认的最小容许时间增量Tmin为10-5乘以总分析步时间。除了总分析步时间的限制外,有限元对最大容许时间增量Tmax没有默认的上界限制。用户也可以根据有限元模拟计算的实际情况,指定不同的最小和/或最大容许的时间增量大小。例如:如果知道模拟计算在所加载荷过大时,求解会出现问题(这也许是因为模型会经历塑性变形),所以可能想减小Tmax,此时就可以指定最大容许时间增量。,2018/12/16,61,几何非线性
22、,局部方向 在几何非线性分析中,局部的材料方向在每个单元中可以随变形而转动。对于壳、梁及桁架单元,局部的材料方向总是随变形而转动的。对于实体单元,仅当单元上定义了非默认的局部材料方向时,它的局部材料方向才随变形而转动;而在默认情况下局部材料方向在整个分析中将始终保持不变。在节点上定义的局部坐标系方向在整个分析中始终保持不变,它们不随变形而转动。,2018/12/16,62,几何非线性,对后继分析步的影响 一旦在一个分析步中包括了几何非线性,所有的后继分析步中都会自动考虑几何非线性。如果在一个后继分析步中没有要求几何非线性,有限元会发出警告信息并声明这一步将被自动强制包括几何非线性。,2018/
23、12/16,63,几何非线性,其它几何非线性效应 模型的大变形并不是唯一要考虑的重要几何非线性效应。有限元中刚度矩阵的计算也包括由于施加荷载引起的单元刚度项(称为载荷刚度)的计算,这些项能改善计算的收敛情况。另外,壳的薄膜荷载、缆索和梁的轴向荷载都对结构对横向荷载响应的刚度产生很大影响。所以在考虑几何非线性时,这种由于薄膜刚度对横向荷载的影响应该考虑在内。,2018/12/16,64,输出控制,非线性分析可能比线性分析得到远为多的结果,原因有两个:1. 许多输出项目只有在非线性分析中才有(如, 每个子步中平衡迭代次数,辐射单元的热流) 2. 多子步结果在每个载荷步中可以使用。可能需要将这些信息
24、存储到输出文件(Jobname.OUT)和/或结果文件(Jobname.RTH) 中用于分析和判断是否有足够的结果可以绘制响应曲线。,在结果文件中可以只存储. . . 自由度解. . (类型) 每五步求解. . . (频率) 单元类型3 (范围),ANSYS 输出控制命令允许用户指定 类型, 频率, 和数据的范围 。,例如,2018/12/16,65,缺省情况下,ANSYS将最后子步的所有单元所有结果写入热分析文件 (Jobname.RTH)。使用下列菜单指定类型,频率和输出范围. . .,1,2,反应热流 单元热流和梯度 单元各种数据,3,3. 选择要控制的结果项目。 例如: 所有求解项 (
25、缺省) 节点温度解,2018/12/16,66,4. 指定结果文件输出频率。如果是I “Every Nth” ,用一个整数填充 “Value of N” 。 如果要使用数组指定输出时间点,使用 “At Time Points” 选项。该选项在第6章中有详细解释。,4,5. 指定预先定义的节点或单元组元。6. 结束请按OK。如果要再设置输出控制时请按APPLY 可以设置多个输出控制。注:缺省情况下,ANSYS同时存储最近的求解结果。,5,6,注: 缺省情况下,ANSYS同时存储最近的求解结果。,2018/12/16,67,求解管理,缺省情况下,ANSYS在工作路径写一个求解管理文件。它包含诸如时
26、间/载荷步,使用的CPU时间,每次迭代子步的最大温度数值等有用的项目。 这一功能在进行时间较长的非线性和瞬态分析时最有效,因为它给分析人员一个完整的数据文件,而不需要在后处理中读取。例如:,Q: 这是非线性分析吗? (见行 4),2018/12/16,68,用户可以定义不多于三个变量来跟踪模型特定节点的温度响应或反应热流。使用下列菜单定义变量:,1,2,3,3. 选择管理的节点。 4. 指定变量号码。 5. 对于节点温度指定TEMP,对于节点反应热流指定HEAT。 6. 点击OK。,2018/12/16,69,收敛验证,许多问题可以造成非线性求解不收敛。在缺省情况下,ANSYS如果发现问题不收
27、敛,求解就会终止,并且最后的不收敛结果会写入结果文件以供分析。用户必须在后处理之前知道求解是不收敛的。ANSYS用下列几种方法指出求解是不收敛的:出错文件(jobname.err)会清楚的判断不收敛的解,并且会对不收敛的可能原因加以说明:,* ERROR * CP= 345.870 TIME= 17:36:15Preconditioned conjugate gradient solver error level 1. Possibly, themodel is unconstrained or additional iterations may be needed. Tryrunning w
28、ith a multiplier MULT 1 in EQSLV command (3 MULT 1).* ERROR * CP= 388.930 TIME= 17:36:58DOF (e.g. Displacement) limit exceeded at time 4(load step 4 substep 1 equilibrium iteration 1)Maximum value= 6.871541204E+14 Limit= 1000000.May be due to an unrestrained or unstable model.* WARNING * CP= 404.140
29、 TIME= 17:37:13The unconverged solution (identified as time 4 substep 999999) isoutput for analysis debug purposes. Results should not be used forany other purpose.,2018/12/16,70,在通用后处理器中查询results summary时,不收敛的求解结果会被指定为子步数目999999 (下面的红色箭头):前面收敛的子步(绿色箭头)可以象平常一样后处理。,Converged Converged Converged Uncon
30、verged,2018/12/16,71,非线性分析后处理考虑因素,非线性分析后处理时,可以使用与线性分析相同的选项来查看结果。但是,有几点必须注意: 误差估计 - 对于线性和非线性热分析,网格误差估计都是有效的。要注意的是ANSYS在计算误差时使用的是应变自由参考温度的随温度变化材料特性而不是实际温度。 结果对于时间的图形显示 - 非线性分析在存储子步数据时会得到比线性分析多许多的数据。为了方便理解这些大量的数据,时间历程后处理器提供了结果相对时间的图形显示功能(如,3号节点的温度随时间的变化)。非线性求解数值如时间步和迭代次数同样可以相对时间绘制。第5章讲解理使用时间历程后处理器的方法。,2018/12/16,72,动画 - 另一个得到非线性系统在载荷历程中的响应的方法是将结果数据动画显示。使用下列菜单进行随时间变化的动画显示:,2,3,6,1,指定做动画的 数据序列,指定做动画的数据项,4,5,注: “Animate Over Time” 选项将结果随固定温度间隔的数值插值并显示。,2018/12/16,73,本章完,