FEM之理论通俗有限元（2）---偏微分方程

先来看一个平面热传导的例子：

该平板左下方的顶点温度设置100度，右上方顶点设置0度，其它边自由。

图1是有限元模型图，图2显示的是温度场的分布，图3是温度梯度向量图。

在上面的例子中，右下方的顶点也加上一个温度为0的约束。

问题来了：

1. 为什么温度场会服从这种分布，而不是别的形式？

2. 温度场的分布是否是唯一的？如果是，为什么？

3. 对于这类问题如何求解？

偏微分方程（Partial Differential Equation），简称PDE

偏微分方程最早是科学家从振动，热传递等一系列的物理现象归纳总结出来的，实际上绝大部分物理现象都可以用PDE来描述。最常用的PDE都是二阶的，即场的导数是二阶，高阶的PDE也可以通过变换降成二阶PDE。

首先介绍导数来预热， 我们知道 把位移对时间求导，就是速度，把速度对时间再求导就是加速度，简写为 a=s(t)'';

a为加速度，S为位移，t为时间；

这个表达式实际上有两层含义：

1. 位移是时间的函数。例 S(t) = t*t*t-5*t 4， S(t) = sin(t)/cos(2*t);

2. 要使加速度不为0， S(t)一次求导以后必须不为常数，因为常数求导为0。

PDE中要求解的是场函数的表达式，即上面中的S(t)，知道了S(t)的表达式，给出t，就能求出任意时刻的S。

下面是一个典型的二阶PDE，其中u为场函数

方程看起来很简单，其实不管是通解，还是定解，u的解析解表达式都是很复杂的，通常都是高阶函数，正弦，余弦，指数，对数的组合表达式。理论上PDE都可以找到解析解，但现实大部分PDE的解析解都无法求出。找解析解是数学家们的任务。

该PDE说明u是含有x，y的函数，当对x求导时，y看做常量，对y求导时，x看做常量。

假设 u(x,y)=x*x-y*y x    ---(1)

u对x求导一次:u'/x=2*x 1; 

u对x求导2次:u''/x=2;

同理对y求导2次:u''/y=-2

说明式1是方程的一个解，但问题是没有给定边界条件，方程的解实际上是有无限多个的，要确定定解，必须要设定边界条件，这也是例子中右下方加上一个新约束，导致温度场发生变化：即使边界一个小的变化也可能导致完全不同的定解，即不同的u表达式。

关于PDE 的边界条件和常用PDE参考：三角单元介绍

------------------------------------------------------------------

上述方程是无热源稳态热传导的PDE，场函数为温度，称为拉普拉斯方程（Laplace，名人得记住），拉普拉斯方程式是电磁，热，流体，天文等领域经常会碰到的一类数学问题。详细介绍可参考PDE方面的书。

从上面两个例子也可以看出，即使很简单的PDE也是很难找到精确的解析解。找到解析解固然最好，但是找不到也没关系。就跟圆周率一样，我们没必要知道它精确多少位，取7位数对于绝大部分应用足矣。找不到解析解，就只能依靠数值解了，数学家们首先用到的是 差分法。

任然以平面正方形温度场为例，假设左边，上，下边温度为100度，右边为0度。 

差分法就是在求解区域画上 一定密度的结构化网格，如下图：

上图中的任意一个点的温度值都可以用其周围临近的四个点的平均值表示，如果是长方形，表达式会略有不同。将网格简化如下：

中心点的温度为 （100 100 100 0）/4 =75 

若将网格加密一些，如上图，中心四个点为未知量，每个点用周围四个点的和求平均，可得四个方程，联立解线性方程组即可，这是差分法的基本思想。

在差分法的基础上，发展起来的是有限元法。

有限元方法的理论基础是变分原理和加权余量法。通常有限元方法，变分法，加权余量，以及编程实现之间的关系很容易让人迷惑。

关于有限元方法的原理，建议直接记住一点：就是伽辽金法(Galerkin)。虽然还有变分，Razi等，但用Galerkin来理解有限元足矣。Galerkin方法是加权余量法的一种。加权余量法的基本原理是：

假设现在有PDE的近似解，而且这个解也满足边界条件，那么将这个解带入PDE，计算结果将产生误差（因为解场函数是近似的），为了在整个求解区域上保持与精确解相同的特性，就必须想办法将这些误差在近似解内部消除掉。

常用加权余量法有：子域法，最小二乘法法，配点法，伽辽金法。使用主要涉及到两个方面：

1. 如何设计近似解。

理论上可以设计出无限多个近似解，指数，幂数，正弦，余弦等。但是教科书上无一例外的都是以多项式和为主。原因只有一个：简单。

2. 如何消除误差。即采用什么方法，强制设置误差为0，以满足方程的解。

伽辽金法要求误差相对于权函数是正交的，即两者的乘积积分为0。通常权函数为假设的解多项式。将求解区域离散后，不同的单元通过节点连接，节点所连接的每一个单元在该点产生的误差与设计的权函数乘积积分之后相加为0（有点拗口）.简单理解就是强制一个节点连接的所有单元在该点产生的误差和为0。由此可见，一个节点就可以产生一个方程（严格来说是一个自由度产生一个方程，因为温度是标量，一个节点只有一个自由度，如果取位移，平面一个节点就是2个自由度，X和Y两个方向），最终的解就是节点个数的多项式。

总之任何形式的PDE，离散成网格之后，都可以通过伽辽金方法得到一组近似的线性方程组K*x=F。而且利用伽辽金方法最大的好处就是K是对称矩阵，方便计算求解。这也是采用伽辽金方法建立有限元方程的主要原因。

说到这里，利用有限元求解PDE原理很简单：设计近似解，产线性方程组，求解方程组。产生线性方程组的方法有很多，最常用的就是伽辽金法。当微分方程存在相应的泛函时，伽辽金法和泛函可以推导出相同的结果。

如果说差分法和有限元法有什么区别的话，从网格角度看：差分法是把网格画好然后铺在求解对象上，而有限元法是直接把对象划分网格。因此差分法适合求解区域比较规则的几何，而有限元可以对任意形状几何求解。

除伽辽金法外，流体中常用有限体积法，有限体积法采用加权余量法中的子域法。

现在实际上已经回答了前面提出的三个问题：

1. 为什么温度会服从这种分布，而不是别的形式的？---服从PDE

2. 温度场的分布是否是唯一的？如果是，为什么？ ---唯一的边界条件对应唯一的分布

3. 对于这类问题如何求解？---差分法，有限元法，有限体积法

想用通俗的语言来描述有限元，发现其实很难，因为很多东西是需要公式的，有限元本来就是利用积分，微分来解决问题，脱离这些东西讲FEM只能算是科普了。

声明：原创文章，欢迎留言与我讨论，如需转载留言