第3回：Gmsh の geometry order は有限要素次数ではない

Introduction

前回は，曲線要素写像を \[ F_K=\Psi_K\circ\Phi_{T_K} \] と分けて考える理由を説明した． ここで，\(\Phi_{T_K}\) は参照三角形から affine core への写像であり，\(\Psi_K\) は affine core から曲線要素への写像である． この分解により，三角形のスケールや異方性と，曲線境界による幾何補正を分けて考えることができる．

今回は，実装上で非常に混同しやすい点を整理する． それは，

\[ \text{geometry order} \neq \text{finite element degree} \]

ということである． Gmsh で曲線メッシュを作るとき，geometry_order という量が出てくる． 一方，FEniCSx で有限要素空間を作るときには，Lagrange 要素の次数 \(k\) を指定する． この二つは似て見えるが，意味はまったく違う．

二つの次数

今回区別したい次数は，次の二つである．

  • 幾何次数：\(q_{\rm geo}\)
  • 有限要素次数：\(k\)

幾何次数 \(q_{\rm geo}\) は，メッシュの形をどの次数で表すかを決める． たとえば，円周を直線で表すのか，二次曲線的に表すのか，三次曲線的に表すのか，という違いである．

一方，有限要素次数 \(k\) は，未知関数 \(u_h\) をどの次数の多項式で近似するかを決める． たとえば，\(k=1\) なら一次 Lagrange 要素，\(k=2\) なら二次 Lagrange 要素である．

したがって，幾何次数は「領域や要素の形」を決める量であり，有限要素次数は「未知関数の近似空間」を決める量である． ここを混同すると，数値結果の解釈を誤りやすい．

幾何次数 \(q_{\rm geo}\) の意味

Gmsh でメッシュを作るとき，幾何次数を変えると，要素の coordinate map が変わる． 簡単に言えば，

  • \(q_{\rm geo}=1\)：直線辺をもつ三角形
  • \(q_{\rm geo}=2\)：二次の幾何表現をもつ曲線三角形
  • \(q_{\rm geo}=3\)：三次の幾何表現をもつ曲線三角形

である．

\(q_{\rm geo}=1\) の場合，円周は直線の chord で近似される． したがって，単位円をメッシュ化しても，実際の計算領域は多角形になる．

一方，\(q_{\rm geo}=2,3\) の場合，境界辺の途中にも幾何ノードが入り，境界が曲線として表現される． そのため，単位円の境界をより正確に表すことができる．

ここで変わっているのは，あくまで要素写像 \[ F_K:\widehat T\to K \] の幾何表現である． 未知関数を近似する多項式次数が変わっているわけではない．

有限要素次数 \(k\) の意味

有限要素次数 \(k\) は，局所有限要素空間を決める． exact-curved Lagrange 要素では，曲線要素 \(K\) 上の局所空間を

\[ V_L^k(K) = \left\{ v:K\to\mathbb R : v=\widehat v\circ F_K^{-1}, \quad \widehat v\in P^k(\widehat T) \right\} \]

と考える．

ここで，\(P^k(\widehat T)\) は参照三角形上の \(k\) 次以下の多項式空間である． したがって，\(k=1\) なら参照三角形上で一次多項式を使う． それを \(F_K^{-1}\) によって物理要素 \(K\) 上へ押し出すことで，曲線要素上の有限要素関数を得る．

重要なのは，\(F_K\) が曲線写像であれば，\(v=\widehat v\circ F_K^{-1}\) は物理座標 \(x\) の多項式とは限らない，という点である． それでも，有限要素次数 \(k\) は，参照要素上の多項式次数として定義される．

今回の数値実験では何を固定し，何を変えるのか

今回の数値実験では，有限要素次数は常に

\[ k=1 \]

に固定する． つまり，未知関数 \(u_h\) は一次 Lagrange 要素で近似する．

その上で，幾何次数だけを

\[ q_{\rm geo}=1,2,3 \]

と変える．

したがって，この実験で比較しているのは，

\[ \text{同じ }P^1\text{ 有限要素空間を使ったときに，幾何表現を変えると何が変わるか} \]

という点である． これは，高次有限要素法の実験ではない． あくまで，未知関数の近似次数を \(P^1\) に固定したまま，メッシュの曲線幾何表現だけを変える実験である．

なぜ \(H^1\) 誤差はあまり変わらないのか

数値結果を見ると，\(q_{\rm geo}=1,2,3\) で幾何誤差は大きく変わる． たとえば，単位円の面積誤差や境界半径誤差は，geometry order を上げるとかなり小さくなる．

一方で，有限要素誤差，特に \(H^1\) 誤差は，geometry order を上げても劇的には変わらない． これは自然である． なぜなら，未知関数の近似空間はすべての場合で \(P^1\) に固定されているからである．

\(P^1\) 有限要素法では，十分滑らかな解に対して，標準的には

\[ |u-u_h|_{H^1}=O(h), \qquad \|u-u_h\|_{L^2}=O(h^2) \]

が期待される． 幾何表現を高次にしても，有限要素空間の次数が \(P^1\) のままであれば，主な \(H^1\) 誤差は \(P^1\) 近似によって決まる．

したがって，geometry order を上げることの効果は，主に幾何誤差の改善として現れる． 有限要素近似そのものの次数を上げたいなら，finite element degree \(k\) も上げる必要がある．

Gmsh では何をしているのか

Gmsh では，メッシュを生成した後に geometry order を指定することで，高次幾何のメッシュを作ることができる． 概念的には，次のような指定である． 

gmsh.model.mesh.generate(2)
gmsh.model.mesh.setOrder(geometry_order)

ここで geometry_order が \(q_{\rm geo}\) に対応する．

この指定により，Gmsh は直線三角形の頂点だけでなく，辺上や要素内部に幾何ノードを配置し，曲線的な coordinate map を構成する． この coordinate map が，FEniCSx に読み込まれた後，数値積分や基底関数の評価で使われる．

ただし，この操作は有限要素空間の次数を変えるものではない． FEniCSx 側で 

V = fem.functionspace(mesh, ("Lagrange", 1))

とすれば，未知関数の有限要素次数は \(k=1\) である． メッシュの geometry order が 2 や 3 であっても，この空間は一次 Lagrange 空間である．

組合せとして見る

幾何次数と有限要素次数は，組合せとして考えると分かりやすい． たとえば，次のような組合せがある．

幾何次数 有限要素次数 意味
\(q_{\rm geo}=1\) \(k=1\) 直線メッシュ上の一次 Lagrange FEM
\(q_{\rm geo}=2\) \(k=1\) 曲線幾何上の一次 Lagrange FEM
\(q_{\rm geo}=3\) \(k=1\) より高次の曲線幾何上の一次 Lagrange FEM
\(q_{\rm geo}=2\) \(k=2\) 二次曲線幾何上の二次 Lagrange FEM
\(q_{\rm geo}=3\) \(k=2\) 三次曲線幾何上の二次 Lagrange FEM

今回の実験は，このうち上の三つ，

\[ (q_{\rm geo},k)=(1,1),\quad (2,1),\quad (3,1) \]

を比較している． つまり，幾何だけを変えて，有限要素空間は変えない実験である． この設定により，曲線幾何が幾何誤差に与える影響を見やすくしている．

まとめ

今回は，Gmsh の geometry order と有限要素次数の違いを整理した． geometry order \(q_{\rm geo}\) は，要素や境界をどの次数の幾何写像で表すかを決める． 一方，finite element degree \(k\) は，未知関数をどの次数の多項式で近似するかを決める．

今回の数値実験では，有限要素次数は常に \(k=1\) に固定し，幾何次数だけを \(q_{\rm geo}=1,2,3\) と変える． したがって，geometry order を上げると幾何誤差は小さくなるが，\(P^1\) 近似の主要な \(H^1\) 誤差次数が自動的に高次になるわけではない．

次回は，FEniCSx/UFL において，弱形式を物理座標で書きながら，実際には参照要素上で積分が行われる仕組みを見る． 特に，grad と dx の中で \(DF_K^{-\top}\) と \(|\det DF_K|\) がどのように現れるかを整理する．

参考文献

この第3回では，曲線幾何の次数と有限要素次数の違いを実装の観点から説明した． 背景として，曲線有限要素法，isoparametric finite elements，および筆者の投稿中論文を挙げておく．

  1. C. Bernardi, Optimal finite-element interpolation on curved domains, SIAM Journal on Numerical Analysis, 26(5), 1212--1240, 1989.
  2. P. G. Ciarlet and P.-A. Raviart, Interpolation theory over curved elements, with applications to finite element methods, Computer Methods in Applied Mechanics and Engineering, 1, 217--249, 1972.
  3. M. Lenoir, Optimal isoparametric finite elements and error estimates for domains involving curved boundaries, SIAM Journal on Numerical Analysis, 23(3), 562--580, 1986.
  4. H. Ishizaka, Exact-curved Lagrange finite elements for the Poisson problem in two dimensions, submitted.

ページの先頭へ