不連続ガラーキン法

応用数学において、不連続ガラーキン法(DG法)は微分方程式を解く数値解析法の一種である。有限要素法有限体積法の枠組みの特徴を組み合わせたもので、幅広い応用分野から生じる双曲型、楕円放物型、混合形式の問題に効果的に適用されてきた。DG法は特に、電気力学流体力学プラズマ物理学など、支配的な一次部分を持つ問題で大きな関心を集めている。実際、こうした問題の解には強い勾配(さらには不連続性)が含まれる場合があり、古典的な有限要素法ではうまくいかず、有限体積法では低次の近似しか行えない。

不連続ガラーキン法は、偏微分方程式を数値的に解く手法として、1970年代初頭に初めて提案され、解析されました。1973年、リードとヒルは双曲型中性子輸送方程式を解くためのDG法を導入しました。

楕円型問題に対するDG法の起源は、ジャンプペナルティなどの現代的な機能が徐々に開発されたため、単一の出版物にまで遡ることはできません。しかし、初期の影響力のある貢献者としては、BabuškaJ.-L. Lions、Joachim Nitsche、Miloš Zlámalなどが挙げられます。楕円型問題に対するDG法は、1977年にGarth Bakerが4次方程式の設定で論文を発表し、すでに開発されていました。楕円型問題に対するDG法の歴史的発展と入門については、Arnold、Brezzi、Cockburn、Mariniによる出版物でより詳しく説明されています。DG法に関する研究の方向性と課題の多くは、Cockburn、Karniadakis、Shuが編集した議事録集にまとめられています。

概要

連続ガラーキン法(CG法)と同様に、不連続ガラーキン法(DG法)は、特定のモデル系の弱定式化に基づいて定式化された有限要素法です。従来のCG法が適合するのに対し、DG法は、区分的に連続する関数の試行空間上でのみ作用するため、適合法で用いられる有限次元の内積部分空間よりも 多くの包括関数空間を構成することがよくあります。

例として、「ソース」や「シンク」のない 空間領域におけるスカラーの未知の連続方程式を考えてみましょう。ρ{\displaystyle \rho }Ω{\displaystyle \オメガ}

ρt+J0{\displaystyle {\frac {\partial \rho }{\partial t}}+\nabla \cdot \mathbf {J} =0,}

のフラックスはどこですか? J{\displaystyle \mathbf {J} }ρ{\displaystyle \rho }

ここで、離散三角形分割に制限された空間領域上の不連続な区分多項式関数の有限次元空間を考える。これは次のように表される。 Ω{\displaystyle \オメガ}Ωh{\displaystyle \Omega _{h}}

ShpΩh{v|ΩePpΩe  ΩeΩh}{\displaystyle S_{h}^{p}(\Omega _{h})=\{v_{|\Omega _{e_{i}}}\in P^{p}(\Omega _{e_{i}}),\ \ \forall \Omega _{e_{i}}\in \Omega _{h}\}}

で添え字付けされた要素上の次数以下の多項式空間に対して、有限要素形状関数に対して解は次のように表される 。PpΩe{\displaystyle P^{p}(\オメガ _{e_{i}})}p{\displaystyle p}Ωe{\displaystyle \オメガ _{e_{i}}}{\displaystyle i}jPp{\displaystyle N_{j}\in P^{p}}

ρhj1被写界深度ρjtj××Ωe{\displaystyle \rho _{h}^{i}=\sum _{j=1}^{\text{dofs}}\rho _{j}^{i}(t)N_{j}^{i}({\boldsymbol {x}}),\quad \forall {\boldsymbol {x}}\in \Omega _{e_{i}}.}

同様にテスト関数を選択する

φhi(x)=j=1dofsφjiNji(x),xΩei,{\displaystyle \varphi _{h}^{i}({\boldsymbol {x}})=\sum _{j=1}^{\text{dofs}}\varphi _{j}^{i}N_{j}^{i}({\boldsymbol {x}}),\quad \forall {\boldsymbol {x}}\in \Omega _{e_{i}},}

連続方程式に を掛けて空間の部分で積分すると、離散DG定式は次のようになります。 φhi{\displaystyle \varphi _{h}^{i}}

ddtΩeiρhiφhidx+ΩeiφhiJhndx=ΩeiJhφhidx.{\displaystyle {\frac {d}{dt}}\int _{\Omega _{e_{i}}}\rho _{h}^{i}\varphi _{h}^{i}\,d{\boldsymbol {x}}+\int _{\partial \Omega _{e_{i}}}\varphi _{h}^{i}\mathbf {J} _{h}\cdot {\boldsymbol {n}}\,d{\boldsymbol {x}}=\int _{\Omega _{e_{i}}}\mathbf {J} _{h}\cdot \nabla \varphi _{h}^{i}\,d{\boldsymbol {x}}.}

スカラー双曲型保存則

スカラー双曲型保存則は次の形式をとる。

tu+xf(u)=0fort>0,xRu(0,x)=u0(x),{\displaystyle {\begin{aligned}\partial _{t}u+\partial _{x}f(u)&=0\quad {\text{for}}\quad t>0,\,x\in \mathbb {R} \\u(0,x)&=u_{0}(x)\,,\end{aligned}}}

ここで、未知のスカラー関数 を解こうとしますが、関数は通常与えられています。 uu(t,x){\displaystyle u\equiv u(t,x)}f,u0{\displaystyle f,u_{0}}

空間離散化

空間は次のように離散化される。 x{\displaystyle x}

R=kIk,Ik:=(xk,xk+1)forxk<xk+1.{\displaystyle \mathbb {R} =\bigcup _{k}I_{k}\,,\quad I_{k}:=\left(x_{k},x_{k+1}\right)\quad {\text{for}}\quad x_{k}<x_{k+1}\,.}

さらに、以下の定義が必要です

hk:=|Ik|,h:=supkhk,x^k:=xk+hk2.{\displaystyle h_{k}:=|I_{k}|\,,\quad h:=\sup _{k}h_{k}\,,\quad {\hat {x}}_{k}:=x_{k}+{\frac {h_{k}}{2}}\,.}

関数空間の基礎

解の関数空間の基底表現を導出する。関数空間は次のように定義される。 u{\displaystyle u}

Shp:={vL2(R):v|IkΠp}forpN0,{\displaystyle S_{h}^{p}:=\left\lbrace v\in L^{2}(\mathbb {R} ):v{\Big |}_{I_{k}}\in \Pi _{p}\right\rbrace \quad {\text{for}}\quad p\in \mathbb {N} _{0}\,,}

ここで、 は区間 への の制限を表し、 は最大次数の多項式空間を表します。添え字は、によって与えられる基礎離散化との関係を示す必要があります。ここで、は交差点 において一意に定義されないことに注意してください。 v|Ik{\displaystyle {v|}_{I_{k}}}v{\displaystyle v}Ik{\displaystyle I_{k}}Πp{\displaystyle \Pi _{p}}p{\displaystyle p}h{\displaystyle h}(xk)k{\displaystyle \left(x_{k}\right)_{k}}v{\displaystyle v}(xk)k{\displaystyle (x_{k})_{k}}

まず、区間 上の特定の多項式基底、すなわちルジャンドル多項式、すなわち [1,1]{\displaystyle [-1,1]}(Pn)nN0{\displaystyle (P_{n})_{n\in \mathbb {N} _{0}}}

P0(x)=1,P1(x)=x,P2(x)=12(3x21),{\displaystyle P_{0}(x)=1\,,\quad P_{1}(x)=x\,,\quad P_{2}(x)={\frac {1}{2}}(3x^{2}-1)\,,\quad \dots }

特に直交関係に注目

Pi,PjL2([1,1])=22i+1δiji,jN0.{\displaystyle \left\langle P_{i},P_{j}\right\rangle _{L^{2}([-1,1])}={\frac {2}{2i+1}}\delta _{ij}\quad \forall \,i,j\in \mathbb {N} _{0}\,.}

区間への変換と正規化は関数によって実現される。[0,1]{\displaystyle [0,1]}(φi)i{\displaystyle (\varphi _{i})_{i}}

φi(x):=2i+1Pi(2x1)forx[0,1],{\displaystyle \varphi _{i}(x):={\sqrt {2i+1}}P_{i}(2x-1)\quad {\text{for}}\quad x\in [0,1]\,,}

直交関係を満たす

φi,φjL2([0,1])=δiji,jN0.{\displaystyle \left\langle \varphi _{i},\varphi _{j}\right\rangle _{L^{2}([0,1])}=\delta _{ij}\quad \forall \,i,j\in \mathbb {N} _{0}\,.}

区間への変換は次のように与えられる。Ik{\displaystyle I_{k}}(φ¯ki)i{\displaystyle \left({\bar {\varphi }}_{ki}\right)_{i}}

φ¯ki:=1hkφi(xxkhk)forxIk,{\displaystyle {\bar {\varphi }}_{ki}:={\frac {1}{\sqrt {h_{k}}}}\varphi _{i}\left({\frac {x-x_{k}}{h_{k}}}\right)\quad {\text{for}}\quad x\in I_{k}\,,}

満たす

φ¯ki,φ¯kjL2(Ik)=δiji,jN0k.{\displaystyle \left\langle {\bar {\varphi }}_{ki},{\bar {\varphi }}_{kj}\right\rangle _{L^{2}(I_{k})}=\delta _{ij}\quad \forall \,i,j\in \mathbb {N} _{0}\forall \,k\,.}

-正規化の場合、を定義し、 -正規化の場合、 を定義する。 L{\displaystyle L^{\infty }}φki:=hkφ¯ki{\displaystyle \varphi _{ki}:={\sqrt {h_{k}}}{\bar {\varphi }}_{ki}}L1{\displaystyle L^{1}}φ~ki:=1hkφ¯ki{\displaystyle {\tilde {\varphi }}_{ki}:={\frac {1}{\sqrt {h_{k}}}}{\bar {\varphi }}_{ki}}

φkiL(Ik)=φiL([0,1])=:ci,andφ~kiL1(Ik)=φiL1([0,1])=:ci,1.{\displaystyle \|\varphi _{ki}\|_{L^{\infty }(I_{k})}=\|\varphi _{i}\|_{L^{\infty }([0,1])}=:c_{i,\infty }\quad {\text{and}}\quad \|{\tilde {\varphi }}_{ki}\|_{L^{1}(I_{k})}=\|\varphi _{i}\|_{L^{1}([0,1])}=:c_{i,1}\,.}

最後に、解の基底表現を定義する。uh{\displaystyle u_{h}}

uh(t,x):=i=0puki(t)φki(x)forx(xk,xk+1)uki(t)=uh(t,),φ~kiL2(Ik).{\displaystyle {\begin{aligned}u_{h}(t,x):=&\sum _{i=0}^{p}u_{ki}(t)\varphi _{ki}(x)\quad {\text{for}}\quad x\in (x_{k},x_{k+1})\\u_{ki}(t)=&\left\langle u_{h}(t,\cdot ),{\tilde {\varphi }}_{ki}\right\rangle _{L^{2}(I_{k})}\,.\end{aligned}}}

ここで、インターフェースの位置では定義されていないことに注意してください。 uh{\displaystyle u_{h}}

また、プリズムベースは平面状の構造に使用され、2D/3Dのハイブリッド化が可能です。

DGスキーム

保存則は、テスト関数を乗じてテスト区間にわたって積分することで弱形式に変換される。

tu+xf(u)=0tu,vL2(Ik)+xf(u),vL2(Ik)=0forvShptu,φ~kiL2(Ik)+xf(u),φ~kiL2(Ik)=0forkip.{\displaystyle {\begin{aligned}\partial _{t}u+\partial _{x}f(u)&=0\\\Rightarrow \quad \left\langle \partial _{t}u,v\right\rangle _{L^{2}(I_{k})}+\left\langle \partial _{x}f(u),v\right\rangle _{L^{2}(I_{k})}&=0\quad {\text{for}}\quad \forall \,v\in S_{h}^{p}\\\Leftrightarrow \quad \left\langle \partial _{t}u,{\tilde {\varphi }}_{ki}\right\rangle _{L^{2}(I_{k})}+\left\langle \partial _{x}f(u),{\tilde {\varphi }}_{ki}\right\rangle _{L^{2}(I_{k})}&=0\quad {\text{for}}\quad \forall \,k\;\forall \,i\leq p\,.\end{aligned}}}

部分積分を使うと、

ddtuki(t)+f(u(t,xk+1))φ~ki(xk+1)f(u(t,xk))φ~ki(xk)f(u(t,)),φ~kiL2(Ik)=0forkip.{\displaystyle {\begin{aligned}{\frac {\mathrm {d} }{\mathrm {d} t}}u_{ki}(t)+f(u(t,x_{k+1})){\tilde {\varphi }}_{ki}(x_{k+1})-f(u(t,x_{k})){\tilde {\varphi }}_{ki}(x_{k})-\left\langle f(u(t,\,\cdot \,)),{\tilde {\varphi }}_{ki}'\right\rangle _{L^{2}(I_{k})}=0\quad {\text{for}}\quad \forall \,k\;\forall \,i\leq p\,.\end{aligned}}}

界面におけるフラックスは、数値フラックスによって近似される。 g{\displaystyle g}

gk:=g(uk,uk+),uk±:=u(t,xk±),{\displaystyle g_{k}:=g(u_{k}^{-},u_{k}^{+})\,,\quad u_{k}^{\pm }:=u(t,x_{k}^{\pm })\,,}

ここで、 は左辺と右辺の極限を表す。最終的に、DGスキームは次のように書ける。 uk±{\displaystyle u_{k}^{\pm }}

ddtuki(t)+gk+1φ~ki(xk+1)gkφ~ki(xk)f(u(t,)),φ~kiL2(Ik)=0forkip.{\displaystyle {\begin{aligned}{\frac {\mathrm {d} }{\mathrm {d} t}}u_{ki}(t)+g_{k+1}{\tilde {\varphi }}_{ki}(x_{k+1})-g_{k}{\tilde {\varphi }}_{ki}(x_{k})-\left\langle f(u(t,\,\cdot \,)),{\tilde {\varphi }}_{ki}'\right\rangle _{L^{2}(I_{k})}=0\quad {\text{for}}\quad \forall \,k\;\forall \,i\leq p\,.\end{aligned}}}

スカラー楕円方程式

スカラー楕円方程式は次の形式をとる。

xxu=f(x)forx(a,b)u(x)=g(x)forx=a,b{\displaystyle {\begin{aligned}-\partial _{xx}u&=f(x)\quad {\text{for}}\quad x\in (a,b)\\u(x)&=g(x)\,\quad {\text{for}}\,\quad x=a,b\end{aligned}}}

この方程式は定常熱方程式であり、 は温度です。空間の離散化は上記と同じです。区間は長さ の区間に分割されることを思い出してください。 u{\displaystyle u}(a,b){\displaystyle (a,b)}N+1{\displaystyle N+1}h{\displaystyle h}

ノードで関数のジャンプと平均を導入します。 []{\displaystyle [{}\cdot {}]}{}{\displaystyle \{{}\cdot {}\}}xk{\displaystyle x_{k}}

[v]|xk=v(xk+)v(xk),{v}|xk=0.5(v(xk+)+v(xk)){\displaystyle [v]{\Big |}_{x_{k}}=v(x_{k}^{+})-v(x_{k}^{-}),\quad \{v\}{\Big |}_{x_{k}}=0.5(v(x_{k}^{+})+v(x_{k}^{-}))}

内部ペナルティ不連続ガラーキン法(IPDG法)は、以下の式を満たす uh{\displaystyle u_{h}}

A(uh,vh)+A(uh,vh)=(vh)+(vh){\displaystyle A(u_{h},v_{h})+A_{\partial }(u_{h},v_{h})=\ell (v_{h})+\ell _{\partial }(v_{h})}

ここで、双線型形式 とA{\displaystyle A}A{\displaystyle A_{\partial }}

A(uh,vh)=k=1N+1xk1xkxuhxvhk=1N{xuh}xk[vh]xk+εk=1N{xvh}xk[uh]xk+σhk=1N[uh]xk[vh]xk{\displaystyle A(u_{h},v_{h})=\sum _{k=1}^{N+1}\int _{x_{k-1}}^{x_{k}}\partial _{x}u_{h}\partial _{x}v_{h}-\sum _{k=1}^{N}\{\partial _{x}u_{h}\}_{x_{k}}[v_{h}]_{x_{k}}+\varepsilon \sum _{k=1}^{N}\{\partial _{x}v_{h}\}_{x_{k}}[u_{h}]_{x_{k}}+{\frac {\sigma }{h}}\sum _{k=1}^{N}[u_{h}]_{x_{k}}[v_{h}]_{x_{k}}}

そして

A(uh,vh)=xuh(a)vh(a)xuh(b)vh(b)εxvh(a)uh(a)+εxvh(b)uh(b)+σh(uh(a)vh(a)+uh(b)vh(b)){\displaystyle A_{\partial }(u_{h},v_{h})=\partial _{x}u_{h}(a)v_{h}(a)-\partial _{x}u_{h}(b)v_{h}(b)-\varepsilon \partial _{x}v_{h}(a)u_{h}(a)+\varepsilon \partial _{x}v_{h}(b)u_{h}(b)+{\frac {\sigma }{h}}{\big (}u_{h}(a)v_{h}(a)+u_{h}(b)v_{h}(b){\big )}}

線型形式と線型形式は {\displaystyle \ell }{\displaystyle \ell _{\partial }}

(vh)=abfvh{\displaystyle \ell (v_{h})=\int _{a}^{b}fv_{h}}

そして

(vh)=εxvh(a)g(a)+εxvh(b)g(b)+σh(g(a)vh(a)+g(b)vh(b)){\displaystyle \ell _{\partial }(v_{h})=-\varepsilon \partial _{x}v_{h}(a)g(a)+\varepsilon \partial _{x}v_{h}(b)g(b)+{\frac {\sigma }{h}}{\big (}g(a)v_{h}(a)+g(b)v_{h}(b){\big )}}

ペナルティパラメータは正の定数です。この値を大きくすると、不連続解のジャンプが減少します。この項は、対称内部ペナルティガラーキン法では に等しくなるように選択され、非対称内部ペナルティガラーキン法では に等しくなります。 σ{\displaystyle \sigma }ε{\displaystyle \varepsilon }1{\displaystyle -1}+1{\displaystyle +1}

直接不連続ガラーキン法

直接不連続ガラーキン法(DDG法)は、拡散問題を解くための新しい不連続ガラーキン法です。2009年にLiuとYanは拡散方程式を解くためのDDG法を初めて提案しました。[ 1 ] [ 2 ]この法は、不連続ガラーキン法と比較して、中間変数を導入することなく、関数の数値フラックスと一次微分項を直接取得することで数値形式を導出する点で優れています。この法を用いることで、妥当な数値結果を得ることができ、導出プロセスが簡素化されるため、計算量が大幅に削減されます。

直接不連続有限要素法は、不連続ガラーキン法の一分野です。主に、問題の変分形式への変換、領域単位の分割、基底関数の構築、不連続有限要素方程式の作成と解法、収束および誤差解析などが含まれます。

たとえば、1次元の非線形拡散方程式を考えてみましょう。

Ut(a(U)Ux)x=0  in (0,1)×(0,T){\displaystyle U_{t}-{(a(U)\cdot U_{x})}_{x}=0\ \ in\ (0,1)\times (0,T)}、その中でU(x,0)=U0(x)  on (0,1){\displaystyle U(x,0)=U_{0}(x)\ \ on\ (0,1)}

空間離散化

まず、、 を定義します。したがって、 の空間離散化が完了しました。また、 を定義します。 {Ij=(xj12, xj+12),j=1...N}{\displaystyle \left\{I_{j}=\left(x_{j-{\frac {1}{2}}},\ x_{j+{\frac {1}{2}}}\right),j=1...N\right\}}Δxj=xj+12xj12{\displaystyle \Delta x_{j}=x_{j+{\frac {1}{2}}}-x_{j-{\frac {1}{2}}}}x{\displaystyle x}Δx=max1j<N Δxj{\displaystyle \Delta x=\max _{1\leq j<N}\ \Delta x_{j}}

となるような近似 値を求めたい。u{\displaystyle u}U{\displaystyle U}t[0,T]{\displaystyle \forall t\in \left[0,T\right]}uVΔx{\displaystyle u\in \mathbb {V} _{\Delta x}}

VΔx:={vL2(0,1):v|IjPk(Ij), j=1,...,N}{\displaystyle \mathbb {V} _{\Delta x}:=\left\{v\in L^{2}\left(0,1\right):{v|}_{I_{j}}\in P^{k}\left(I_{j}\right),\ j=1,...,N\right\}}は、次以下のにおける多項式空間です。 Pk(Ij){\displaystyle P^{k}\left(I_{j}\right)}Ij{\displaystyle I_{j}}k{\displaystyle k}

計画の策定

フラックス: . h:=h(U,Ux)=a(U)Ux{\displaystyle h:=h\left(U,U_{x}\right)=a\left(U\right)U_{x}}

U{\displaystyle U}: 方程式の正確な解。

この方程式に滑らかな関数を掛けると、次の方程式が得られます。 vH1(0,1){\displaystyle v\in H^{1}\left(0,1\right)}

IjUtvdxhj+12vj+12+hj12vj12+a(U)Uxvxdx=0{\displaystyle \int _{I_{j}}U_{t}vdx-h_{j+{\frac {1}{2}}}v_{j+{\frac {1}{2}}}+h_{j-{\frac {1}{2}}}v_{j-{\frac {1}{2}}}+\int a\left(U\right)U_{x}v_{x}dx=0}

IjU(x,0)v(x)dx=IjU0(x)v(x)dx{\displaystyle \int _{I_{j}}U\left(x,0\right)v\left(x\right)dx=\int _{I_{j}}U_{0}\left(x\right)v\left(x\right)dx}

ここでは任意であり、方程式の正確な解は近似解に置き換えられます。つまり、必要な数値解は微分方程式を解くことによって得られます。 v{\displaystyle v}U{\displaystyle U}u{\displaystyle u}

数値フラックス

DDG 法の精度を上げるには、適切な数値フラックスを選択することが重要です。

数値フラックスは次の条件を満たす必要があります。

♦ これは、h=b(u)x=a(u)ux{\displaystyle h={b\left(u\right)}_{x}=a\left(u\right)u_{x}}

♦ 数値フラックスは 上の単一の値で保存されます。 xj+12{\displaystyle x_{j+{\frac {1}{2}}}}

♦安定性があります。 L2{\displaystyle L^{2}}

♦ 方法の精度を向上させることができます。

したがって、数値フラックスの一般的なスキームは次のようになります。

h^=Dxb(u)=β0[b(u)]Δx+b(u)x¯+m=1k2βm(Δx)2m1[x2mb(u)]{\displaystyle {\widehat {h}}=D_{x}b(u)=\beta _{0}{\frac {\left[b\left(u\right)\right]}{\Delta x}}+{\overline {{b\left(u\right)}_{x}}}+\sum _{m=1}^{\frac {k}{2}}\beta _{m}{\left(\Delta x\right)}^{2m-1}\left[\partial _{x}^{2m}b\left(u\right)\right]}

このフラックスにおいて、は隣接する2つの計算ユニットにおける多項式の最大次数です。は関数のジャンプです。非一様グリッドでは が、一様グリッドでは となることに注意してください。 k{\displaystyle k}[]{\displaystyle \left[\cdot \right]}Δx{\displaystyle \Delta x}(Δxj+Δxj+12){\displaystyle \left({\frac {\Delta x_{j}+\Delta x_{j+1}}{2}}\right)}1N{\displaystyle {\frac {1}{N}}}

誤差推定

正確な解と数値解の間の誤差がであることを示します。 U{\displaystyle U}u{\displaystyle u}e=uU{\displaystyle e=u-U}

次の基準でエラーを測定します。

|||v(,t)|||=(01v2dx+(1γ)0tj=1NIjvx2dxdτ+α0tj=1N[v]2/Δxdτ)0.5{\displaystyle \left|\left|\left|v(\cdot ,t)\right|\right|\right|={\left(\int _{0}^{1}v^{2}dx+\left(1-\gamma \right)\int _{0}^{t}\sum _{j=1}^{N}\int _{I_{j}}v_{x}^{2}dxd\tau +\alpha \int _{0}^{t}\sum _{j=1}^{N}{\left[v\right]}^{2}/\Delta x\cdot d\tau \right)}^{0.5}}

そして、|||U(,T)||||||U(,0)|||{\displaystyle \left|\left|\left|U(\cdot ,T)\right|\right|\right|\leq \left|\left|\left|U(\cdot ,0)\right|\right|\right|}|||u(,T)||||||U(,0)|||{\displaystyle \left|\left|\left|u(\cdot ,T)\right|\right|\right|\leq \left|\left|\left|U(\cdot ,0)\right|\right|\right|}

参照

参考文献

  1. ^ Hailiang Liu, Jue Yan,拡散問題に対する直接不連続ガラーキン法(DDG法) , SIAM J. NUMER. ANAL. Vol. 47, No. 1, pp. 675–698.
  2. ^ Hailiang Liu, Jue Yan、「界面補正機能を備えた拡散のための直接不連続ガラーキン法(DDG)」、「Commun. Comput. Phys. Vol. 8, No. 3, pp. 541-564」