本系列文章研究成熟的有限元理论基础及在商用有限元软件的实现方式,通过
1) 基础理论
2) 商软操作
3) 自编程序
三者结合的方式将复杂繁琐的结构有限元理论通过简单直观的方式展现出来,同时深层次的学习有限元理论和商业软件的内部实现原理。
有限元的理论发展了几十年已经相当成熟,商用有限元软件同样也是采用这些成熟的有限元理论,只是在实际应用过程中,商用CAE软件在传统的理论基础上会做相应的修正以解决工程中遇到的不同问题,且各家软件的修正方法都不一样,每个主流商用软件手册中都会注明各个单元的理论采用了哪种理论公式,但都只是提一下用什么方法修正,很多没有具体的实现公式。商用软件对外就是一个黑盒子,除了开发人员,使用人员只能在黑盒子外猜测内部实现方式。
一方面我们查阅各个主流商用软件的理论手册并通过进行大量的资料查阅猜测内部修正方法,另一方面我们自己编程实现结构有限元求解器,通过自研求解器和商软的结果比较来验证我们的猜测,如同管中窥豹一般来研究的修正方法,从而猜测商用有限元软件的内部计算方法。我们关注CAE中的结构有限元,所以主要选择了商用结构有限元软件中文档相对较完备的Abaqus来研究内部实现方式,同时对某些问题也会涉及其它的Nastran/Ansys等商软。为了理解方便有很多问题在数学上其实并不严谨,同时由于水平有限可能有许多的理论错误,欢迎交流讨论,也期待有更多的合作机会。
自主结构有限元求解器iSolver介绍视频:
http://www.jishulink.com/college/video/c12884
==显式求解 Step By Step==
动力学问题是将力的方程和运动学方程耦合在一起的理论,在实际问题中,运动过程无时不刻的与力相关,因此数学上两个方程联立求解是最理想的。但数值计算中却难以实现,也没必要联立求解,只需先求一个方程,然后再求另一个方程,只要时间步长足够短,那么精度依然可以保证。有限元亦是如此解耦,有限元两个时刻点之间会认为力的状态不变,那么可以分两步求解。
1) 在每个时刻点只求力学方程,得到运行学的初始条件。
2) 在两个时刻点之间由于力的状态不变,那么可以只按运动学来求的运动学量。
本文首先研究商用有限元中最常用的显式求解算法中心差分法的理论,并给出了一个弹簧显式动力学分析的Step by Step例子,通过这个例子猜测了Abaqus中采用中心差分法求解显式动力学问题的过程,然后在自编有限元程序iSolver采用同样的算法,验证iSolver的结果和Abaqus完全一致,从而证明Abaqus的内部算法和我们设想的一致。
中心差分法的标准理论可查看相应的论文,由于和Abaqus的中心差分法比较接近,所以在此不累述。
式中,
表示速度,
表示加速度,
表示时间步大小。
带入i时刻的位移,则有,
其中,i时刻的加速度可根据牛顿定律计算得出,即,
式中,M表示集中质量阵,F表示外力,I表示内力。
针对初始化、处理某些约束条件以及在后处理过程中,Abaqus对速度有特殊的处理,如下式,
在初始时,即(t=0),除非用户自定义,一般情况下速度和加速度都设为0。假定有下式,
带入式(1),则有
将i=-1带入式(4)中,也可得到式(6),故前后是一致的。
初始化流程如下,
以某一时刻t为例,整体流程如下:
Example 1:弹簧的显式动力学分析
(模型详见附件Job-ExplicitSpring-AddF-NoBC.inp)
创建一根线弹簧,在弹簧两侧各加两个点质量,无约束,右端X方向拉力。
参数如下:
尺寸:X方向长度L=1;
质量:各向同性,大小为100;
刚度:弹簧刚度为1000;
力:分析时间内恒定大小,为800;
时间设置:总时间为1,时间增量固定为0.2。
与第一个增量步类似,不再累述。
可以发现iSolver和Abaqus完全一致。
本文概要性地介绍了Abaqus中心差分法的理论以及算法实现的整体流程,并通过简单弹簧显式动力学分析算例与Abaqus计算结果进行对比,验证了算法和整体流程的正确性。后续文章,我们会逐步深入显式动力学的一些细节,敬请期待。
如果有任何其它疑问或者项目合作意向,也欢迎联系我们:
snowwave02 From www.jishulink.com
email: snowwave02@qq.com
免责声明:本文系网络转载或改编,未找到原创作者,版权归原作者所有。如涉及版权,请联系删