第6回：Exact-curved FEMs と isoparametric FEMs は何が違うのか
Introduction
このシリーズでは，曲線領域を有限要素法で扱うために，exact-curved FEM という見方を使っている． しかし，曲線メッシュや曲線要素に少し詳しい人なら，おそらく次の疑問を持つと思う．
それは，普通の isoparametric FEM と何が違うのか？
これは非常に自然な疑問である． 実際，Gmsh で曲線メッシュを作り，FEniCSx に読み込んで，参照要素上で数値積分するという実装だけを見ると，isoparametric FEM とかなり似ている．
だから，まず最初に結論を言っておく．
実装の見た目はかなり似ている．しかし，解析で何を分けて見ているかが違う．
本シリーズで強調している exact-curved FEM の特徴は，要素写像を
\[ F_K=\Psi_K\circ\Phi_{T_K} \]と分けて見る点にある．ここで，\(\Phi_{T_K}\) は直線三角形である affine core を作る写像であり，\(\Psi_K\) はその affine core を曲線要素へ移す写像である．
つまり，曲線要素を一つの写像としてまとめて扱うのではなく，
- 要素の大きさ・向き・細長さを表す部分
- 境界の曲がりや曲率を表す部分
を分けて考える．この分け方が，isoparametric FEM との大きな違いである．
まず isoparametric FEM とは何か
isoparametric FEM は，曲線境界をもつ領域を有限要素法で扱うための標準的な方法である．考え方はとても自然である．
まず，参照要素 \(\widehat K\) を用意する． そして，参照要素上の形状関数を使って，物理要素 \(K\) への写像を作る．たとえば，幾何ノードを \(a_i\)，参照要素上の形状関数を \(\widehat\varphi_i\) とすると，典型的には
\[ F_K(\widehat x) = \sum_i a_i\,\widehat\varphi_i(\widehat x) \]のように書く．この写像 \(F_K\) によって，参照要素が物理要素へ移される．
このとき，未知関数も同じ参照要素上の形状関数を使って近似する．したがって，幾何の写像と未知関数の近似が，同じような形状関数で作られる．これが isoparametric という名前の由来である．
たとえば，二次 Lagrange 要素を使うなら，幾何も二次的に曲がり，未知関数も二次の形状関数で近似する，というのが典型的な考え方である．
もちろん，現代のソフトウェアでは，幾何次数と有限要素次数を別々に指定できることも多い．しかし，isoparametric FEM の基本的な発想は，
参照要素から物理要素への写像 \(F_K\) を作り，その写像を使って計算する
というものである．
実装を見ると，たしかに似ている
曲線要素上で Poisson 問題を解くとき，局所剛性行列は，形式的には次のように書ける．
\[ A^K_{ij} = \int_{\widehat T} \left( DF_K(\widehat x)^{-\top} \nabla_{\widehat x}\widehat\varphi_j(\widehat x) \right) \cdot \left( DF_K(\widehat x)^{-\top} \nabla_{\widehat x}\widehat\varphi_i(\widehat x) \right) |\det DF_K(\widehat x)|\,d\widehat x. \]ここには，有限要素法でおなじみの二つの量が出てくる．
- \(DF_K^{-\top}\)：勾配を物理要素上の勾配へ変換するための因子
- \(|\det DF_K|\)：面積要素を変換するための因子
これは isoparametric FEM でも現れる式である．そのため，実装だけを見ると，
これは isoparametric FEM と同じでは？
と思うのは自然である．
実際，本シリーズでも，実装では Gmsh が作る curved coordinate map を使っている．FEniCSx/UFL は，弱形式を物理座標で書いても，内部では参照要素に引き戻して積分してくれる．この意味では，実装の表面は isoparametric FEM とかなり近い．
では，どこが違うのか
違いは，写像 \(F_K\) をどう見るかにある．
isoparametric FEM では，通常，参照要素から物理要素への写像 \(F_K\) を一つの写像として扱う．もちろん，そのヤコビアン \(DF_K\) や行列式 \(\det DF_K\) を評価する．しかし，解析の中心にあるのは，あくまで一つの写像 \(F_K\) である．
一方，本シリーズの exact-curved FEM では，その写像を
\[ F_K=\Psi_K\circ\Phi_{T_K} \]と分けて見る．
この分解は，次の二段階を表している．
\[ \widehat T \overset{\Phi_{T_K}}{\longrightarrow} T_K \overset{\Psi_K}{\longrightarrow} K. \]ここで，
- \(\widehat T\)：参照三角形
- \(T_K\)：曲げる前の直線三角形，つまり affine core
- \(K\)：物理的な曲線要素
である．
つまり，まず参照三角形を直線三角形 \(T_K\) に写し，そのあと \(T_K\) を曲げて曲線要素 \(K\) にする．この「いったん直線三角形を経由する」見方が大切である．
たとえるなら，形を作る工程を分ける
少し直感的に言うと，isoparametric FEM では，最終的な曲線要素を一つの写像で作る．一方，exact-curved FEM の見方では，次の二段階に分ける．
- まず，参照三角形を伸ばしたり回したりして，直線三角形 \(T_K\) を作る．
- 次に，その直線三角形を曲げて，曲線要素 \(K\) を作る．
この二段階に分けることで，
- どれくらい細長い三角形なのか
- どの方向に伸びているのか
- どれくらい境界が曲がっているのか
を別々に見ることができる．特に異方性メッシュを考えるときには，「細長さ」と「曲がり」を混ぜてしまうと，誤差評価が非常に見えにくくなる．そこで，affine core と curved correction を分けるのである．
比較表
ここまでの違いを表にまとめると，次のようになる．
| 観点 | isoparametric FEM | exact-curved FEM in this series |
|---|---|---|
| 基本的な見方 | 参照要素から物理要素への写像 \(F_K\) を作る | \(F_K=\Psi_K\circ\Phi_{T_K}\) と分解して見る |
| 幾何の扱い | 幾何ノードと形状関数で曲線要素を表す | affine core と curved correction に分ける |
| 解析で見たいもの | 写像全体の性質 | affine scaling と curvature effect の分離 |
| 異方性との相性 | 写像全体のヤコビアンに混ざりやすい | 異方性は \(\Phi_{T_K} \)，曲率は \(\Psi_K\) で追跡する |
| 実装上の見え方 | 標準 FEM ソフトウェアで自然に実装される | 実装は isoparametric 的に見えることが多い |
| 目的 | 曲線要素を実装する標準的枠組み | 曲線幾何をもつ有限要素法を解析しやすく整理する枠組み |
この表で重要なのは，exact-curved FEM が isoparametric FEM を否定しているわけではない，という点である．むしろ，実装では isoparametric 的な仕組みを使うことが自然である． 違いは，その写像を解析の中でどう分解して見るかにある．
Bernardi の exact curved element との関係
曲線要素上の補間理論として重要なのが，Bernardi の exact curved element の理論である．そこでは，曲線要素の写像を，アフィン写像と曲線的な摂動の和として考える．形式的には，
\[ G_K=\widetilde G_K+\varphi_K \]のように表される．ここで，\(\widetilde G_K\) はアフィン写像，\(\varphi_K\) は曲線的な摂動である．
この見方は，本シリーズの見方と非常に近い．実際，形式的には
\[ G_K = \left( \mathrm{id} + \varphi_K\circ\widetilde G_K^{-1} \right) \circ \widetilde G_K \]と書ける．つまり，加法的な摂動として書かれた写像も，合成写像として見ることができる．
ただし，本シリーズでは，この合成構造を最初から解析の中心に置く．すなわち，
\[ F_K=\Psi_K\circ\Phi_{T_K} \]として，affine scaling と curved correction を明示的に分ける．
異方性メッシュでなぜ効くのか
この分解が特に効いてくるのは，異方性メッシュを考えるときである．異方性メッシュでは，要素が細長くなることがある．そのため，単に要素サイズ \(h\) だけで誤差を測ると，方向ごとの情報が消えてしまう．
そこで，affine core \(T_K\) 上で方向を決める．たとえば，\(T_K\) 上の方向を \(r_i\) とする．この方向は，曲線写像 \(\Psi_K\) によって，物理要素 \(K\) 上の方向へ運ばれる．
\[ \tau_{K,i}(x) = D\Psi_K(\Psi_K^{-1}(x))r_i. \]この \(\tau_{K,i}\) が transported direction である．
こうすると，
- 方向ごとの長さや異方性は affine core \(T_K\) 上で見る
- その方向が曲線要素上でどう変わるかは \(\Psi_K\) で見る
という整理ができる．この分離により，異方性の影響と曲率の影響を，混ぜずに追跡できる． ここが，今回の exact-curved decomposition の一番大きな利点である．
では exact-curved FEM は新しい実装方法なのか
ここも誤解しやすい点である．本シリーズでいう exact-curved FEM は，「まったく新しい実装ソフトウェアの作り方」を意味しているわけではない．
実装では，Gmsh や FEniCSx のような既存の有限要素ソフトウェアを使ってよい．実際，曲線メッシュを生成し，それを FEniCSx に読み込むと，ソフトウェア側は coordinate map を使って数値積分を処理してくれる．
本シリーズの主張は，実装そのものよりも，解析の見方にある．すなわち，ソフトウェアが使っている coordinate map を，理論的には
\[ F_K=\Psi_K\circ\Phi_{T_K} \]という分解をもつものとして理解し，その上で補間誤差や有限要素誤差を評価する，という立場である．
この意味で，exact-curved FEM は，isoparametric FEM と競合するものではない．むしろ，既存の曲線要素実装を，異方性や曲率の観点から見直すための理論的な整理である．
まとめ
今回は，exact-curved FEMs と isoparametric FEMs の違いを整理した． 実装だけを見ると，両者はかなり似ている．どちらも参照要素から物理要素への写像を使い，ヤコビアンを通じて勾配や積分測度を変換する．
しかし，解析の立場では違いがある．isoparametric FEM では，要素写像 \(F_K\) を一つの写像として扱うことが多い．一方，本シリーズで扱う exact-curved FEM では，
\[ F_K=\Psi_K\circ\Phi_{T_K} \]と分けて，affine core \(T_K\) と curved correction \(\Psi_K\) を別々に扱う．
この分解により，
- 異方性は affine core 側で見る
- 曲率は curved correction 側で見る
- 実装で使われる coordinate map を，解析しやすい形に整理できる
という利点が得られる．
次回以降は，この視点を使って，Ritz 射影や有限要素解の \(L^\infty\) 評価，多様体上の FEM への展開についても考えていく．
参考文献
この回では，isoparametric FEM と exact-curved FEM の関係を整理した． 背景となる文献として，曲線要素，isoparametric FEM，exact curved interpolation，および本シリーズの基礎となる筆者の投稿中論文を挙げておく．
- C. Bernardi, Optimal finite-element interpolation on curved domains, SIAM Journal on Numerical Analysis, 26(5), 1212--1240, 1989.
- 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.
- M. Zlámal, Curved elements in the finite element method. I, SIAM Journal on Numerical Analysis, 10(1), 229--240, 1973.
- M. Lenoir, Optimal isoparametric finite elements and error estimates for domains involving curved boundaries, SIAM Journal on Numerical Analysis, 23(3), 562--580, 1986.
- H. Ishizaka, K. Kobayashi and T. Tsuchiya, Anisotropic interpolation error estimates using a new geometric parameter, Japan Journal of Industrial and Applied Mathematics, 40, 475--512, 2023.
- H. Ishizaka, Exact-curved Lagrange finite elements for the Poisson problem in two dimensions, submitted.
