本系列文章研究成熟的有限元理论基础及在商用有限元软件的实现方式,通过
1) 基础理论
2) 商软操作
3) 自编程序
三者结合的方式将复杂繁琐的结构有限元理论通过简单直观的方式展现出来,同时深层次的学习有限元理论和商业软件的内部实现原理。
有限元的理论发展了几十年已经相当成熟,商用有限元软件同样也是采用这些成熟的有限元理论,只是在实际应用过程中,商用CAE软件在传统的理论基础上会做相应的修正以解决工程中遇到的不同问题,且各家软件的修正方法都不一样,每个主流商用软件手册中都会注明各个单元的理论采用了哪种理论公式,但都只是提一下用什么方法修正,很多没有具体的实现公式。商用软件对外就是一个黑盒子,除了开发人员,使用人员只能在黑盒子外猜测内部实现方式。
一方面我们查阅各个主流商用软件的理论手册并通过进行大量的资料查阅猜测内部修正方法,另一方面我们自己编程实现结构有限元求解器,通过自研求解器和商软的结果比较来验证我们的猜测,如同管中窥豹一般来研究的修正方法,从而猜测商用有限元软件的内部计算方法。我们关注CAE中的结构有限元,所以主要选择了商用结构有限元软件中文档相对较完备的Abaqus来研究内部实现方式,同时对某些问题也会涉及其它的Nastran/Ansys等商软。为了理解方便有很多问题在数学上其实并不严谨,同时由于水平有限可能有许多的理论错误,欢迎交流讨论,也期待有更多的合作机会。
从NASA在1964年推出Nastran算起,结构CAE软件已经发展了将近60年,虽然在60年内商用结构软件层出不穷,但最流行的通用结构CAE软件现在依然只有Nastran、Ansys和Abaqus三款软件,iSolver以前的结果精度完全对标Abaqus,线性和材料非线性与Abaqus计算结果对比精度已经没有误差,但在推广过程中,发现在线性问题上,客户只相信Nastran的结果,因为很多工程的校核规范和经验方法都以Nastran的结果为基准,如果换成其它软件,那么那些经验系数和校核公式都可能会修改,而这些修改还需要做大量的工程验证,这是客户承受不了的代价。
客户改不了的习惯只能是自主软件改,我们的自主软件只是一个后来者,想要推广只能按客户要求修改精度,因此,iSolver开始兼容Abaqus和Nastran的精度,计算出的结果既可以对标Abaqus,修改设置后又可以对标Nastran,因此需要对Abaqus和Nastran的各方面的算法差异进行深入研究。如果是线性问题,那么Nastran和Abaqus的精度误差主要体现在单元算法、边界处理、MPC约束关系等,在2017年第二篇:S4壳单元质量矩阵研究文章中我们就曾经分析过Abaqus的S4壳单元和Nastran的Quad4壳单元质量矩阵的内部实现方式和差异,在这里主要研究Abaqus、iSolver与Nastran梁单元差异,由于这三款软件的梁单元的差异较多,我们分几篇文章来说明,
有限元是求受力情况下的位移等变形情况,也就是位移等未知量和外力存在一定关系。对于一根三维实体梁,梁实际受到的外力是三维全局空间的,如果直接用全局坐标系下三维的力来求梁的受力分析,那么就需要对梁划分为三维的体单元求解,网格数目和计算效率比较差,一种简单方法是对那些细长的梁(Abaqus认为是细长比>8),此时可以用简单的等效为线单元的形式来表达位移和外力的关系,这样只要用一个线单元就可以表示这个三维实体梁了,大大简化了求解矩阵。
实际的加载是多个力的组合,譬如下方采用手轮加载的力、弯矩和扭矩外载荷
但梁的有限元中可以把这个线单元受力关系分为:
(1) 轴向拉伸力
(2) 轴向扭转
(3) 横向弯曲力,可以加力载荷或者弯矩
三部分,此时每部分都有简单的位移和外力的公式,也就是存在一个局部坐标系,简化梁理论总是先求出梁单元局部坐标系的刚度和质量阵,然后再用三维变换直接转到全局坐标系下。
对(1)(2)轴向的受力,沿梁的轴向方向,而对(3)弯曲力,沿截面方向。
梁的轴向方向是定死的,就是两个节点连线就行,但截面的方向是不定的,对一个平面来说,只需要两个互相垂直的方向1和2。梁的弯曲效应主要由I11(1方向惯性矩)、I12(1-2方向惯性矩)、I22(2方向惯性矩)、形心位置等参数决定,对于圆形或者圆管,这些参数和1、2取哪个方向无关,但其它型材就相关了。对同一个截面,理论上1、2方向可以任意取,实际上基本原则总是取1、2方向后这些截面参数可以很容易的计算,譬如对L、T等肯定是沿翼板和沿腹板来取方向。但到底是1方向沿翼板还是2方向沿腹板,或者沿腹板的方向每个软件并不相同。
一般软件都有几种方式设置梁的截面方向,有几种:
(1) 指定一个三维方向矢量v
(2) 指定一个三维节点
(3) 指定一个转角,这个转角为默认的梁的方向沿轴线旋转
这三种方式后台最终都是先求出一个截面方向,另外一个方向只要满足右手定则即可,在Abaqus、Nastran、iSolver中我们都以第一种方式直接指定一个三维矢量v来说明截面方向。同时以下方的L型材来说明同一个型材在不同软件中的截面方向(算例名称为OneBeam.cae/Beam-OneMesh-L)。
在Abaqus中,创建一根梁由节点1和2组成,同时,节点1、2分别是梁的第一和第二节点。此例子中x=0为节点1。
梁局部坐标系的x方向永远都是1->2节点,abaqus中称为t方向,如下方:
设置梁方向时,输入n1的一个三维方向矢量,简单起见,n1和t直接取垂直,取默认的0,0,-1:
Abaqus后台得到
局部坐标系的z方向(即截面2的方向)Abaqus.2= t×n1
局部坐标系的y方向(即截面1的方向)Abaqus.1=Abaqus.2×t
最终t,S.1,S.2满足右手定则,得到局部坐标系方向
三维显示为:
很简单,梁截面几何尺寸的设置方向的1、2就是Abaqus的局部坐标系的y、z轴。
Nastran的局部坐标系的x方向和Abaqus完全一致,都是节点1->2方向t。和Abaqus一样的模型如下:
Patran中可以设置Bar Orientation,在Nastran中称为v方向,同样设置为0,0,-1:
按照Nastran帮助手册,Patran设置的局部坐标系和Abaqus完全一致,Nastran由t和v首先确定一个局部坐标系Patran的y,z方向,此时
Patran.z= t×v
Patran.y=Patran.z×t
即下图表示:
那既然局部坐标系和Abaqus完全一致,那么Abaqus的L型定义的参数输入到Patran中是否在三维全局坐标系下也完全一致呢?
可惜不是的,把上面的L型几何参数四个值原封不动输入到Patran的Section中:
Patran打开三维显示梁的方式,转到Abaqus的同一个角度,显然实体和Abaqus完全不同,Nastran的后台计算的刚度矩阵等必然也和Abaqus不同了。
所以型材几何尺寸的设置方向和Abaqus不同
Nastran后台计算时局部坐标系的Iyy和Izz分别采用梁截面几何尺寸设置的I22和I11。
梁截面几何尺寸的方向的向上(即1方向)是Abaqus局部坐标系的y,截面方向的向右(即2方向)是Abaqus局部坐标系的z方向。
想要Nastran结果和Abaqus一致,只需要把yz颠倒就行,譬如按这个原则输入上面L型材的Patran的Section的四个参数,把1、2方向颠倒:
在Patran全局坐标系下显示三维模型,可发现和Abaqus完全一致:
iSolver的梁截面方向采用Abaqus的形式,不过后台也支持了Nastran的梁截面按Nastran形式的自动转换,使得iSolver能同时处理Abaqus和Nastran的梁模型定义问题。
譬如对上面的截面形式,L型梁的创建参数如下:
在iSolver中三维显示如下,和Abaqus完全一致:
如果觉得上面的文字太复杂,也可以看一下视频的简要讲解和软件演示,地址如下:
https://www.jishulink.com/college/video/c12884 20理论系列文章38-梁单元差异(2)-梁截面方向
在Nastran和Abaqus(iSolver)的梁截面的方向上,区别表示如下,第二个截面几何尺寸的设置方向其实和局部坐标系完全无关,但个人觉得这两个方向很容易搞混,如果你想在Abaqus(iSolver)和Nastran针对同一对象建模,为了结果的一致性,建议直接忽略属性关联时的局部坐标系,直接认为几何尺寸的设置方向就是局部坐标系来建模就行了。
免责声明:本文系网络转载或改编,未找到原创作者,版权归原作者所有。如涉及版权,请联系删