0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

ニュートン法と有限要素法(FEM)

Posted at

1. 序論 ― 連続体力学を最適化問題として解く意義

最適化の視点からは、未知ベクトル $\mathbf z$ の最小化問題

\min_{\mathbf z}\; \Phi(\mathbf z)

に対し、ニュートン法は

\Delta\mathbf z
  = -\bigl[\nabla^{2}\Phi(\mathbf z_n)\bigr]^{-1}
     \nabla\Phi(\mathbf z_n), \qquad
\mathbf z_{n+1}=\mathbf z_n+\Delta\mathbf z

という 2 次収束アルゴリズムを提供する。本稿では、連続体力学の平衡条件を同一形式の最小化として書き下し、連続最適化の既存実装で直接解けることを示す。


2. 取り扱う変数と仮定

記号 定義 役割
$\Omega_0$ 静止した 3 次元領域 積分領域
$\boldsymbol\chi:\Omega_0!\to!\mathbb R^3$ 変形写像 変位の基礎
$\mathbf u = \boldsymbol\chi - \mathbf X$ 変位ベクトル場 最適化変数
$\mathbf F = \partial\boldsymbol\chi/\partial\mathbf X$ 変形勾配 エネルギー $\Psi$ の入力
$\mathbf E = \tfrac12(\mathbf F^{!T}\mathbf F - \mathbf I)$ Green 歪 小変形極限で使用
補足 A ― 材料モデルと枠の無関係性
  • 材料は履歴を持たない超弾性体とし,エネルギー密度 $\Psi = \Psi(\mathbf F)$ のみで記述する。
  • 剛体回転 $\mathbf R!\in!\mathrm{SO}(3)$ に対して $\Psi(\mathbf R\mathbf F) = \Psi(\mathbf F)$ を仮定する(objectivity)。

3. 無限次元を有限次元へ写像する手順

3.1 形状関数による空間離散化

\mathbf u(\mathbf X) = \sum_{I=1}^{n_{\text{dof}}} N_I(\mathbf X)\,\mathbf u_I,

により未知場を係数ベクトル

\mathbf u = \{\mathbf u_I\}\in\mathbb R^{n_{\text{dof}}}

へ縮約する。

3.2 ガウス–ルジャンドル積分による数値積分

\int_{\Omega_e} f\,dV
  \;\approx\;
\sum_{g=1}^{n_g} w_g\,f(\boldsymbol\xi_g).
補足 B ― 積分点数とロッキング

次数 $p$ の形状関数には
$n_g \ge \lceil (2p+1)/2\rceil^3$
が経験的最小点数である。一次要素で大変形解析を行う場合、体積ロッキングを抑制する目的で選択積分や混合変数法の導入が推奨される。


4. 全ポテンシャルエネルギーとニュートン法

4.1 連続体としてのエネルギー

\Pi[\mathbf u] =
\int_{\Omega_0} \Psi(\mathbf F)\,dV
- \int_{\Omega_0} \mathbf B\!\cdot\!\mathbf u\,dV
- \int_{\Gamma_0^t} \mathbf T\!\cdot\!\mathbf u\,dA.

ここで $\mathbf B$ は体積力,$\mathbf T$ は表面力である。

4.2 離散スカラー関数

\Phi(\mathbf u)=
\sum_e \sum_g w_g\,
  \Psi\bigl(\mathbf F(\boldsymbol\xi_g;\mathbf u)\bigr)
-
\mathbf u^{\mathsf T}\mathbf f_{\text{ext}}.

4.3 ニュートン反復

\Delta\mathbf u
  = -\bigl[\nabla^{2}\Phi(\mathbf u_n)\bigr]^{-1}
     \nabla\Phi(\mathbf u_n), \qquad
\mathbf u_{n+1}=\mathbf u_n+\Delta\mathbf u.
補足 C ― 収束対策と分岐検出
  • ヘッセ行列が半正定値の場合や荷重増分が大きい場合,反復が発散し得る。線形検索・ダンピング・アークレングス法が有効である。
  • $\det\nabla^{2}\Phi \to 0$ は分岐点接近の指標となるため,最小固有値を監視してクリティカルロードを推定できる。

5. 勾配とヘッセ行列の具体形

5.1 残差ベクトル

\mathbf r(\mathbf u)=
\nabla\Phi(\mathbf u) =
\sum_e \sum_g w_g\,
  \mathbf B_g^{\mathsf T}\,\mathbf P_g
- \mathbf f_{\text{ext}},
\mathbf P = \frac{\partial\Psi}{\partial\mathbf F}\bigl(\mathbf F\bigr)

第 1 Piola 応力である。

5.2 一貫接線剛性

\mathbf K(\mathbf u) =
\nabla^{2}\Phi(\mathbf u) =
\sum_e \sum_g w_g\,
  \mathbf B_g^{\mathsf T}\,
  \mathbf D_g\,
  \mathbf B_g\,
  |\det\mathbf J_g|,
\mathbf D = \frac{\partial^2\Psi}{\partial\mathbf F^2}\bigl(\mathbf F\bigr)

材料剛性テンソルである。$\mathbf K$ はニュートン反復におけるヘッセ行列に相当し、完全微分を用いることで 2 次収束が保証される。


6. 小変形理論への極限

変位が微小なとき

\mathbf F \approx \mathbf I + \nabla\mathbf u, \quad
\Psi(\mathbf F) \approx \tfrac12\,\mathbf E:\mathbf C:\mathbf E,

と 2 次まで展開すると、要素剛性は

\mathbf K_e =
\int_{\Omega_e}
  \mathbf B^{\mathsf T}\,
  \mathbf C\,
  \mathbf B
  \,d\Omega.

ここで $\mathbf C$ は定数 4 階テンソルであり、線形有限要素法における剛性行列そのものとなる。なお得られる応力は 第 2 Piola–Kirchhoff 応力の 1 次項 に対応する。


7. まとめ

  1. 形状関数と数値積分により、連続体力学の平衡を有限次元スカラー関数 $\Phi(\mathbf u)$ に帰着した。
  2. ニュートン法をそのまま適用することで効率的に平衡解を得られる。
  3. 勾配・ヘッセ行列は、それぞれ応力 $\mathbf P$ と材料剛性 $\mathbf D$ の形で物理的意味をもつ。
  4. 小変形極限では標準的な線形 FEM への収束が確認できる。
0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?