我与有限元


触摸有限元时刻不算常也不短,从03年到如今吧,中心断断续续的学习有限元,大大小小软件也用了不少,牛人的帖子大师的作品读了一些,近来常来论坛逛逛,发现还有许多学弟学妹初学fem,提出的疑问形形色色但大多不对正路,不竟回想起了当年自个学习有限元时情形,孤立无助处处瞎碰,通宵熬夜啃手册也得不出一个所以然来,如今回头想想假如有几个严师诤友会是多么美好的工作,写此文仅仅只为了让后生少走些当年自个走过的弯路,有更多的精力去研讨更深层次的疑问,而不是在门槛上打转,一同也算是报答论坛多年来陪同。
这篇文章如有遗漏的当地,还请各位高手弥补纠正。

 

Chapter 1
啥是有限元 by caoer
1.1 PDE
有限元是一种求解疑问的数值办法,求解啥疑问呢?--求解PDE(偏微分方程).那么PDE是做啥用的呢?--描绘客观物理国际。我想假如这两个疑问搞清楚了也就了解了为啥要用fem,fem能够做那些东西。 PDE能够描绘许多物理现象,电磁,流体,换热,diffusion,力学,河槽变迁,物种竞争,股市金融,等等等等。。。。乃至全部宇宙,当然也不是一切的物理现象都能够用PDE描绘,如微观国际分子原子的运动等等,所以我历来都不建议用有限元办法仿真微观物质现象的因素,但也有PDE运用于微观方位如possion 方程来解析plasma的物理现象,这在量子物理里边用核算的办法过于巨大,泼松方程反而使疑问简略并且能吻合试验,这些都是题外话就不多说了。除了PDE以外,ODE一样也能够描绘客观国际,但ODE多用于操控体系,很有一些线性PDE的解法也都是将PDE转化为ODE来做解析解的
1.2 求解PDE
有了PDE今后,疑问是怎样求解并得到成果,首先要阐明的是不是一切的PDE都有解的,通常有解的PDE才有实习工程含义。关于数值解法,常用的是有限差分,有限元和谱办法,还有蒙特卡罗法。有限差分呈现的较早,核算精度相对较高,可是费时,且模型形状有必要规矩,鸿沟条件处理艰难,长处是能够对比便利的操控核算精度,适用于流体类的仿真。有限元办法功率高且满意精度请求,鸿沟条件简略处理,得到了广阔的运用,尤其是固体范畴。谱办法因为能够选用FFT办法的来求解,使得程序有着精度高,收敛快的特色,也克服了有限元条件下运用高阶插值方程核算费时的缺陷,常常运用periodic boundary condition,但也有不断添加的算法使得一类二类鸿沟变成可能,适合微观规范的PDE解,谱办法和有限元联络发生的谱元法取两者之长处,使得运用远景十分好。蒙特卡罗法不是根据弱解办法的,随机数的多维采样终究得到核算上的成果,多用于金融剖析。咱这儿仍是着重有限元解PDE,望文生义,有限元将全部核算几许模型区分为许多小的单元(element),每个单元的富含必定数量的节点(node),详细单个单元有多少节点,有对应的不一样算法与差值方程,拿一个简略的线性4节点平面单元来说,每个单元包含4个节点,每个节点有对应的variable值,比方简略固体力学疑问,每个节点就有对应的位移值,热力学疑问每个节点就有对应的温度值,等等。然后单元内部的variables就经过差值办法核算得出弱解(weak formation)是树立在变分法根底上的,经过这个办法将strong form的PDE变换为weak form,使得有限元的求解变成可能,详细怎样推导weak form能够参阅一些有限元书本,假如一本根底的有限元书本没有介绍怎样推导weak form,那么能够思考挑选换一本了。推导所得的弱解办法仍需求经过核算机来核算,Galerkin办法推导出富含求和符号的公式,在核算机中八成以loop的办法来核算这个量,通常这个量即是stiffness matrix中的component,这个loop八成是环绕gauss积分点进行的。公式中还会存在积分核算,有限元办法多用gauss quaradure的办法来核算,精度通常能够满意。也就说通常简略的限元核算中存在两次approximation,一次是Galerkin一次是gauss,这也是许多人在核算完今后需求进行validation的因素,一同关于非线性疑问newton法自身即是经过核算误差来判别收敛的,固体有限元常常运用能量增量作为newton求解器的判据。单个单元的stiffness matrix核算完成后,还需求将一切单元的矩阵安装为一个大型的矩阵,然后进行线性代数核算。这个安装是很有窍门的,因为通常状况下stiffness matrix是一个很大的稀少矩阵,0值通常能够省去以节约核算量。一个由10x10x10个8节点矩形单元组成的模型会有11x11x11个node, 假如是热力学疑问变量是温度T, 每个节点上有一个自由度,拼装后的的stiffness matrix会有(11x11x11)x(11x11x11)=1771561之大,跟着单元数或变量的添加,核算是惊人的,这么一个大的矩阵求解显然不能用惯例的办法(gaussian elimination),这即是各种这么求解稀少矩阵算法存在的因素了,有用且迅速求解线代方程组是一个好的solver的规范。
1.4 后处理
本来关于最基本的有限元办法,求解得到的仅仅是variable的值,如力学即是节点位移值,thermal problem即是温度值,流体即是位移速度加压力值,假如咱们想知道应力或许应变怎样办呢?后处理体系里边个都会添加相应子程序来核算stress, strain, flux等等。这也即是为何咱们能看到各式各样的contour的因素了,当然读者也能够自个参加核算各种量的子程序,如应变能密度啥的。关于啥是有限元就介绍到这儿,仅仅是一些漫笔和主意,详细的理论和推导需求自身实习与探究,这篇文章行文匆促仅仅论述自个对有限元的粗浅了解。有不对的当地还请纠正。
下一章会谈及一些我从前用过得软件。

 

Chapter 2 有限元软件的介绍与对比
有限元软件许多,有商业的,开源的,免费的,并行的,多物理场的,专业于某范畴
【这儿仅仅介绍一些笔者从前用过的,遗漏的还请高手弥补。】
2.1ANSYS
第一个学习的有限元软件,也是上世纪末本世纪初国内最盛行的,为啥盛行呢?应为ansys的教程漫山遍野,所以我们都学习ansys因为有教程上手快,那个时abaqus的教程可谓是百里挑一,天然学习的人也就不多,其时组里导师都是用ansys,天然我也就用了。假如一个有限元软件推行的好,那么他的教程必定要推行的好,fem是一个专业性很强的软件,教程推行不够,出售天然不可。题外话,接下来说说对用他的感触。
其时学习ansys是直奔apdl去的,加上理论布景很弱所以学习进程反常艰苦。总得说来ansys是一个很中规中具的软件,功用很强壮,核算也很可靠,速度快,基本上固体力学的疑问都能处理。作为一个工程软件仍是很不错的。帮助文件很丰厚,我也挑不出他啥缺陷。仅仅感觉国内学习ansys现已陷入了一个怪圈,用户喜爱就某一个运用窍门钻的很深,关于其根底的理论却疏忽不视,这无不于许多ansys教材的编写误导有关。是指算例,却不关怀原理,悲痛。
2.2COMSOL(femlab)
首次用comsol的时分它叫femlab2.0, 无法的是仅仅作为一种爱好学习的,因为时刻精力疑问最终也没有很深入,偶尔时机1年后又开端认真的学习,逐步体会到他的强壮,我想说的是comsol关于各种物理场的research的确是一把利器,他能够恣意的由用户输入PDE,核算成果并出contour,这关于懂得有限元理论却又不想花许多的时刻来编程的人来说,即是种bliss,可是缺陷即是速度太慢,这也是能够了解的,究竟comsol是将单元与pde信息进行重组后核算,而不是许多商业有限元现已将Pde固化在单元信息里了,可是其超卓的weak form PDE模块是稀少难得的有限元剖析模块。还有其建模与后处理模块适当便利,比ansys便利许多。至于核算的准确性笔者还抱有置疑态度,细节就不提了。国内如今这个软件推的极好,大有赶超ansys之势。值得一提的是comsol的并行核算才干,选用的umfpack+mpich的构架,全体并行功率仍是十分令人满意的。
2.3FEAP
这个本来是老大哥,尽管许多人不知道它。开源,不免费,可是很便宜,也有免费的学生版别,我们能够下载学习用。据说abaqus和ansys都是源于这套程序,abaqus技术部里边的人都是feap搞这个组里出来的,尽管是用fortran编的,可是编程思路清晰,注解丰厚,参阅文档附近,实在是稀少难得学习fem的神器。不供给gui的输入办法,指令的输入和apdl和umat办法很像,都是文本macro办法。核算速度超快,带有后处理功用,供给微弱的二次开发接口,能够自个编写子程序输出别的后处理软件如tecplot.具有openmpi并行才干,供给丰厚的elastic, plastic,hyperelastic 模块。缺陷是之用于固体力学与热核算,当然能够自个开发别的单元,如电磁的。还有即是feap只运转于linux下,运用者需求知道怎样make程序。
2.4libmesh 和 deal II
之所以把这两个程序放在一同,是因为都出自于同一个大师。c++开源免费程序,十分的robust,可是还有一些小bug存在,能够运用于更广泛的范畴如流体,微弱在于adaptive meshing,和remesh,然后能够处理许多惯例有限元无法处理的疑问,如singularity疑问。libmesh尽管是后起之秀,可是开展反常迅猛,无法只要那20多个比方作为教程,关于初学FEM的兄弟那即是灾难了,并且好像固体不是libmesh的主攻方向,假如要在libmesh上开发一个壳单元,就有必要对全部内部插值环境有极好的了解,只能说难度很大。需求有c++布景以及有限元理论知识,也是笔者遇到门槛最高的一个有限元软件,引荐给专业人事。一样是在Linux下运转的软件。
2.5 ALGOR_pipeline
这个软件在也挺有名,可是国内用的少,首要是固流耦合这快做得对比好,就我所知一些有名的飞机和流体机械厂都用他,我首要是用pipeline这个插件,或许姓名记住不对,即是用来核算全部管路体系的轰动的,前处理中设定管路,就能够核算出各阶振型与频率,很便利,据说准确性挺好。其主模块用的不多,不做评论了。
基本上笔者首要运用的即是这几种有限元软件,我想假如能通晓上面任何一种软件,学习别的的软件都是一个很快的进程。假如你是博士生且有一个对比长得研讨学习时刻,引荐用feap和libmesh, 或许是根据如PETSc/blas求解器这么的自编程序, 假如是硕士,引荐comsol,本科毕业规划,引荐ansys.当然每自己水平周围环境都不一样,不能混为一谈,总之多学一点东西没有坏处,每个软件都有其深沉的布景。

 

 

Chapter 3 有限元的数值办法
这儿先略过有限元的几许建模与网格区分有些,因为一来不是数值办法的首要有些,二来我对这块也不是很通晓,所以直接从数值解法下手。这有些能够说是ch4与ch5的全体概览有些。也是fem的中心。
3.1 单元离散化与jacobian
区分好网格后,就意味着单元编号以及节点都已断定下来了,本来这步在有限元剖析的一开端就设定好了。拿平面单元为例,假如选用平面8节点矩形单元来核算热力学疑问,那么每个单元有8个节点温度值,要知道每个节点的方位并不是规矩均匀分布在实习有限单元上的,我们通常将节点的global坐标变换到local坐标上,以便利核算推导,今后再经过jacobian变换到global上,这点和连续介质力学的reference变换是同一道理。可是不是jacobian的树立是不是有必要的呢?不是,假如一切的插值方程都是根据global的,那么local空间能够省去,jacobian天然不需求了,这种做法不被引荐,应当编出来的程序可扩展性十分差,并且核算量较大,所以在大型通用有限元里边都是核算jacobian的。
3.2 运动方程与各种矩阵
得到了单个单元的8节点test function 或 grediant of test function的积分值今后(积分的核算用gauss法),就需求联络一切的单元,安装变成一个全体的矩阵,这个矩阵即是stiffness matrix,假如是一个4单元简略正方形区域,那么安装好的stiffness matrix即是一个21x21矩阵,以k标记,除此以外,假如是瞬态疑问,也会有质量矩阵存在,八成是一个对角矩阵,其值通常是test function的积分值,假如存在damping那么还会存在damping matrix,一个有限元体系通常只存在这3个矩阵,然后经过一个运动方程将其连立起来,即 M*u_tt+D*u_t+K*u=F,F是source。这即是将pde转化为弱解再转化为有限元方程的最终办法。这个矩阵肯定是一个稀少矩阵,关于稀少矩阵传统的gaussian elimination办法好像就显得十分笨拙了,所以各式各样求解稀少矩阵的求解器就呈现了。
3.3 Ax=b求解.
有了 M*u_tt+D*u_t+K*u=F 方程,接下来要做得即是求解他,假如是一个2x2的矩阵,手艺即可核算得到,实习运用中,通常都是上万阶的矩阵。或许有人会问了u_tt和u_t是怎样解得的?瞬态疑问因为有了时刻变量的参加,一切有必要要有对应的时刻积分求解器,二阶的有newmarks, HHT,energy conserved 办法,一阶的有crank nicolson, 向前向后euler法等等,这点在求解器有些会胪陈,经过时刻求解器的核算后,运动方程还能够化简为 K_tilde*u=F的办法,也即是线代求解器需求解的最终方程,求解后,各节点的u值得到。结束。
3.4 多场耦合办法
假如有多个场存在怎样办呢?道理很简略,比方一个固体场一个温度场,两个pde出来的运动方程是这么的 M1*u1_tt+D1*u1_t+K1*u1=F1和 M2*u2_t+K2*u2=F2,将两个场的运动方程线性进行叠加(线性耦合),得到 M1*u3_tt+(D1+M2)*u3_t+(K1+K2)*u3=F1+F2 ,再化简,得到 K_tilde*u3=F3,运用求解器求解既得成果。这儿仅仅简略的线性耦合,关于杂乱物理现象的非线性耦合都是在这个根底之上进行改变。并且要留意的是,最终所得的K_tilde纷歧定是对称矩阵,也就意味着求解器有必要要有解Unsymmetric matrix的才干。
3.5 鸿沟条件加载
鸿沟条件的加载,都是在一切单元矩阵安装结束今后进行的,通常是variable或是gradiant of variable的值被限制住了,关于低阶鸿沟条件,如Dirichlet boundary condition, 常用的办法是添加一个赏罚量,这个赏罚量适当大>1e10,这么便使得u被固定在人为设定的u0方位。关于natrual boundary condition, 也即是q_n的值怎样断定,通常于变量的梯度值有关,假如q_n的表达式已知,那么只需求核算出q_n的准确值带入source相即可,需求留意的是natural boundary condition请求对比精密的网格,否则会致使部分节点的核算失真。关于固体疑问,q_n代表的是压力,详细鸿沟上的力(或弯矩在板壳杆单元中)能够经过修正源相{F}来施加载荷。

 

Chapter 4 常用固体单元
记住自个刚学ansys的时分就被自带的100多个单元模型搞的晕头转向,到底应当挑选啥样的单元?这些单元有啥不一样?这些疑问一向困惑着我,本来把每一种单元都搞透是一件不太可能的工作,但有件工作能够做得即是关于几个大类的单元有个底层数值办法上的了解,这关于单元的挑选与开发是至关重要的。这儿仅以3d单元作为比方。
有限元办法中以单个单元为根底,一切的单元都有一样的特点,这些单元组成了全部domain,能够说了解了单元也就了解了全部模型,现代有限元软件都将pde整合在了单元信息里边,还句话说,Pde所推导出的弱解办法中非鸿沟条件term都体如今了单元中。
4.1 固体单元, 变量:位移u1,u2,u3
工程力学中最常用的单元,常用于连续介质的力学核算,通常分为small deformation 和 finite deformation 两个区域,记住ansys里边有个大变形开关,一打开后核算变得不忍目睹,所以挑选small仍是finite要看详细的疑问, 此外关于不一样的资料,如超弹,粘弹,弹塑资料都有不一样的数学模型,这儿挑选的时分必定要稳重,不要想当然,最佳要看看手册和理论,核算办法上常用的有displacement,mixed,enhanced strain办法,displacement是最初始的,收敛效果不好,能够用于处理通常根底疑问,引荐后边两种来解非惯例资料。此外,关于一样的pde模型,单个单元节点越多核算越准确,单元数目一样状况下,6面体单元比4面体单元精度高。商业软件通常都思考的很周全,假如是自编软件的话,就需求思考许多疑问了,如locking, buckling, convergence 等等,这些东西很风趣也很艰深,既需求深沉的理论功底又要有实习的编程才干,假如要做实习可用的大规模疑问,对核算机cpu的构架,cache的大小,data layout 和编译器都需求极好的了解,才干编写出高质量的代码。所以编写一个small deformation的固体单元关于科技开展的今日,现已是本科生的课堂规划试验课了。但关于大变形,多场耦合,粘弹性资料,多规范,仍是需求对比深沉的张量推导内力的,所以搞力学理论布景的要多了解了解核算机有关原理及树立高效编程习气,搞核算机出世的要树立张量运算推导和物理现象的概念,那样即是有限元开发的有用人才,这句话我和广阔有限元爱好者共勉。
4.2 热单元,变量:温度T.
是相对简略的一种单元,pde自身也很简略,即是换热方程,节点自由度即是T,因为温度T没有方向性,所以T只要1维量,核算量小。别的heat flux = gradient of T, 要点参看一下Fourier heat conduction equation.
4.3 壳单元,变量:u1 , u2 , u3 , θ1 , θ2 , (θ3)
有些具有曲面外壳的资料,其一个方向的变形要远远小于别的方向的。 也有一种大变形的壳单元,只富含5个变量,少一个θ或是u,可是有可能在核算的时分不收敛,所以在运用不一样的shell单元时,必定要多看看阐明,遇到不收敛的状况也不用紧张,多从手册剖析处理疑问。壳单元许多,最著名的要属MITC了,bathe在84年开发的,影响至今还在,基本上如今一切的壳单元(除了degenerated solid method)都是根据MITC的。壳单元(理论)能够说是近100年来固体力学工作者最心旷神往的当地,因为起广泛的运用使得我们不得不对器感爱好,小到人的眼球,大到房顶,重要到导弹潜艇,处处都是壳啊。。。。壳太杂乱了,以至于有许多人尽量避开他去研讨板,要说如今最高深的板壳理论仍是钱伟长大师的非线性板壳理论(自己了解),可惜在如今西方的文献里,我们都只提von karman nonlinear plate/shell theory, 钱老的研讨能够说是提早的人类文明最少50年(仍是自己了解),因为就如今来说nonlinear板壳的运用仍是百里挑一。将nonlinear shell/plate理论变成有限元单元的,如今笔者还没看到,等候某位大牛呈现将其编出。
4.4 构造单元,变量:u1,u2,u3,θ1,θ2,θ3
常用于构造核算,如弯曲,剪切变形。也即是资料力学里边各种梁的核算,这儿pde方程是关于u的4阶grediant,解的办法是将二阶gradiant of u 看作是Moment,3阶看作是shear stress,这类单元通常是2节点单元,一头一尾相互连接。
索单元和板单元就不赘述了,能够参看有关资料,与构造单元迥然不一样。索单元乃至更简略,和热单元相似。还有一些薄膜单元点单元啥的运用很少,就不多说了。
看到论坛上许多人开口闭口二次开发,本来所谓二次开发很大程度上都集中于新单元的开发,假如你有新的pde,要想核算完成它,就有必要写出新的单元代码,然后完成核算。别的方面的二次开发无外乎一些后处理或是mesh上面的新功用,而这些商用有限元基本老练,所以再提二次开发之前,先想想现有的单元库里边有没有我需求的单元,假如没有再二次开发它,而这个进程是绵长而检苦的。

 

 

Chapter 6 有限元的未来
基本上有限元数值办法的最中心内容现已差不多说完了,对应不一样的疑问我们题出了许多极好的办法,但这些办法都是根据我之前所述的理论根底之上的。接下来谈谈我对有限元将来开展趋势的观点,或许能和计划或正在从事fem的兄弟有些共识。
6.1 并行核算
跟着多核cpu进入自己电脑商场,并行软件现已是大势所趋,那么并行fem软件作为大型科学核算软件,并行趋势极为明朗,各大fem公司分分研发推出并行版,可是并行版的核算准确性与功率还需求经受检测,因为并行机制的设定,每个cpu返回的核算成果不一样步,会致使latency,矩阵的blocks分化怎样树立行至有用的ghost node/halos, 削减communication time, 关于各种并行构架的适应性,都是大型有限元并行化要处理的疑问,不过当Nvidia公司推出根据Fermi构架的GPU今后,科学核算界为之振奋,但其是否能催生出有用的并行通用有限元程序仍是一个挑战。并且如今自己电脑上并行核算的speedups并不能进步许多,有2-3倍就现已很不错了。所以并行fem仍需阅历检测,但必定是大势所趋。
6.2 多规范模型核算
高性能显微设备的诞生,我们关于微观规范的物质越发感爱好,fem能否描绘微观国际呢?纳米级资料,微流场,原子的diffusion,dislocation等等,量子力学因为思考的电子的效果,使得核算量过于巨大,分子动力学办法尽管已获得了必定成功,可是终归不如fem来的有用和便利。描绘不一样规范下的资料,俨然以变成很多专家们研讨的首要方向,在坚凯围棋的第6版有限元书中添加了此方面内容,也有不少专家联络md和mc的办法,来树立多规范仿真体系,都有一些开展。
动网格,鸿沟元,无网格等等都是fem的开展或联络方向,本人不是很了解,不赘述,期待弥补纠正。