Differential Forms with Applications to the Physical Sciencesのp.48の6番でz=f(x,y)と表せるような曲面についてガウス曲率Kと平均曲率Hを求める問題がある。これをsympyでやってみる。
続きを読むメッシュ上の数値積分
- 実際にはメッシュの頂点上に与えられた離散的な関数を球面調和関数で展開する
- 球面調和関数は連続なのでできる限りその情報を失わないようにしたい
- 連続関数のメッシュ上での積分はガウス求積が一般的に用いられる
- ガウス求積 - Wikipedia
- ガウス=クロンロッド求積法 - Wikipedia: 誤差見積もりが必要な場合
- Simpson法だと3種類の係数で2次の精度→一般にN種類の係数でN-1次の精度(積分点が少なくていい訳ではない)
- ガウス求積だとN個の重みとN個の積分点で2N-1次の精度(積分点はN個だけ)
- 重みには球面調和関数に用いたルジャンドルの陪関数を使う
- 積分点はSimpson法だと均一に分布するのに対し、ガウス求積だと境界部分により多く分布する
sympyで行列計算
物理数学の教科書を読んでいて文字入りの逆行列の計算が合わなかったのでsympyを使った。sympyがインストールされていなければSymPy Liveが少し遅いが便利。
from sympy import * l, m, s = symbols("λ μ σ") A = Matrix([[l + 2 * m, l, l], [l, l + 2 * m, l], [l, l, l + 2 * m]]) y = Matrix([s, 0, 0]) x = A.inv() * y print( x.symplify() )
Matrix([ [σ*(λ + μ)/(μ*(3*λ + 2*μ))], [ -λ*σ/(2*μ*(3*λ + 2*μ))], [ -λ*σ/(2*μ*(3*λ + 2*μ))]])
球面メッシュ
2015-04-06にあるRのコードをC++に翻訳した。
続きを読む