Newmark-β法による過渡応答解析
理論と物理
概要
先生、地震応答解析とか爆発荷重の解析って「過渡応答解析」ですよね？ ソルバーの設定にNewmark-βってよく出てくるんですが…
過渡応答解析は時間変動する荷重に対する構造の動的応答を求める解析だよ。地震だけじゃなく、衝撃荷重・風荷重・機械振動の起動停止なんかも対象になる。Newmark-β法は1959年にNewmarkが提案した陰的時間積分スキームで、構造解析の過渡応答解析の世界標準アルゴリズムだ。パラメータ β と γ を調整することで安定性と精度のトレードオフを制御できる。
「陰的」と「陽的」の時間積分の違いって何ですか？ LS-DYNAは「陽解法」ですよね？
陽解法（Explicit）は次ステップの加速度を現在の情報だけで直接計算できるが、安定性のために時間刻みが非常に小さくないといけない（CFL条件）。衝突・爆発などの超短時間現象向け。陰解法（Implicit）のNewmark-β法は現在と次ステップの両方を使うため、連立方程式を解く必要があるが、大きな時間刻みで安定して計算できる。地震（数十秒）など比較的ゆっくりの現象に向いている。
「陰解法は連立方程式を解く必要がある」とのことですが、毎タイムステップで剛性マトリクスを逆行列にするんですか？ コストが心配です。
線形問題（線形弾性）であれば、実効剛性マトリクス $\hat{K}$ を最初に一度だけLU分解しておけば、毎ステップはForward-Backward置換（$O(n)$相当）だけで済む。だから線形過渡応答では陰解法でも非常に高速だよ。問題は非線形解析（大変形・材料非線形・接触）で、この場合はステップごとにNewton-Raphson反復が必要になってコストが増加する。LS-DYNAがExplicitを使う「衝突解析」は非線形かつ超短時間のため、陽解法の方が有利なんだ。
運動方程式の時間積分
過渡応答解析の支配方程式を教えてください。
FEMによる多自由度系の運動方程式：
$[M]$ は質量マトリクス、$[C]$ は減衰マトリクス（Rayleigh減衰: $C = \alpha_M M + \beta_K K$）、$[K]$ は剛性マトリクス、$\{F(t)\}$ は時間依存外力ベクトル。変位 $\{u\}$、速度 $\{\dot{u}\}$、加速度 $\{\ddot{u}\}$ が未知数。
Newmark-β公式
Newmark法の具体的な更新式を教えてください。
Newmarkの予測式はこうなる：
$\{\ddot{u}\}_{n+1}$ は未知数だから、$\{u\}_{n+1}$ の表現を運動方程式に代入して連立方程式を解く。実効剛性マトリクス：
を構築して $\hat{K} \{u\}_{n+1} = \hat{F}_{n+1}$ を解く（$\hat{F}$ は既知量で構成される右辺ベクトル）。
$\hat{K}$ にはM、C、Kすべてが入っているんですね。これを一度因子分解しておけば毎ステップ高速に解けるということですか？
その通り。線形問題では $\hat{K}$ は時間に依存しない定行列だから、最初にLU（またはCholesky）分解しておけば毎タイムステップの解法コストは $O(n)$ の前後代入だけ。ただし $\Delta t$ が変化したり（アダプティブ時間積分）、非線形問題だと $\hat{K}$ が変化するため毎ステップ再因子分解が必要になる。これが非線形過渡応答でコストが跳ね上がる理由だ。
数値解法と実装
安定性と精度：βとγの選択
βとγをどう設定すればいいですか？ ソルバーのデフォルト値を使えばいいんですか？
代表的な設定を示すよ。
| β | γ | 手法名 | 安定性 | 精度 | 推奨用途 |
|---|---|---|---|---|---|
| 0 | 0 | 前進差分（陽的） | 条件付き | 1次 | —（Explicit用） |
| 1/6 | 1/2 | 線形加速度法 | 条件付き | 2次 | —（不安定リスク） |
| 1/4 | 1/2 | 定加速度法（台形則） | 無条件 | 2次 | 構造一般（最も一般的） |
| 1/4+ε | 1/2+ε | 数値減衰付き | 無条件 | 1次へ劣化 | —（HHT-αを使う方が合理的） |
実務では $\beta = 1/4, \gamma = 1/2$（定加速度法）が標準。ただし数値減衰ゼロのため高次モードが人工的に振動し続ける場合がある。その場合は HHT-α 法（次節）を使う。
HHT-α法への発展
HHT-α法ってどういうものですか？
Hilber-Hughes-Taylor（HHT）法はNewmark法に $\alpha$ パラメータを追加した改良版だよ。運動方程式の右辺を現ステップと次ステップの加重平均にすることで、精度を保ちながら高次モードに対して人工的な数値減衰を導入できる：
$\alpha = 0$ はNewmark定加速度法と同等。$\alpha = -1/3$ は最大の数値減衰。Abaqus・COMSOL・Ansysではこのデフォルトがよく使われている（$\alpha \approx -0.05$〜$-0.1$）。
「高次モードに数値減衰を導入」というのが重要なのはなぜですか？
FEMモデルには解析対象周波数をはるかに超える高次モードが無数にある。これらは物理的に意味がなく、モデルの数値的なアーティファクトに近い。Newmark定加速度法（α=0）ではこれらが数値的に全く減衰せず、振動ノイズとして結果に残る。HHT-αの数値減衰は「高次モードほど強く減衰させ、低次モードには影響しない」という理想的な特性を持つ。2次精度も維持されるから、ほぼすべての実務でNewmarkより優れている。
Rayleigh減衰の設定
Rayleigh減衰の係数 $\alpha_M$ と $\beta_K$ ってどう決めればいいんですか？
2つの基準周波数 $\omega_1$（低次モード）、$\omega_2$（解析帯域上限モード）での減衰定数 $h$（構造減衰比）を指定して逆算する：
$\omega_1 = 2\pi f_1$、$\omega_2 = 2\pi f_2$（Hz換算）。鋼構造で典型的には $h = 0.02$〜$0.05$（2〜5%）。注意点として $\omega_1$ と $\omega_2$ の間のモードだけが目標減衰比 $h$ に近くなり、この範囲外では減衰が過小または過大になる。広い周波数帯が必要な場合は、帯域ごとに別の減衰を設定するモーダル減衰の方が合理的。
実践ガイド
地震応答解析の実践
地震応答解析のタイムステップはどう決めればいいですか？
実践的なガイドラインを示すよ。
- $\Delta t$ の基準：解析したい最高振動数の1/10以下。例：100Hzまで精度が必要なら $\Delta t \le 0.001$s。一般的には最小固有周期 $T_{min}$ の1/10が目安。
- 地震波の入力：加速度時刻歴（観測記録またはAI生成波形）を基盤または支点変位として与える。サンプリング周波数に注意（100Hz以上推奨）。
- 減衰の設定：Rayleigh減衰係数 $\alpha_M, \beta_K$ を1次・2次モードの減衰定数（h = 2〜5%）から逆算。$\alpha_M = 2h\omega_1\omega_2/(\omega_1+\omega_2)$、$\beta_K = 2h/(\omega_1+\omega_2)$。
- 応答確認：変位・速度・加速度応答スペクトルを作成し、設計規準との整合性を確認。
日本の耐震設計で使われる応答スペクトルとNewmark法の結果をどう照合すればいいですか？
応答スペクトルはモーダル解析（SRSS/CQC法）で求めるのが標準だが、Newmark法の時刻歴解析との対応は以下のように行う：
- Newmark時刻歴解析で計算した加速度応答 $a(t)$ から、1DOF系応答スペクトル $S_a(T, h)$ を逆算して照合
- 日本建築基準法（告示波形）を入力地震波として使う場合、3波以上の解析の最大値で評価するのが要求
- 時刻歴解析はモーダルスペクトル解析の「モードの組み合わせ誤差」がない点が利点
特に高次モードが効いてくる超高層建築（60m以上）や不整形建築では、モーダルスペクトル解析よりNewmark時刻歴解析の方が精度が高い。
実践チェックリスト
- タイムステップ Δt ≤ 最高解析周波数の1/10 に設定
- Rayleigh減衰係数 αM、βK を目標周波数帯と減衰定数から逆算
- 数値減衰が必要な場合はHHT-αを使用（α = -0.05〜-0.1）
- 初期条件（変位・速度）が正しく設定されているか確認
- 地震波入力の場合：支点加速度入力か等価節点力かを明確に区別
- 線形問題：有効剛性行列の因子分解が一度で完了しているか確認
- 非線形問題：各ステップのNewton反復収束を監視（残差基準≤10⁻⁴推奨）
- 応答スペクトルを計算して設計基準との照合を実施
- 解析時間全体のエネルギー収支を確認（入力エネルギー＝ひずみ＋運動＋散逸）
ソフトウェア比較
主要ツールでのNewmark法の実装を比較してください。
主要ツールの比較を見てみよう。
| ツール | 過渡応答ソルバー | デフォルト時間積分 | 特徴 |
|---|---|---|---|
| MSC Nastran（SOL 109/129） | Newmark, HHT | Newmark β=1/4 | モーダル過渡応答も搭載、Sturm確認機能 |
| Ansys Mechanical（Transient） | HHT-α（Newmark系） | α=-0.1 | 自動時間ステップ制御、収束問題の自動対処 |
| Abaqus/Standard（*DYNAMIC） | HHT-α | α=-0.05 | 非線形構造との統合、接触・大変形対応 |
| COMSOL（Solid Mechanics） | Generalized-α | α=-0.05 | マルチフィジックス連成が容易 |
| OpenSees | Newmark, HHT, Collocation | Newmark | 地震工学専用、オープンソース、RC非線形に強い |
先端技術
Newmark法の発展形や最新の時間積分法はどんなものがありますか？
いくつかの重要な発展がある。
- Generalized-α法（Chung & Hulbert, 1993）：HHT-αをさらに一般化。Spectral radiusで無条件安定かつ高次モードを選択的に減衰。現代ソルバーのデフォルト。
- GSSSS法（Tamma et al.）：さらに一般的な枠組みで設計自由度を増やしたスキーム。
- 陰的-陽的分離（IMEX）法：剛性の高い部分を陰的、その他を陽的に扱うハイブリッド。FEM-DEM連成などに有効。
- 機械学習ベースの時間刻み適応：PINNを用いた適応時間積分法が研究段階にある。
デジタルツインや構造ヘルスモニタリングとNewmark法はどう関係しますか？
構造ヘルスモニタリング（SHM）との組み合わせが近年急速に進展している。センサーで計測した加速度データをNewmark法のFEMモデルに同化（データアシミレーション）することで、リアルタイムに構造状態を更新する「デジタルツイン」が実用化されつつある。具体的には：
- カルマンフィルタ × Newmark法：加速度センサーと有限要素法を融合して未計測点の変位推定
- モデル更新（Model Updating）：実測固有振動数と解析固有振動数の差を最小化してFEMパラメータを修正
- 損傷識別：Newmark時刻歴応答の差異から亀裂・ゆるみの位置を同定
橋梁・大型構造物での適用事例が増えており、今後の主要応用分野になる。
Newmark法と1971年San Fernando地震
Nathan Newmarkは1959年に本手法を提案したが、それが実際の設計基準に影響したのは1971年のSan Fernando地震（M6.6）の後だ。この地震でOlive View病院が崩壊し、動的解析の重要性が一気に高まった。Newmarkはその後、確率論的地震ハザード解析（PSHA）の基礎も作り、現代の耐震工学の父と呼ばれる。
トラブルシューティング
過渡応答解析でよくあるトラブルはどんなものですか？
数値減衰不足による高次モード振動、タイムステップが大きすぎる精度低下、減衰設定の間違い、境界条件の不整合など様々なトラブルがある。
| 症状 | 主な原因 | 対策 |
|---|---|---|
| 高次モードの数値振動ノイズ | Newmark β=1/4の数値減衰ゼロ | HHT-α法（α=-0.05〜-0.1）に変更 |
| 応答が過大評価される | Δtが大きすぎ（精度不足） | Δtを1/2〜1/5に縮小してメッシュ収束確認 |
| 非線形解析が収束しない | Newton-Raphson反復が発散 | Δtを小さくする、アークレング法（弧長制御法）に変更 |
| 地震応答が過小評価される | Rayleigh減衰の設定ミス（過減衰） | αM・βK を再計算、モーダル減衰の確認 |
| 初期から大きな変動が出る | 初期条件の不整合（静的釣り合いが未適用） | 解析前に静的自重荷重での釣り合いステップを追加 |
数値減衰問題、タイムステップ精度、Rayleigh減衰設定、収束困難など詳細解説