系统主要思想和做法
计算机发展至今,在带动社会不断发展和变化的同时,其本身也在不断的发展与进步,从机器代码到高级语言(FORTRAN/C)再到面向对象的程序设计方法的发展,已经引领社会发生了巨大的变革,目前人们都在期待着更为方便的新一代语言的出现。针对有限元研究领域,FEPG为各专业的研究人员创造了一门独特的有限元语言,为各学科领域的数值计算提供了极大的方便,这种语言可以节省90%以上的软件编写代码量。例如针对一个三维稳态热传导问题,采用FEPG提供的有限元语言,用户只需填写以下VDE/PDE、GCN、GIO等几个简单的文件即可生成全部生成有限元计算源代码。
稳态热传导控制方程:
假设求解边界条件为:
按照有限元方法,转化为虚功方程的弱形式:
对应的偏微分方程PDE文件书写如下:
计算流程由以下几行语言规定:
以上为采用FEPG独特的有限元语言来求解一个三维稳态热传导问题所需要的文件,它的书写形式接近于通常习惯的写法,采用有限元语言描述有限元方程的代码量远远低于直接采用FRORTRAN/C编程,因此可以把我们从庞大复杂的编程劳动中解脱出来,将更多精力投入到寻找合理的物理模型与更准确的算法中。
为什么会有如此效果?原因就在于FEPG采用了元件化的程序设计方法和程序自动生成技术,有限元程序的元件化,使用户在保证程序完整性的同时可以根据需要更改核心的算法和公理论式,大大减少了工作量,又保证了程序的灵活性和适用性。
1、元件化思想
有限元程序总体可分为三个组成部分:前处理;有限元计算;后处理。
FEPG系统采用元件化程序设计方法,把庞大的有限元计算程序分解为若干个元件程序。每个元件程序都是一个完整的FORTRAN程序,可以单独进行编译、连接与运行,它们之间的通信完全通过磁盘文件。由元件程序组成的程序系统是通过命令流产生的批命令方式实现的,即把所要执行的元件程序按照一定的顺序以操作系统所能接受的命令流方式写在一个命令流文件中,由FEPG系统自动生成相应的批命令文件,然后运行这个批命令文件,由操作系统解释执行。
有限元计算程序主要由START、BFT、SOLV和E、U五个元件程序组成,这五个元件程序的功能如下:
START程序:该程序给出每个节点的各个自由度与将来要形成的代数方程组的变量(即方程号)的对应关系(即哪个节点的哪个自由度将要对应于方程组的哪个变量)以及解的初值。
BFT元件程序:此元件程序的主要功能是给出每一时刻解的边值,即指定节点位移和载荷,以及对时间的更新和生成保存计算结果的批命令文件“post.bat”。
E元件程序:E元件程序用于计算单元刚度、质量和载荷等,并把它们由节点各自由度表示转换成由代数方程组的变量表示,同时处理边界约束条件,并形成代数方程组的右端项。
SOLV求解器:求解器用于叠加形成总体刚度矩阵及求解线性代数方程组。
U元件程序:U元件程序用于把求解器求出的变量位移转换为节点各自由度的位移以及其它的后处理计算。
2、有限元程序自动生成原理
采用元件化程序设计方法把完整的有限元程序分解成若干个元件程序后,进一步将这些元件程序分解为可变部分和不可变部分。不变部分系统直接给出,可变部分是根据方程和算法用有限元语言描述,通过生成系统自动产生,然后可变部分和不变部分组成完整的有限元程序。
下面的流程图说明,通过有机组合元件程序,即可实现静态或动态、线性或非线性等各种有限元问题的求解。
3、应用FEPG进行有限元计算的一般过程
在FEPG系统中,用户可以以两种方式得到计算所需的全部程序,如下图所示,
第一种方法是使用公式库生成程序,我们把常见的物理问题,如固体力学、电磁场、传热、渗流、流体等的描述方程用有限元语言描述好放在公式库中。用户只需点击公式库菜单即可生成用户所需的全部有限元计算程序。
另一种方法便是由用户公式生成程序。用户根据自己研究的物理问题,用有限元语言将控制方程写成VDE、PDE文件,将计算方法写成gcn和gio文件,然后用FEPG系统命令(Gio命令)产生全部有限元程序。
4、多场耦合问题
FEPG最善于解决多场耦合问题。多物理场耦合问题有限元求解的主要困难是多物理场往往不能直接联立求解,因为不同物理场的微分方程的性质差异比较大,因此需要分开迭代求解;由此引起的耦合问题算法因问题而异,数据结构复杂,编程工作量大;现实世界中耦合问题是千变万化的,不同的场耦合,不同的耦合方式组合成无数种模式。在此情况下按通常的方式编制程序很难跟上对耦合场求解的需要。
由于FEPG采用了元件化程序设计和代码生成技术,在求解多物理场耦合问题时,有自己独特的处理方法。因为有多个控制方程,用户需准备多个VDE/PDE文件,然后在GCN文件中描述各场之间的耦合关系,最后以和单场同样的方式生成多场耦合有限元程序。采用FEPG生成耦合问题的有限元程序,既可以让您免去大量繁琐的编程劳动,保证了程序的正确性和统一性,又可以在新的算法出现时即时更新,重新生成计算程序,满足您在计算精度、稳定性、计算时间等各方面的要求。
FEPG系统的特点
通过以上对系统结构和原理的介绍,我们就不难理解FEPG系统与众不同的三大特点:
1、采用元件化设计方法
因此大大降低程序的复杂性,大大提高程序的可读性和再用性,减少代码量90%以上。为有限元程序的维护和发展创造了前所未有的前景。
2、自动生成源程序
FEPG完全基于有限元方法的基本原理(虚位移原理),不受专业领域的限制,各种有限元问题和有限元方法均可由用户填写微分方程描述文件和算法文件自动生成全部计算程序。
3、开放源代码
用户可以看到由公式库或自己的方程生成的所有源程序,并可进行修改和编译。同时系统允许用户插入自己的FORTRAN源程序,很容易与其它FORTRAN程序融合。这就意味着用户可以很容易采用本系统开发自己的软件,也可以将程序生成的代码嵌入到自己的软件中。这是目前所有商业软件都不能实现的。
正由于具有以上这些特点,使得使用FEPG可以方便地生成求解各种多学科、多物理场的非线性耦合问题的有限元程序,并可获得有限元问题的全部FORTRAN源程序。使得以往需要数月甚至数年来完成的编程劳动在数天甚至数小时内即可完成。
应用案例
前面第一节所描述的是一个简单的传热学算例,针对具体的复杂工程和科学问题,FEPG也能胜任,解决他们所不能解决的问题,下面介绍三个不同领域的应用案例,从中我们可以了解到我们应用FEPG解决复杂工程和科学问题的能力。
一、FEPG在大变形金属平板透缝结构中的应用
1 问题说明
一金属结构,其上设有6条贯穿厚度的透缝,初始为平板形状,受压面一侧有不透气薄膜贴合。周边被固定约束。在逐步增加的侧向均布面载荷(准静态)作用下,结构不断产生挠曲变形,直至结构达到拉伸失稳条件而失去继续承载的能力。欲得到极限承载压力以及变形过程中,位移、应变和应力的变化及分布规律。材料为多线性的(等效)应力应变关系
2、方程描述
平衡方程
本构方程
上面的方程可简写为
3、文件的填写
1)、计算位移的文件
将上面的公式写成PDE文件并且加上计算等效应变、读取应力应变曲线以及修改E的程序的nesa.pde文件如下:
disp u v w
coor x y z
自定义变量,由上一迭代步得到的位移计算等效应变,由应力应变关系得到弹性模量E
$cv enxx={un/x}
$cv ee=dsqrt(ee)*dsqrt(2.0d0)/3
$c6 do i=1,nee
$cv if(ee .ge. ees(i) .and. ee .lt. ees(i+1))then
$cv sgm=(ee-ees(i))*(sigmas(i+1)-sigmas(i))/(ees(i+1)-ees(i))
$cc +sigmas(i)
$cv if(ee .le. ees(2))then
$cv pe=sigmas(2)/ees(2)
$c6 else
$cv pe=sgm/ee
$cv endif
$c6 goto 1111
$c6 else if(ee .ge. ees(nee+1))then
$cv write(*,*)’超出所给数据表—-out data file’
$cv pe=sigmas(nee+1)/ee
$cv endif
$c6 enddo
$c0 1111 continue
自定义变量,并书写刚度矩阵和载荷项
$c6 fact = pe/(1.+pv)/(1.-2.*pv)*vol
$c6 shear = (0.5-pv)
stif
dist=+[exx;exx]*(1.-pv)*fact+[exx;eyy]*pv*fact+[exx;ezz]*pv*fact
+[eyy;exx]*pv*fact+[eyy;eyy]*(1.-pv)*fact+[eyy;ezz]*pv*fact
+[ezz;exx]*pv*fact+[ezz;eyy]*pv*fact+[ezz;ezz]*(1.-pv)*fact
+[eyz;eyz]*shear*fact+[exz;exz]*shear*fact+[exy;exy]*shear*fact
load=+[u]*fu*vol+[v]*fv*vol+[w]*fw*vol+[exx]*envx*fact
+[eyy]*envy*fact+[ezz]*envz*fact+[eyz]*enn(2,3)*shear*fact
+[exz]*enn(1,3)*shear*fact+[exy]*enn(1,2)*shear*fact
end
2)、计算应力的文件及算法文件
对于计算应力的单元程序,是弹塑性本构出发,采用最小二乘法。由Green-Lagrange应变得到Kirchhoff应力,再变换到Euler应力。我们也可以写出相应的pde文件。
准备好上述文件以及算法文件gcn、gio文件后,只需运行GIO命令生成程序,即可生成计算所需的全部程序
4、计算结果