[有限元方法基础理论] 质量集中有限元（谱元）

[有限元方法基础理论] 质量集中有限元（谱元）

集中质量有限元，也称谱元，要求具有对角质量矩阵，且对角线上的值等于行和。可以通过用数值积分公式以及其他方法来构造谱元。

不妨考虑一个初边值问题，
$$\begin{array}{rl}& u_t-\Delta u=f\text{ in }\Omega ,t>0\\ & u=0\text{ on }\partial \Omega ,t>0,\text{ with }u(\cdot ,0)=v\text{ in }\Omega \end{array}$$
它的半离散方法可以写为：寻找$u_h:[0,\infty )\to S_h$，使得
$$\left(u_{h,t},\chi \right)+\left(\nabla u_h,\nabla \chi \right)=(f,\chi ),\forall \chi \in S_h,t>0,u_h(0)=v_h$$
写成矩阵的形式，就是
$$\mathcal{B}{\alpha }^{\prime }(t)+\mathcal{A}\alpha (t)=\tilde{f}(t),\text{ for }t>0,\text{ with }\alpha (0)=\gamma$$
这里的$\mathcal{B}$是质量矩阵，$\mathcal{A}$是刚度矩阵。以上就是标准的有限元离散到代数方程组，用到的一些符号的具体的含义不言自明，就不细说了。

所谓的质量集中有限元，我们就是希望，这里得到的质量矩阵$\mathcal{B}$，满足条件：它是个**对角矩阵，对角线上的值刚好等于行和。**不妨用$\overline{\mathcal{B}}$来表示。

那么问题来了，通过什么样的方式，上述的要求才能被满足呢？我们这里提两种方法，一个是通过数值积分的方法构造，另外一个是通过构造分片常数有限元空间的方法来构造。实现的时候，总是方便的，但是在收敛性分析的时候，我们希望有它的一种构造方式作为桥梁。

采用数值积分公式计算单元质量矩阵时，当求积点与插值点一致时，就可以得到谱元。积分公式要求：达到指定阶，且权重全正。

对于分片 $k$ 次多项式有限元，精确计算质量矩阵的积分 公式需安达到 $2k$ 阶。对常见二阶 PDE，积分公式只需要达到 $2k-2$ 阶便可保 证数值解精度。

考虑最简单的线性元 $P_1$插值点为三个顶点 { v 0 , v 1 , v 2 } \left\{v_{0}, v_{1}, v_{2}\right\} {v0​,v1​,v2​}， 将其作为求积点，积分公式可达到 1 阶且权重 $\left(\frac{1}{3}\right)$ 为正，满足要求。下面我们来验证一下，它满不满足质量集中的条件。

在单元$\tau$上，它的以顶点为积分点的数值积分公式为，
$$Q_{\tau ,h}(f)=\frac{1}{3}\text{ area }(\tau )\sum _{j=1}^{3}f\left(P_{\tau ,j}\right)\approx {\int }_{\tau }fdx$$
刚好是三个顶点的平均值再乘以面积。为了方便，我们接下来用如下符号来表示数值积分意义下的内积，
$$(\psi ,\chi {)}_h=\sum _{\tau \in {\mathcal{T}}_h}{Q}_{\tau ,h}(\psi \chi )$$
下面只要证明，对于
$${\left(u_{h,t},\chi \right)}_h+\left(\nabla u_h,\nabla \chi \right)=(f,\chi ),\forall \chi \in S_h,t>0,u_h(0)=v_h$$
代数化得到的质量矩阵，满足集中质量的要求即可。事实上，若令$u_h(t)={\sum }_{j=1}^{N_h}{\alpha }_j(t){\Phi }_j$，则代数方程组可以写为：
$$\sum _{j=1}^{N_h}{\alpha }_j^{\prime }(t){\left({\Phi }_j,{\Phi }_k\right)}_h+\sum _{j=1}^{N_h}{\alpha }_j(t)\left(\nabla {\Phi }_j,\nabla {\Phi }_k\right)=\left(f,{\Phi }_k\right),k=1,\dots ,N_h$$
只要证明两点：

• ${\left({\Phi }_j,{\Phi }_k\right)}_h=0$ for $j\ne k$，既然是0狄利克雷边界。
• ${\Vert {\Phi }_j\Vert }_h^2={\left({\Phi }_j,{\Phi }_j\right)}_h={\sum }_{k=1}^{N_h}\left({\Phi }_j,{\Phi }_k\right)$

第一点是自然的，因为我们可以看到，不同的形函数的乘积在顶点处的值都等于 0，那么数值积分值就等于 0。下面来说明第二点。

首先，需要说明的是，在同一个单元上，两个不同的形函数的之间的乘积积分等于 1/12 的三角形面积，而形函数自身和自身的乘积的积分等于 1/6 的单元面积，即
$${\int }_{\tau }{\Phi }_j{\Phi }_{k}dx=\frac{1}{12}\text{ area }(\tau )\text{ 和}{\int }_{\tau }{\Phi }_j^{2}dx=\frac{1}{6}\text{ area }(\tau )$$
这个可以通过标准单元变换和简单计算得到。可以看到，若用$D_j$表示和顶点$P_j$有关的三角形区域的并集，那么原来的那个质量矩阵的行和就表达为
$$\sum _{k=1}^{N_h}\left({\Phi }_j,{\Phi }_k\right)=\frac{1}{3}\text{ area }\left(D_j\right)$$
这个不清楚的话，画一下图，就可以看明白。一个基函数在相关的某个三角形单元上，自身和自身作用了一次，和同一三角形上其他顶点对应的基函数各作用了一次。

我们再看利用数值积分公式计算内积得到的质量矩阵对角线上的值，
$${\Vert {\Phi }_j\Vert }_h^2=\sum Q_{\tau ,h}\left({\Phi }_j^2\right)=\frac{1}{3}\text{ area }\left(D_j\right)$$

这个刚好就等于原始质量矩阵的行和，这就说明了，用这样一个数值积分格式得到的质量矩阵刚好满足质量集中的要求。

二次元 $P_2$插值点为三个顶点 { v 0 , v 1 , v 2 } \left\{v_{0}, v_{1}, v_{2}\right\} {v0​,v1​,v2​} 和三条边中点 { e 0 , e 1 , e 2 } \left\{e_{0}, e_{1}, e_{2}\right\} {e0​,e1​,e2​} ，相应的积分公式可 以达到 2 阶精度。但是，三个顶点的权重为 $0$。一般采取增加自由度或者扩展多项式空间的方式构造谱元。
三次元 $P_3$可以找到由三个顶点、每条边上两个对称的内点以及三 角形重心构成的 4 阶积分公式。但是，三个顶点的权重是负的!

下面我们看用另外一种方法来构造谱元。

对于一个三角形单元$\tau$，我们把每个顶点和它对边的中点用直线段连接，这些直线段就交于三角形的重心，并且把三角形分成了等面积的 6 块。我们用 $B_{j,\tau }$来表示六块中和顶点$P_j$相邻的两块，它的面积刚好是三角形面积的1/3。我们用$B_j$来表示顶点$P_j$周围的$B_{j,\tau }$的并区域。这个自己画个图就可以很容易理解。

那么，我们就可以以$B_j$为分割单元，构造一个分片常数有限元空间${\bar{S}}_h$，里面的元素可以表示为，
$$\bar{\chi }(x)=\sum _{j=1}^{N_h}\bar{\chi }\left(P_j\right){\bar{\Phi }}_j(x)$$
这里的${\bar{\Phi }}_j=1$是关于区域$B_j$的特征函数。

那么只针对原来的质量矩阵这一块进行处理，
$$\left({\bar{u}}_{h,t},\bar{\chi }\right)+\left(\nabla u_h,\nabla \chi \right)=(f,\chi ),\forall \chi \in S_h,t>0,u_h(0)=v_h$$
很容易看到，这样离散得到的质量矩阵，是满足集中质量的要求的。事实上，$\left({\bar{\Phi }}_j,{\bar{\Phi }}_k\right)=0$ 当 $j\ne k$。并且${\Vert {\bar{\Phi }}_j\Vert }^2=$ area $\left(B_j\right)=$ area $\left(D_j\right)/3$。

关于其误差分析，可以利用方法一提到数值积分的构造，参考Thomée V. Galerkin finite element methods for parabolic problems[M]. Berlin: Springer-Verlag, 1984.。