对流扩散特征值问题在流体力学、环境科学、能源开发和电子科学等多个领域的应用使它引起了科学界的广泛关注, 例如文献[1, 2]研究了对流扩散方程的格式及其解的估计. 利用有限元方法求解对流扩散特征值问题也引起了越来越多学者的研究兴趣. 例如讨论对流扩散特征值问题的后验误差估计、讨论对流扩散特征值问题的多水平校正方法, 讨论自适应算法等等. 因此, 有限元方法解决这一问题成为引起数学和物理领域关注的重要课题. 自适应有限元方法是科学计算的主流, 近年来, 这种方法已得到了深入研究, 并广泛地应用于许多问题上, 如文献[3, 4]讨论了后验误差估计和自适应算法. 使用后验误差估计和自适应有限元算法的想法是在1978年由Babuska和Rheinbolt首次在文献[5]提出的. Reed和Hill首次针对文献[6]中的线性双曲问题引入并分析了的间断有限元方法, 间断有限元方法的主要特点是测试函数沿网格中的面(或边)不连续, 具有局部质量守恒、易于与其他方法相结合、$ hp $自适应、处理多边形网格等优点. 因此, 间断有限元方法被用于解决各种特征值问题, 如拉普拉斯特征值问题、经典的自伴随Steklov特征值问题、双调和特征值问题等. 本文首次使用间断有限元方法计算对流扩散特征值问题, 建立了一个后验误差估计, 并验证间断有限元方法特征函数后验误差估计的可靠性和有效性, 实现自适应计算. 结果表明, 该自适应算法能达到最优收敛阶数, 从误差曲线也可以看出, 在相同自由度下, 自适应算法得到的近似比均匀网格计算的近似更准确.
设$ \Omega \subset {R^2} $是一个Lipshitz边界$ \partial \Omega $的有界域, 设$ \textbf{n} $是$ \partial \Omega $的单位外法向量, 考虑Dirichlet边界条件特征值问题: 求$ \lambda \in C $和$ u \in H_0^1(\Omega ) $, 使得
令$ (u, v) = \int_\Omega {u\bar vdx, } $并定义连续的双线性形式
假设$ \textbf{b} $和$ c $是$ \Omega $上的有界函数, $ \nabla \cdot \textbf{b} $存在且满足$ - \frac{1}{2}\nabla \cdot \textbf{b} + c \geqslant 0, \; \; in\; \; \Omega . $在这些假设下, 存在与$ u $、$ v $无关的两个正常数$ A $和$ B $, 使得双线性形式$ a( \cdot , \cdot ) $满足
$ (2.1) $的弱形式是求$ (\lambda , u) \in C \times {H_0^1}(\Omega ), \; \; u \ne 0, \; \; $使得下面等式成立
设$ {{\mathcal {T}}_h} = \left\{ \kappa \right\} $为$ \Omega $的形状规则网格, 单元$ \kappa $中边的长度用$ {h_e} $表示, 单元$ \kappa $的直径用$ {h_\kappa} $表示, 并且$ h = \mathop {\max }\limits_{_{\kappa \in {{\cal T}_h}}} {h_\kappa } $. $ {\Gamma _h} = \Gamma _h^i \cup \Gamma _h^b $, 其中$ \Gamma _h^i $表示内部边, $ \Gamma _h^b $表示边界$ \partial\Omega $上的边. 定义$ v $在$ e $上的均值和跳跃$ \left\{ v \right\} = \frac{1}{2}({v^ + } + {v^ - }), \left[\kern-0.15em\left[ v \right]\kern-0.15em\right] = \left( {{v^ + } - {v^ - }} \right){\mathbf{n}}, $其中$ e = \partial {\kappa ^ + } \cap \partial {\kappa ^ - }, \; \; {v^ + } = v{|_{{\kappa ^ + }}}, {v^ - }{\text{ = }}v{|_{{\kappa ^ - }}} $, $ \textbf{n} $是从$ {\kappa ^ + } $到$ {\kappa ^ - } $的单位外法向量. 如果$ e \in \Gamma _h^b $, 定义$ v $在$ e $上的均值和跳跃$ \left\{ v \right\} = v, \; \; \left[\kern-0.15em\left[ v \right]\kern-0.15em\right] = v\textbf{n}. $定义
其中$ \sigma $为惩罚参数, $ \theta $的公共值为$ 0 $, $ \pm 1 $. 定义间断有限元空间
其中$ {{\mathbb{P}_m}(\kappa )} $是$ \kappa $上的$ m $次多项式空间. 引入剖分$ {\mathcal {T} _h} $上的分片函数空间$ {H^s}({\mathcal {T} _h}) = \left\{ {v \in {L^2}(\Omega ):v{|_\kappa } \in {H^s}(\kappa ) , \; \; \forall \kappa \in {\mathcal {T} _h}} \right\}. $ $ (2.4) $的有限元近似是求$ ({\lambda _h}, {u_h}) \in C \times {S^h} $, $ {u_h} \ne 0 $, 使得
$ (2.4) $的源问题为: 求$ w \in H_0^1(\Omega ) $, 使得
$ (2.7) $的间断有限元近似是求$ {w_h} \in {S^h} $, 使得
定义线性有界算子$ T:{L^2}(\Omega ) \to H_0^1(\Omega ) $满足
由$ (2.6) $可定义对应的离散解算子$ {T_h}:{L^2}(\Omega ) \to {S^h} $满足
引入赋予间断有限元范数的和空间$ V(h) = {S^h} + H_0^1(\Omega ) $, 其中间断有限元范数为:
并且在分片函数空间$ {H^{1 + s}}({{\cal T}_h})(s > \frac{1}{2}) $上定义$ h $范数为:
注意在间断有限元空间$ {S^h} $上, $ {\left\| \cdot \right\|_G} $与$ {\left\| \cdot \right\|_h} $是等价的.
由文献[3]和格林公式可以推导出间断有限元方法的一致性, 即设$ w $为$ (2.7) $的解, $ f \in {L^2}(\Omega ) $, 有
由上式和$ (2.8) $式, 可以得到误差公式为:
不难看出, 如下连续性和椭圆性成立:
令$ w $是$ (2.7) $的解, 且$ f \in {L^2}\left( \Omega \right) $, 则成立如下正则性估计$ {{{\left\| w \right\|}_{1 + r}} \lesssim {{\left\| f \right\|}_{0, \Omega }}}{\left( {\tfrac{1}{2} < r \leqslant 1} \right)} . $由文献[7]中$ 2.1.4 $节, 可知$ {w^I} \in {S^h} $表示$ w $在$ {{\cal T}_h} $上的插值, 且有如下插值估计$ {\left\| {w - {w^I}} \right\|_h} \lesssim {h^r}{\left\| w \right\|_{1 + r}}, $注意$ \left[\kern-0.15em\left[ {{w^I} - w} \right]\kern-0.15em\right] = 0. $
定理2.1 设$ w $和$ {w_h} $分别是$ (2.7) $式和$ (2.8) $式的解, 那么有以下不等式成立
证 首先, 我们证明$ (2.17) $式, 通过利用$ (2.16) $式, $ (2.14) $式和$ (2.15) $式, 可以推导出
利用三角不等式, $ (2.19) $式, 可以得到
因此, 对于足够小的$ h $, 可以得到$ (2.17) $式.
下面, 我们证明$ (2.18) $式. 由$ (2.16) $式和$ (2.5) $式, 得到
由插值估计, 间断有限元范数的定义和正则性估计, 可得
因此, 结合$ (2.21) $式和$ (2.22) $式, 可得到
由上式, 可得
利用三角不等式和插值误差估计得$ (2.18) $式, 证明完成.
定理2.2 设$ w $和$ {w_h} $分别是$ (2.7) $式和$ (2.8) $式的解, 那么有以下不等式成立
证 考虑$ (2.4) $式对偶问题的源问题$ a(v, {w^*}) = (v, g), \; \; \forall v \in H_0^1(\Omega ) $, 对于任意固定的$ g \in {L^2}(\Omega ) $, 有正则性估计$ {\left\| {{w^*}} \right\|_{1 + r}}{\text{ }} \lesssim {\left\| g \right\|_{0, \Omega }} $成立.
由$ \left[\kern-0.15em\left[ {{w^*} - {w^{{\text{*}}I}}}\right]\kern-0.15em\right]{\text{ = }}0 $和误差公式$ (2.14) $推导出
利用插值估计, 正则性估计和$ h $范数的定义, 有
将$ (2.28) $代入$ (2.27) $, 利用Riesz表示定理可得到$ (2.25) $式.
下面, 我们证明$ (2.26) $, 由插值估计, 间断有限元范数的定义和正则性估计, 得
将$ (2.29) $代入$ (2.27) $, 利用Riesz表示定理, 可以得到$ (2.26) $, 证明完成.
设$ \lambda $是$ (2.4) $的第$ j $个特征值, 具有代数重数$ q $和陡度$ \alpha $, 其中$ {\lambda _j} = {\lambda _{j + 1}} = \cdots = {\lambda _{j + q - 1}} $. 当$ {\left\| {{T_h} - T} \right\|_{0, \Omega }} \to 0 $, $ (2.6) $的$ q $个特征值$ {\lambda _{j, h}}, \cdots {\lambda _{j + q - 1, h}} $将收敛到$ \lambda $. 设$ M(\lambda) $是与$ \lambda $相关的$ (2.4) $广义特征向量空间, $ {M_h}(\lambda ) $是与$ {\lambda _h} $相关的$ (2.6) $广义特征向量空间的直接和, $ {\lambda _h} $收敛于$ \lambda $.
类似文献[3]中定理3.1的证明方法可以证明下面定理.
定理2.3 设$ M(\lambda ) \subset {H^{1 + s}}(\Omega )(m \geqslant s > \frac{1}{2}) $, 那么有以下不等式成立
设$ {u_h} \in {M_h}(\lambda ) $是$ (2.6) $的广义特征向量空间的直接和, 那么存在$ (2.3) $的特征值$ u $使得
如果设$ \alpha = 1 $, 那么
设$ \left( {{\lambda _h}, {u_h}} \right) $为$ (2.6) $的特征对, 在每个单元$ \kappa \in {{\cal T}_h} $和$ e \in {\Gamma _h} $上分别定义元素残差和面残差为$ {R_\kappa } = \Delta {u_h} + \left( {{\lambda _h} - c} \right){u_h} - {\mathbf{b}} \cdot \nabla {u_h};\; \; {J_{F, 1}} = {\text{ }}\left[\kern-0.15em\left[ {\nabla {u_h}} \right]\kern-0.15em\right], \; \; \forall e \in {\Gamma _h^i};\; \; {J_{F, 2}} = \left[\kern-0.15em\left[ {{u_h}} \right]\kern-0.15em\right], \; \; \forall e \in {\Gamma _h^i}. $
定义每个单元$ \kappa \in {{\cal T}_h} $上的局部误差指示子
全局误差指示子为
下面, 我们将证明这个误差估计是可靠的.
引进提升算子$ {\cal L}:V(h) \to {\left[ {{S^h}} \right]^2} $如下
由文献[8]我们知道提升算子具有稳定性, 即$ \left\| {{\cal L}(v)} \right\|_{0, \Omega }^2 \lesssim \sum\limits_{e \in {\Gamma _h}} {\left\| {h_e^{ - 1/2}\left[\kern-0.15em\left[ v \right]\kern-0.15em\right]} \right\|_{0, e}^2}. $利用这个算子, 我们可以定义一个辅助双线性型式$ {\tilde a_h}( \cdot , \cdot ):V(h) \times V(h) \to C $,
由文献[3]中定理$ 4.1 $可知下列定理成立.
定理3.1 设$ (\lambda , u) $和$ ({\lambda _h}, {u_h}) $分别是$ (2.4) $和$ (2.6) $的特征对, $ u \in {H^{1 + s}}(\Omega )(s > \frac{1}{2}) $, 那么对于任意的$ v \in H_0^1(\Omega ) $, 有下式成立
证 注意在$ H_0^1\left( \Omega \right) \times H_0^1\left( \Omega \right) $上$ a = {\tilde a_h} $. 设$ w \in H_0^1\left( \Omega \right) $, 由双线性形式的椭圆性和连续性, 可推导出
取$ v = u - w $, 可得
由三角不等式, 得
由$ w $的任意性, 可得定理成立.
由文献[9]中Scott-Zhang插值, 我们可以得到以下引理.
引理3.1 对于任意的$ \varphi \in H_0^1(\Omega ) $, 存在一个分片线性插值$ {I^h}\varphi \in {S_h} $满足
其中$ {{U_\kappa }} $是与$ \kappa $共享至少一个节点的所有单元的并集, $ {{U_e}} $是与边$ e $共享至少一个节点的所有边的并集.
引理3.2 设$ (\lambda , u) $和$ ({\lambda _h}, {u_h}) $分别是$ (2.4) $和$ (2.6) $的特征对, 对于任意的$ v \in H_0^1(\Omega ) $, 成立
证 由插值性质, 我们得到$ \left[\kern-0.15em\left[ {v - {I^h}v}\right]\kern-0.15em\right] = 0 $, 利用格林公式, 可得到
由柯西-斯瓦兹不等式、$ (3.9) $式、$ (3.10) $式和提升算子的稳定性, 有
结合$ {B_1} - {B_4} $上述四个不等式, 可以得到
由$ (2.6) $和$ (3.15) $我们有
故$ (3.11) $成立, 证明完成.
定理3.2 在定理$ 3.1 $的条件下, 有下式成立
证 由文献[10]可知, 对于任意的$ v \in {S^h} $, 存在一个丰富算子$ {E_h}:{S^h} \to {S^h} \cap H_0^1(\Omega ) $, 使得
对$ (3.5) $右边第二项利用$ (2.11) $和$ (3.18) $式, 并注意到$ \left[\kern-0.15em\left[ {{E_h}{u_h}}\right]\kern-0.15em\right] = 0 $, 有
将$ (3.19) $和$ (3.11) $代入$ (3.5) $, 可得到$ (3.17) $, 即证明完成.
由定理2.3, 我们知道陡度$ \alpha=1 $时, $ {\left\| {\lambda u - {\lambda _h}{u_h}} \right\|_{0, \Omega }} $和$ {\left\| {u - {u_h}} \right\|_{0, \Omega }} $都是$ {\left\| {u - {u_h}} \right\|_G} $的高阶小量, 因此$ (3.17) $告诉我们误差估计指示子$ \eta ({u_h}) $是间断有限元能量范数的上界之一, 因此误差估计是可靠的.
为了保证我们的估计方法对于实际的自适应改进是有效的, 我们的下一个目标是证明局部误差估计指示子$ {\eta _\kappa } $提供了$ \kappa $上误差的局部下界. 标记$ {b_\kappa } \in H_0^1(\kappa ) $为标准单元气泡函数, $ {b_e} \in H_0^1({U_e}) $为面上的气泡函数, 其中$ {U_e} $是两个单元$ {\kappa ^ + } $和$ {\kappa ^ - } $共享$ e $的并集, 通过利用Verf$ \ddot u $rth开发的气泡函数技术, 我们引入并介绍以下知识.
引理3.3 对于所有多项式函数$ v \in {P_k}(\kappa ) $,
对于所有多项式函数$ w \in {P_k }(e) $, 我们有
对于每一个$ {b_e}w $, 存在扩展的$ {W_b} $满足$ {W_b}{|_e} = {b_e}w $, $ {W_b} \in H_0^1({U_e}) $
根据以上引理, 并使用标准参数(见文献[11]中引理3.13), 我们可以证明以下有局部下界.
引理3.4 设$ (\lambda , u) $和$ ({\lambda _h}, {u_h}) $分别是$ (2.4) $和$ (2.6) $的第$ j $个特征对, 然后我们有以下局部下界:
$ (i) $对于任意的$ \kappa \in {{\cal T}_h} $, $ {h_\kappa }{\left\| {\Delta {u_h} + \left( {{\lambda _h} - c} \right){u_h} - {\mathbf{b}} \cdot \nabla {u_h}} \right\|_{0, \kappa }} \lesssim {\left\| {\nabla (u - {u_h})} \right\|_{0, \kappa }} + {h_\kappa }{\left\| {u - {u_h}} \right\|_{0, \kappa }}\\ + {h_\kappa }{\left\| {\lambda u - {\lambda _h}{u_h}} \right\|_{0, \kappa }}. $
$ (ii) $设$ e \in \Gamma _h $是$ {\kappa ^ + } $和$ {\kappa ^ - } $共享的内部边, 有
其中$ { U _e} = \{ {\kappa ^ + }, {\kappa ^ - }\} $.
$ (iii) $对于每条边$ e \in \Gamma _h $, $ h_e^{ - 1}\left\| {{J_{F, 2}}} \right\|_{0, e}^2 = h_e^{ - 1}\left\| {\left[\kern-0.15em\left[ {{u_h}} \right]\kern-0.15em\right]} \right\|_{0, e}^2 = h_e^{ - 1}\left\| {\left[\kern-0.15em\left[ {u - {u_h}} \right]\kern-0.15em\right]} \right\|_{0, e}^2. $
证 $ (i) $设$ {v_h} = \Delta {u_h} + \left( {{\lambda _h} - c} \right){u_h} - {\mathbf{b}} \cdot \nabla {u_h} $和$ {v_b} = {b_\kappa }{v_h} $. 注意在$ {L^2}(\kappa ) $中$ \Delta u + \left( {\lambda - c} \right)u - {\mathbf{b}} \cdot \nabla u = 0 $, 在$ \partial\kappa $上, $ {v_b} = 0 $, 利用分部积分, 有
利用$ (3.20) $式和柯西-斯瓦兹不等式, 可得
则$ (i) $得证.
$ (ii) $对于任意的$ e \in \Gamma _h $, 设$ {w_h} = \left[\kern-0.15em\left[ {\nabla {u_h}}\right]\kern-0.15em\right] $, $ {w_b} = {b_e}{w_h} $. 设$ {W_b} \in H_0^1({U_e}) $是$ {w_b} $满足$ (3.22) $和$ (3.23) $的扩展. 注意$ \left[\kern-0.15em\left[ {\nabla u}\right]\kern-0.15em\right] = 0 $, 利用格林公式, $ (3.21) $, $ (3.22) $, $ (3.23) $, 可以推导出
结合$ {\left\| {\Delta {u_h} + \left( {{\lambda _h} - c} \right){u_h} - {\mathbf{b}} \cdot \nabla {u_h}} \right\|_{0, \kappa }} $在$ (i) $中的界和网格形状规则性得
由此可得到$ (ii) $.
$ (iii) $对于任意的$ e \in \Gamma _h $, 我们有$ \left[\kern-0.15em\left[ u\right]\kern-0.15em\right] = 0 $, 从而得到$ (iii) $.
定理3.3 在定理$ 3.1 $的条件下, 有下面等式成立
证 通过$ {\eta _\kappa } $的定义和引理3.4, 可得到$ (3.26) $, 再利用间断有限元范数$ {\left\| \cdot \right\|_G} $的定义, 可以得到$ (3.27) $.
定理3.3表明误差估计指示子$ \eta ({u_h}) $是有效的.
引理3.5(文献[3]中引理4.6) 设$ \left( {\lambda , u} \right) $和$ \left( {{\lambda _h}, {u_h}} \right) $分别是$ (2.4) $和$ (2.6) $的特征对, 设$ \left( {{\lambda ^*}, {u^*}} \right) $和$ \left( {\lambda _h^*, u_h^*} \right) $分别是$ (2.4) $的对偶问题及离散变分形式的特征对, $ \left( {{u_h}, u_h^*} \right) \ne 0 $, 那么
定理3.4 在引理$ 3.5 $的的条件下, 设特征函数空间$ M(\lambda ), \; \; M({\lambda ^*}) \subset {H^{1 + s}}(\Omega )(s > \frac{1}{2}) $, 那么
其中, $ \varepsilon $是一个足够小的常数, $ \sigma>\frac{1}{2} $充分接近$ \frac{1}{2} $.
证 定理$ 2.3 $表明$ {\left\| {u - {u_h}} \right\|_{0, \Omega }} $比$ {\left\| {u - {u_h}} \right\|_h} $更高阶, $ {\left\| {{u^*} - u_h^*} \right\|_{0, \Omega }} $也比$ {\left\| {{u^*} - u_h^*} \right\|_h} $更高阶. 因此, 由$ (3.28) $, $ {u_h} $的估计式$ (3.17) $和$ {u_h^*} $的估计式, 可得
由逆估计和迹不等式, 得
通过以上三个公式, 可证明该定理.
从定理$ 3.2 $和定理$ 3.3 $, 可以知道特征函数误差$ \left\| {u - {u_h}} \right\|_G^2 + \left\| {{u^*} - u_h^*} \right\|_G^2 $的估计指示子$ \eta {({u_h})^2} + \eta {(u_h^*)^2} $是可靠和高效的, 因此, 基于该估计指示子的自适应算法可以生成良好的梯度网格. 使近似特征函数在$ \left\| \cdot \right\|_G^2 $中达到最优收敛阶$ O(do{f^{ - m}}) $. 因此, 我们期望得到:
因此, 从$ (3.29) $可得$ \left| {\lambda - {\lambda _h}} \right| \le do{f^{ - m}} $. 所以$ \eta {({u_h})^2} + \eta {(u_h^*)^2} $可看作$ {\lambda _h} $的误差估计指示子, 第4节的数值实验表明$ \eta {({u_h})^2} + \eta {(u_h^*)^2} $作为$ {\lambda _h} $的误差估计指示子是可靠的和高效的.
在本节中, 将报告一些数值实验, 以此来证明我们方法的有效性. 考虑问题$ (2.1) $, 其中取$ \textbf{b}=(0, 0)^{T} $, $ (1, 1)^{T} $, $ (3, 0)^{T} $, $ c=0 $. 我们的程序是在iFEM软件包下编译的, 我们使用SIPG方法$ (\theta=-1) $来进行计算. 考虑以下两个测试域: L形域$ {\Omega _L} = {( - 1, 1)^2}\backslash ([0, 1) \times ( - 1, 0]) $, 裂缝结构域$ {\Omega _{SL}} = {( - 1, 1)^2}\backslash {\left\{ {0 \leqslant x \leqslant 1, y = 0} \right\}} $. 由于确切的特征值是未知的, 我们在L形域中取参考特征值$ \lambda_1={\left| \textbf{b}\right|^2}/4 +9.63972384472 $; 在裂缝结构域$ {\Omega _{SL}} $中取参考特征值$ \lambda_1={\left| \textbf{b}\right|^2}/4 +8.3713297112 $. 这些参考特征值是通过自适应计算尽可能精确得到的.
我们在表 1到表 3中列出了通过自适应计算得到的特征值数值解结果, 并在图中描述了自适应网格和误差曲线. 从图 1到图 3中, 我们可以看出, 当$ \textbf{b}=(0, 0)^{T} $, $ (1, 1)^{T} $, $ (3, 0)^{T} $时, 一次间断元的误差曲线近似平行于斜率为$ -1 $的直线. 结果表明, 该自适应算法能达到最优收敛阶数, 从误差曲线也可以看出, 在相同自由度$ (dof) $下, 自适应算法得到的近似比均匀网格计算的近似更准确.