弱Galerkin有限元方法^[1-3]是一种新型的有限元方法.由于引入了弱梯度的概念，该离散方法具有很高的灵活性^[1-2].目前，弱Galerkin有限元方法已被成功地用于求解二阶和四阶椭圆问题及Stokes问题等^{[2, 4-5]}，但其离散所得代数方程组的数值解法研究还比较有限.

另一方面，瀑布型多重网格方法是一类易于实施的单向多重网格方法.对于二阶椭圆问题和Stokes问题，文献[6-10]构造了基于有限元或有限体积法离散的瀑布型多重网格方法，并证明当采用适合的光滑迭代子时，方法在计算精度和计算复杂度方面均是最优的.

本文首先对于弱Galerkin有限元空间，提出并分析了一类新的网格转移算子，然后基于文献[8]的瀑布型多重网格的理论框架，证明了当采用共轭梯度法作为迭代算子时，基于弱Galerkin有限元离散的瀑布型多重网格算法具有最优的计算精度和计算复杂度.

1 模型问题和弱Galerkin有限元方法

设Ω是R²中有界凸多边形区域，其边界为∂Ω.考虑如下椭圆问题：

$ \left\{ \begin{array}{l} -\nabla \cdot (a(x)\nabla u) = f, \;\;在\mathit{\Omega 内}\\ u = 0, \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;在\partial \mathit{\Omega }上 \end{array} \right. $

(1)

其中，f∈L²(Ω)，a(x)∈W^{1, ∞}(Ω) 且存在正数α和β，满足α≤a(x)≤β.这里及后文采用与文献[11-12]相同的Lebesgue空间记号、Sobolev空间记号和相关范数记号.对于L²空间，本文省略内积和范数的下标.用c或C表示与网格参数无关的正常数，且在不同地方可取不同值.

设T_h是正则且拟一致的有限元网格剖分，并用下标h表示其网格尺度. j、j₁和r为非负整数.对于任意K∈T_h, ∂K和$\mathop K\limits^ \circ $分别表示其边界和内部. P_j($\mathop K\limits^ \circ $) 表示$\mathop K\limits^ \circ $上次数小于等于j的多项式全体，P_j1(∂K) 表示∂K上次数小于等于j₁的多项式.令W(K):={{v₀, v_b}:v₀∈L²(K), ${v_b} \in {H^{\frac{1}{2}}}(\partial K)$}H(div, Ω):={q:q∈(L²(Ω))², ∇·q∈L²(Ω)}

对任意v∈W(K)，其弱梯度∇_wv定义为

$ {({\nabla _w}v, \mathit{\boldsymbol{q}})_K} =-{({v_0}, \nabla \cdot \mathit{\boldsymbol{q}})_K} + {\langle {v_b}, \mathit{\boldsymbol{q}} \cdot \mathit{\boldsymbol{n}}\rangle _{\partial K}} $

对任意q∈H(div, Ω) 成立.这里，n是∂K的单位外法线向量，(·, ·)_D表示区域D⊂Ω上的L²内积，〈·, ·〉_∂D表示边界∂D上的L²内积.进一步，定义其离散弱梯度∇_{d, r, K}v为

$ {({\nabla _{d, r, K}}v, \mathit{\boldsymbol{q}})_K} =-{({v_0}, \nabla \cdot \mathit{\boldsymbol{q}})_K} + {\langle {v_b}, \mathit{\boldsymbol{q}} \cdot \mathit{\boldsymbol{n}}\rangle _{\partial K}} $

对任意q∈V(K, r) 成立.这里V(K, r) 表示r次多项式向量空间的子空间，简单起见，后文将省略离散弱梯度的下标r、K.

定义空间

$ \begin{array}{l} {V_h}: = \{ \{ {v_0}, {v_b}\} :{v_0}\left| {_k} \right. \in {P_j}(\mathop K\limits^ \circ ), \\ {v_b}\left| {_{\partial K} \in {P_{{j_1}}}(\partial K), \forall K \in {T_h}\} } \right.\\ {V_{h, 0}}: = \{ v:v \in {V_h}, {v_b} = 0在\partial \mathit{\Omega }上\} \end{array} $

则弱Galerkin有限元方法可写为：找u_h={u₀, u_b}∈V_{h, 0}，使得

$ {a_h}({u_h}, v) = (f, {v_0}), \forall v = \{ {v_0}, {v_b}\} \in {V_{h, 0}} $

(2)

其中a_h(w, v)=(a(x)∇_dw, ∇_dv)_Ω，∀v, w∈V_h.

对任意v∈V_h，引入离散范数^[3]

$ \begin{array}{l} \left\| v \right\|_{0, h, K}^2 = \left\| {{v_0}} \right\|_{0, K}^2 + h\left\| {{v_0}-{v_b}} \right\|_{\partial K}^2, \\ \;\;\;{\left\| v \right\|_{0, h}} = {\left( {\sum\limits_{K \in {T_h}} {\left\| v \right\|_{0, h, K}^2} } \right)^{\frac{1}{2}}} \end{array} $

令Q_hv={Q₀v, Q_bv}表示从H₀¹(Ω) 到V_{h, 0}的L²投影，其中Q₀v|_K是到P_j($\mathop K\limits^ \circ $) 的L²投影，Q_bv|_∂K是到P_j1(∂K) 的L²投影.

令j₁=j，r=j+1，V(K, r) 为Raviart-Thomas元空间^[13-14](简称RT元)，即

$ V(K, r) = {[{P_j}(K)]^2} + {\hat P_j}(K)\boldsymbol{x} $

其中${\hat P_j}(K)$是关于向量x的齐次j次多项式.

根据文献[1]中定理8.3及8.4，有如下估计：

定理1 设u_h是弱Galerkin有限元方法 (2) 的解，u是问题 (1) 的弱解，且满足u∈H^m+1(Ω)，则成立

$ h\left\| {{\nabla _d}{u_h}-\nabla u} \right\| + {\left\| {{u_h}-u} \right\|_{0, h}} \le C{h^{m + 1}}\left\| f \right\| $

2 瀑布型多重网格算法

T_l=T_hl(l≥0) 是一系列嵌套的三角形网格剖分，且网格参数h_l满足h_l=2^-lh₀. V_l=V_{hl, 0}是相应的弱Galerkin有限元空间，V_l^(FE)是相应的连续j次拉格朗日有限元空间，记Q_l=Q_hl，‖·‖_{0, hl}=‖·‖_{0, l}.

第l层上的弱Galerkin有限元方法为：找u_l={u₀, u_b}∈V_l，使得

$ {a_l}({u_l}, {v_l}) = (f, {v_0}){\rm{ }}\forall {v_l} = \{ {v_0}, {v_b}\} \in {V_l} $

(3)

定义A_l算子

$ (({A_l}{u_l}, {v_l})) = {a_l}({u_l}, {v_l}){\rm{ }}\forall {u_l}, {v_l} \in {V_l} $

和能量范数

$ {\left| {\left| {\left| v \right|} \right|} \right|_l} = {(({A_l}v, v))^{1/2}}, \forall v \in {V_l} $

其中，对任意v, w∈V_l，((v, w))=$\sum\limits_{K \in {T_h}} {((} {v_0}, {w_0}{)_K} + h{\langle {v_0}-{v_b}, {w_0}-{w_b}\rangle _{\partial K}})$

基于弱Galerkin有限元的瀑布型多重网格算法的描述如图 1所示.

图 1中，I_l表示网格转移算子，C_l表示迭代算子，m_l是第l层上的迭代次数.

对任意v={v₀, v_b}∈V_l-1，网格转移算子I_l:V_l-1→V_l定义如下：

(1) 对单元内部，即x∈${\mathop K\limits^ \circ _l}$.对任意K_l∈T_l，存在K_l-1∈T_l-1，使得K_l⊂K_l-1，则

$ {({I_l}v)_0}(x) = {v_0}(x) $

(2) 对单元边界，若x∈e_l，e_l落在K_l-1∈T_l-1的内部，则

$ {({I_l}v)_b}(x) = {v_0}(x) $

(3) 对单元边界，若x∈e_l，e_l是边界∂K_l-1的一部分的内部，则

$ {({I_l}v)_b}(x) = {v_b}(x) $

对于上述定义的网格转移算子I_l，有

引理1 对任意v∈V_l-1，成立

$ {\left\| {{I_l}v} \right\|_{0, l}} \le {\left\| v \right\|_{0, l-1}}, \left\| {{\nabla _d}({I_l}v)} \right\| \le C{\rm{|||v||}}{{\rm{|}}_{l-1}} $

证明：由网格转移算子I_l的定义易知第一个不等式成立.

任意K∈T_l-1都由T_l上的若干个单元组成，不妨设K=K₁∪…∪K_m，其中K_i∈T_l.根据文献[3]中式 (3.9)，对任意v∈V_l-1，有

$ {\left\| {{\nabla _d}({I_l}v)\cdot\boldsymbol{n}} \right\|_{\partial {K_i}}} \le Ch_l^{-\frac{1}{2}}{\left\| {{\nabla _d}({I_l}v)} \right\|_{{K_i}}} $

(4)

由文献[3]中的引理3.2的证明，可知

$ {\left\| {\nabla {v_0}} \right\|_K} \le {\left\| {{\nabla _d}v} \right\|_K} $

(5)

$ {\left\| {{v_0}-{v_b}} \right\|_e} \le {h^{\frac{1}{2}}}{\left\| {{\nabla _d}v} \right\|_K} $

(6)

根据离散弱梯度的定义及以上三个式子，可知

$ \begin{array}{l} {({\nabla _d}({I_l}v), {\nabla _d}({I_l}v))_{{K_i}}} = {\rm{ }}{(\nabla {({I_l}v)_0}, {\nabla _d}({I_l}v))_{{K_i}}}-\\ {\langle {({I_l}v)_0}-{({I_l}v)_b}, {\nabla _d}({I_l}v)\cdot\mathit{\boldsymbol{n}}\rangle _{\partial {K_i}}} \le \\ {\left\| {\nabla {{({I_l}v)}_0}} \right\|_{{K_i}}}{\left\| {{\nabla _d}({I_l}v)} \right\|_{{K_i}}} + \\ {\left\| {{{({I_l}v)}_0}-{{({I_l}v)}_b}} \right\|_{\partial {K_i}}}{\left\| {{\nabla _d}({I_l}v)\cdot\mathit{\boldsymbol{n}}} \right\|_{\partial {K_i}}} \le \\ {\left\| {\nabla {v_0}} \right\|_K}{\left\| {{\nabla _d}({I_l}v)} \right\|_{{K_i}}} + \\ {\left\| {{v_0} - {v_b}} \right\|_{\partial K}}{\left\| {{\nabla _d}({I_l}v)\mathit{\boldsymbol{n}}} \right\|_{\partial {K_i}}} \le \\ {\left\| {\nabla {v_0}} \right\|_K}{\left\| {{\nabla _d}({I_l}v)} \right\|_{{K_i}}} + \\ h_{l - 1}^{\frac{1}{2}}{\left\| {{\nabla _d}v} \right\|_K}h_{l - 1}^{\frac{1}{2}}{\left\| {{\nabla _d}({I_l}v)} \right\|_{{K_i}}} \le \\ C{\left\| {{\nabla _d}v} \right\|_K}{\left\| {{\nabla _d}({I_l}v)} \right\|_{{K_i}}} \end{array} $

即第二个不等式成立.

引理2 设u_l为第l层网格上的弱Galerkin有限元问题 (3) 的解，则

$ {\left\| {{u_l}-{I_l}{u_{l-1}}} \right\|_{0, l}} \le Ch_l^2\left\| f \right\| $

证明：假设u_l^(FE)和u_l-1^(FE)分别是T_l和T_l-1剖分上的有限元解.由三角不等式和引理1可得

$ \begin{array}{l} {\left\| {{u_l}-{I_l}{u_{l-1}}} \right\|_{0, l}} \le {\left\| {{u_l}-u} \right\|_{0, l}} + {\left\| {u - u_{l - 1}^{({\rm{FE}})}} \right\|_{0, l}} + \\ {\left\| {{I_l}u_{l - 1}^{({\rm{FE}})} - {I_l}{u_{l - 1}}} \right\|_{0, l}} \le {\left\| {{u_l} - u} \right\|_{0, l}} + \\ \left\| {u - u_{l - 1}^{({\rm{FE}})}} \right\| + {\rm{ }}{\left\| {u_{l - 1}^{({\rm{FE}})} - {u_{l - 1}}} \right\|_{0, l - 1}} \le \\ {C_1}h_l^2\left\| f \right\| + {C_2}h_l^2{\left\| u \right\|_2} + {C_3}h_{l - 1}^2{\left\| u \right\|_2} \le \\ Ch_l^2\left\| f \right\| \end{array} $

引理得证.

设C_l:V_l→V_l是第l层上Richardson、Jacobi、Gauss-Seidel或Conjugate Gradient (CG) 迭代方法对应的迭代算子.利用类似文献[4, 6, 14]中的方法，可以证明存在线性算子R_l:V_l→V_l，使得

$ {u_l}-C_l^{{m_l}}u_l^{(0)} = R_l^{{m_l}}({u_l}-u_l^{(0)}) $

且有

$ \begin{array}{l} |||R_l^{{m_l}}v||{|_l} \le C\frac{{h_l^{-1}}}{{m_l^\gamma }}{\left\| v \right\|_{0, l}}{\rm{, }}\forall v \in {V_l}\\ |||R_l^{{m_l}}v||{|_l} \le {\left| {\left| {\left| v \right|} \right|} \right|_l}{\rm{, }}\forall v \in {V_l} \end{array} $

对于Richardson、Jacobi、Gauss-Seidel迭代方法，γ=1/2.对于CG迭代方法，γ=1.

定义投影算子P_l-1:V_l-1→V_l-1^(FE)为

$ {a_l}({P_{l-1}}u, v) = {a_l}({I_l}u, v), \forall v \in V_{l-1}^{({\rm{FE}})} $

由投影算子的性质，易知

$ {a_l}({P_{l-1}}u, v) = {a_{l-1}}({P_{l-1}}u, v), \forall v \in V_{l - 1}^{{\rm{(FE)}}} $

(7)

$ {a_l}({I_l}u, v) = {a_{l-1}}\left( {u, v} \right), \forall v \in V_{l-1}^{{\rm{(FE)}}} $

(8)

引理3 对于上述定义的投影算子P_l-1，成立

$ |||{P_{l-1}}v||{|_l} \le {\left| {\left| {\left| v \right|} \right|} \right|_{l-1}}, \forall v \in {V_{l-1}} $

(9)

$ {\left\| {{I_l}v-{P_{l-1}}v} \right\|_{0, l}} \le C{h_l}{\left| {\left| {\left| v \right|} \right|} \right|_{l-1}}, \forall v \in {V_{l - 1}} $

(10)

证明：由投影算子P_l的定义及式 (7) 和 (8) 可知，对任意v∈V_l-1，有

$ {a_{l-1}}({P_{l-1}}v, {P_{l-1}}v) = {a_{l - 1}}(v, {P_{l - 1}}v) $

于是

$ \begin{array}{l} {a_{l-1}}\left( {v, v} \right)-{a_l}({P_{l-1}}v, {P_{l - 1}}v) = \\ {a_{l - 1}}(v - {P_{l - 1}}v, v - {P_{l - 1}}v) \ge 0 \end{array} $

这表明|||v|||_l-1≥|||P_l-1v|||_l，即式 (9) 得证.

下面证明式 (10).对任意给定函数v∈V_l-1，引入如下辅助问题：

$ \left\{ \begin{array}{l} -\nabla \cdot (a(x)\nabla \xi ) = {({I_l}v-{P_{l-1}}v)_0}, \;\;\;\;在\mathit{\Omega 内}\\ \xi = 0, \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;在\partial \mathit{\Omega 上} \end{array} \right. $

由于Ω是有界凸多边形区域，所以该问题存在唯一解ξ∈H²(Ω)∩H₀¹(Ω)，且满足正则性估计

$ {\left\| \xi \right\|_2} \le \left\| {{{({I_l}v-{P_{l-1}}v)}_0}} \right\| $

(11)

对任意q∈H(div, Ω)，任意K∈T_l，令Π_lq∈(P_r(K))²满足

$ {(\nabla \cdot\mathit{\boldsymbol{q}}, {v_0})_K} = {(\nabla \cdot{\mathit{\Pi }_l}\mathit{\boldsymbol{q}}, {v_0})_K}, {\rm{ }}\forall {v_0} \in {P_j}(\mathop K\limits^ \circ ) $

由文献[1]中引理7.2和7.3可知

$ \begin{array}{l} {\left\| {{{({I_l}v-{P_{l-1}}v)}_0}} \right\|^2} = \\ -(\nabla \cdot(a\left( x \right)\nabla \xi ), {({I_l}v - {\rm{ }}{P_{l - 1}}v)_0}) = \\ - ({\mathit{\Pi }_l}(a\left( x \right)\nabla \xi ), {\nabla _d}({I_l}v - {P_{l - 1}}v)){\rm{ }}\\ \left\| {{\mathit{\Pi }_l}(a\left( x \right)\nabla \xi ) - {\rm{ }}a\left( x \right){\nabla _d}({Q_l}\xi )} \right\| \le Ch{\left\| \xi \right\|_2}\\ \left\| {\nabla \xi - {\nabla _d}({Q_l}\xi )} \right\| \le Ch{\left\| \xi \right\|_2} \end{array} $

于是

$ \begin{array}{l} {\left\| {{{({I_l}v-{P_{l-1}}v)}_0}} \right\|^2} = \\ ({\mathit{\Pi }_l}(a\nabla \xi )-a{\nabla _d}({Q_l}\xi ), {\nabla _w}({I_l}v - {P_{l - 1}}v)) + \\ (a{\nabla _d}({Q_l}\xi ), {\nabla _d}({I_l}v - {P_{l - 1}}v)) \le \\ \left\| {{\mathit{\Pi }_l}(a\nabla \xi ) - a{\nabla _d}({Q_l}\xi )} \right\|\left\| {{\nabla _d}({I_l}v - {P_{l - 1}}v)} \right\| + \\ {a_l}({Q_l}\xi, {I_l}v - {P_{l - 1}}v) \le \\ C{h_l}{\left\| \xi \right\|_2}\left\| {{\nabla _d}({I_l}v - {P_{l - 1}}v)} \right\| + {a_l}({Q_l}\xi, {I_l}v - {P_{l - 1}}v) \end{array} $

对任意ξ_l-1^(FE)∈V_l-1^(FE)，由

$ {a_l}({I_l}v-{P_{l-1}}v, \xi _{l-1}^{{\rm{(FE)}}}) = 0 $

可得

$ \begin{array}{l} {a_l}({Q_l}\xi, {I_l}v-{P_{l-1}}v) = {a_l}({Q_l}\xi-\xi _{l - 1}^{{\rm{(FE)}}}, {I_l}v - {P_{l - 1}}v) \le \\ C\left\| {{\nabla _d}({Q_l}\xi ) - \Delta \xi _{l - 1}^{{\rm{(FE)}}}} \right\|\left\| {{\nabla _d}({I_l}v - {P_{l - 1}}v)} \right\| \le \\ C{h_l}{\left\| \xi \right\|_2}\left\| {{\nabla _d}({I_l}v - {P_{l - 1}}v)} \right\| \end{array} $

再由式 (11)，得

$ \left\| {{{({I_l}v-{P_{l-1}}v)}_0}} \right\| \le \left\| {{\nabla _d}({I_l}v-{P_{l - 1}}v)} \right\| $

由文献[3]中的引理3.2可知

$ {\left\| {{v_0}-{v_b}} \right\|_e} \le {h^{\frac{1}{2}}}{\left\| {{\nabla _d}v} \right\|_K} $

因此，可得

$ \begin{array}{l} \left\| {{I_l}v-{P_{l-1}}v} \right\|_{0, l}^2 = {\left\| {{{({I_l}v-{P_{l - 1}}v)}_0}} \right\|^2} + \\ \sum\limits_{K \in {T_l}} {{h_l}} \left\| {{{({I_l}v - {P_{l - 1}}v)}_0} - {{({I_l}v - {P_{l - 1}}v)}_b}} \right\|_{\partial K}^2 \le \\ {\left\| {{{({I_l}v - {P_{l - 1}}v)}_0}} \right\|^2} + h_l^2{\left\| {{\nabla _d}({I_l}v - {P_{l - 1}}v)} \right\|^2} \le \\ Ch_l^2{\left\| {{\nabla _d}({I_l}v - {P_{l - 1}}v)} \right\|^2} \end{array} $

即

$ {\left\| {{I_l}v-{P_{l-1}}v} \right\|_{0, l}} \le C{h_l}(\left\| {{\nabla _d}({I_l}v)} \right\| + \left\| {{\nabla _d}({P_{l-1}}v)} \right\|) $

(12)

综合式 (12)、式 (9) 和引理1，即可得到式 (10).

设m_l(0≤l≤L) 是满足m_l≥β^L-1m_L的最小整数，其中β>1，m_L是第L层上的迭代次数.根据文献[7]的理论框架，可得到如下结论：

定理2 对于CG迭代方法，若2≤β≤4，则瀑布型多重网格方法最优，即

$ |||{u_L}-u_L^*||{|_L} \approx |||u-{u_L}||{|_L} $

且计算复杂度为O(n_L).

对于Richardson、Jacobi、Gauss-Seidel迭代方法，若β=4，第l层上的迭代次数满足m_l≥m_*l²(给定m_*≥1)，则有

$ |||{u_L}-u_L^*||{|_L} \le C\frac{{{h_l}}}{{m_*^{\frac{1}{2}}}}\left\| f \right\| $

且计算复杂度为

$ \sum\limits_{k = 1}^L {{m_k}{n_k}} \le C{m_L}{n_L}{(1 + {\rm{log}}{n_L})^3} $

这里n_L表示L层未知量个数.

3 数值实验

本节考虑模型问题 (1)，选取如下算例：

算例1 a(x₁, x₂)=e^x₁+x₂，选取右端函数使得真解u=sin (πx₁) sin (πx₂).

算例2 a(x₁, x₂)=$\frac{1}{{1 + x_1^2}}$，选取右端函数使得真解u=x₁(1-x₁)x₂(1-x₂) e^x₁+x₂.

本文的实验环境为Intel (R) Core (TM) i5-2500 CPU 3.30GH, 12.0 GB内存，64位Windows 7操作系统，实现算法的软件为Matlab R2015b.采用CG迭代算子，最细网格迭代次数m_L=35，参数β=3.

表 1和2分别是算例1和2的计算结果.采用CG迭代，得到基于弱Galerkin有限元离散的瀑布型多重网格算法的计算结果，其中L=1表示仅用CG迭代求解的结果.结果表明：① 由瀑布型多重网格得到的数值解，其能量范数误差与弱Galerkin有限元方法是同阶的，均为O(h)；② 该多重网格方法在计算复杂度方面是最优的，即求解时间与最细层网格上未知量个数n_L成正比.这与前面的理论分析一致.

4 结论与展望

本文应用瀑布型多重网格方法求解弱Galerkin有限元离散二阶椭圆型偏微分方程所得到的代数方程组，证明了当采用共轭梯度法作为光滑迭代子时，本文方法所得数值解的能量误差与弱Galerkin有限元解同阶，且该方法具有最优的计算复杂度.数值算例验证了理论分析结果.本文所提出和分析的网格转移算子，可进一步用于构造其他类型偏微分方程的弱Galerkin有限元瀑布型多重网格算法.