有限元杂谈之一 -- 有限元基础(转)
前言
写本文的初衷 接触有限元时间不算常也不短,从05年到现在吧,中间断断续续的学习有限元,大大小小 软件也用了不少,牛人的帖子和大师们的著作读了一些,最近常来simwe 逛逛,发现还有 很多初学有限元的朋友们, 提出的问题五花八门但大多不对正路, 不竟回想起了当年自己学 习有限元时情景, 孤立无助到处瞎碰, 通宵熬夜啃手册也得不出一个所以然来, 问导师问师 兄还是常常一头雾水, 现在回头想想如果有几个良师益友会是多么幸福的事情, 写此文仅仅 只为了让对有限元感兴趣的朋友们少走一些弯路,有更多的精力和时间去研究更深层次的问 题,同时也算是报答simwer多年来陪伴。
本文并不是一套完整的教科式学习材料, 对于正真想在有限元上做一些工作的, 还是建议多 读大师的著作和程序, 少看此文。 另外此文仅仅是自己对有限元方法粗浅的理解, 如有疏漏 的地方,还请各位高手大师补充指正。
什么是有限元
1.1 PDE
有限元是一种求解问题的数值方法,求解什么问题呢?--求解 PDE(偏微分方程)。 那么PDE是做什么用的呢?--描述客观物理世界。我想如果这两个问题搞清楚了也就明白了 为什么要用有限元,有限元可以做那些东西。 PDE可以描述很多物理现象,电磁,流体, 换热, 声学,扩散,相变,各种力学,河床变迁,物种竞争,股票金融,等等等等。。。。乃 至整个宇宙,当然也不是所有的物理现象都可以用PDE等,所以我一般不建议用有限元方法仿真微观物质现象的原因,但也有PDE应用于微观 物质并得到很好的结果,如泊松方程来解析plasma的物理现象,这在量子物理里用统计的 方法过于庞大, 泊松方程反而使问题简单而且能吻合实验,这些都是题外话就不多说了。除 了PDE以外, ODE同样也可以描述客观世界, ODE相对简单, 多用于控制系统,很有一 些线性PDE的解法也都是将PDE转化为ODE来求解的
1.2 求解PDE
有了PDE以后,问题是如何求解并得到结果,首先要说明的是不是所有的PDE都有解的, 往往有解的PDE才有实际工程意义。对于数值解法,常用的是有限差分,有限元,有限体 积和谱方法, 还有蒙特卡罗法。 有限差分出现的较早, 计算精度相对较高且可控, 但模型形 状必须规则, 边界条件处理困难, 好处是可以比较方便的控制计算精度, 编程简便, 固定节 点的网格划分形式适用于流体类的仿真。有限元方法效率高且满足精度要求, 边界条件容易 处理,得到了广大的应用,尤其是在固体领域。谱方法由于可以采用高阶差值方程和FFT 方法的来求解, 使得程序有着精度高, 收敛快的特点, 也克服了有限元条件下使用高阶插值 方程计算费时的缺点, 常常使用periodicboundarycondition, 但也有越来越多的算法使得一 类二类边界成为可能,适合微观尺度的PDE解,谱方法和有限元结合产生的谱元法取两者 之优点, 使得应用前景非常好。 有限体积类似是有限差分和有限元的结合, 特点是计算点在 单元中心,雷曼边界处理单元之间的交接,处理流体中的激波有着特殊优势,计算方法边界, 速度要比同等级的有限元方法快。蒙特卡罗法不是基于弱解形式的, 随机数的多维采样最终 得到统计上的结果,多用于金融分析。咱这里还是着重有限元解PDE,顾名思义,有限元 将整个计算几何模型划分为很多小的单元( element), 每个单元的含有一定数量的节点(node), 具体单个单元有多少节点, 有对应的不同算法与差值方程, 拿一个简单的线性4节 点平面单元来说, 每个单元包含4个节点, 每个节点有对应的自变量值, 比如简单固体力学 问题(位移法),每个节点就有对应的位移值,热力学问题每个节点就有对应的温度值,等 等。然后单元内部的结果就通过差值方法显示出来
1.3Galerkinapproximation 与 weakformation(变分/虚功)
弱解(weakformation)是建立在变分法基础上的,通过这个方法将 strongform的PDE转换为weakform, 使得有限元的求解成为可能,具体如何推导weakform可以参考一些有限元书 籍, 如果一本基础的有限元书籍没有介绍如何推导weakform, 那么可以考虑选择换一本了。 推导所得的弱解形式仍需要通过计算机来计算, Galerkin方法推导出含有求和符号的公式, 在计算机中多半以 loop的形式来计算这个量, 往往这个量就是stiffnessmatrix中的各个量。 公式中还会存在积分计算, 有限元方法多用gaussquaradure的方法来计算, 精度一般可以满 足。也就说一般简单的限元计算中存在两次 approximation, 一次是 Galerkin一次是 gauss 法求积分, 这也是很多人在计算完以后需要进行验证的原因,同时对于非线性问题 newton 法本身就是通过计算误差来判断收敛的, 固体有限元常常使用能量增量活residual增量作为newton求解器的判据。 单个单元的stiffnessmatrix计算完成后, 还需要将所有单元的矩阵装 配为一个大型的矩阵,然后进行线性代数计算。这个装配是很有技巧的,因为一般情况下stiffnessmatrix是一个很大的稀疏矩阵, 0值往往可以省去以节省计算量。 一个由10x10x10 个8节点矩形单元组成的模型会有11x11x11个 node, 如果是热力学问题自变量是温度T, 每个节点上有一个自由度,组装后的的 stiffnessmatrix会有(11x11x11)x(11x11x11)=1771561 之大, 随着单元数或变量的增加, 计算是惊人的, 这样一个大的矩阵求解显然不能用常规的 方法法(如gaussianelimination法),这就是各种这样求解稀疏矩阵算法存在的原因了, 稀疏矩阵存储转换方法可以参见 Bathe的《 FiniteElementProcedure》, 有效且快速求解线代方程组 是一个好的solver的标准
1.4 后处理
其实对于最基本的有限元方法, 求解得到的仅仅是各节点的自变量的值, 对于基本的位移法 如力学就是节点位移值, 热力问题就是温度值, 流体就是位移速度加压力值, 如果我们想知 道应力或者应变怎么办呢?后处理系统里面个都会增加相应子程序来计算 stress, strain,flux等等。 这也就是为何我们能看到各种各样的云图的原因了, 当然读者也可以自己加入计 算各种量的子程序, 如应变能密度什么的。 关于什么是有限元就介绍到这里, 仅仅是一些随和想法, 具体的理论和推导需要自身实践与探索, 本文行文仓促只是阐述自己对有限元的 粗浅理解。有不对的地方还请指正。下一章会谈及一些我曾经用过得软件。
查看更多评论 >