モード法過渡応答解析(モーダルトランジェント)
理論と物理
モード法過渡応答とは何か
先生、「モード法の過渡応答解析」って周波数応答のモード法とどう違うんですか? 名前は似てるけど、やってることが違うのか同じなのかわからなくて。
いい質問だね。共通点は「固有モードに展開して問題を縮約する」ところ。違いは入力と出力の領域だ。周波数応答は「一定周波数で加振し続けたときの定常応答」を求める。過渡応答は時間変化する任意の外力——衝撃、地震波、ステップ荷重——に対する時刻歴の応答を求める。
具体的にはどんな場面で使うんですか?
典型例を3つ挙げよう。(1) ロケット打上げ時の衛星振動——エンジンの推力変動が数秒間にわたって構造を加振する。(2) 自動車のドア閉め衝撃——ドアを閉めた瞬間のインパルスが車体全体に伝わる過渡応答。(3) 建築の地震応答——弾性範囲なら数十モードで免震建物の挙動を予測できる。いずれも「荷重が時間とともに変化する線形問題」で、モード法が威力を発揮する。
支配方程式とモード展開
出発点は多自由度系の運動方程式だ:
$[M]$:質量行列、$[C]$:減衰行列、$[K]$:剛性行列、$\{F(t)\}$:時間変化する外力ベクトル。$N$自由度なら、毎時間ステップで$N \times N$の連立方程式を解くことになる。これが直接法だ。
100万自由度で1万ステップだと、毎回100万×100万の行列を解くんですか? それは大変ですね……。
そこでモード重ね合わせ法(Mode Superposition Method)の出番だ。変位を固有モードで展開する:
$\{\phi_i\}$は第$i$固有モード形、$q_i(t)$はモード座標(一般化座標)、$m$は採用モード数($m \ll N$)。この展開を運動方程式に代入して、左から$[\Phi]^T$を掛けると——
モードの直交性を使って、連立方程式が分離されるんですよね?
その通り。質量正規化されたモードに対して$[\Phi]^T[M][\Phi] = [I]$、$[\Phi]^T[K][\Phi] = \mathrm{diag}(\omega_i^2)$が成り立つ。減衰が比例減衰(Rayleigh型)なら$[\Phi]^T[C][\Phi] = \mathrm{diag}(2\zeta_i\omega_i)$も対角化される。結果として、各モードが独立な1自由度系になる:
$f_i(t) = \{\phi_i\}^T\{F(t)\}$がモード荷重(一般化力)だ。100万自由度が30個の1自由度系に縮約される——計算時間は劇的に短縮される。
モード参加係数と有効質量
「何モード使えば十分か」はどうやって判断するんですか?
鍵になるのがモード参加係数と有効質量だ。ある方向$\{d\}$(例:X方向=[1,0,0,…])に対する第$i$モードの参加係数は:
有効質量は$M_{\mathrm{eff},i} = \Gamma_i^2$。全モードの有効質量を合計すると全質量になる:$\sum_{i=1}^{N} M_{\mathrm{eff},i} = M_{\mathrm{total}}$。実務では有効質量の累積が全質量の90%以上になるモード数を採用するのが基本ルール。原子力規格(NQA-1)や宇宙規格(ECSS-E-ST-32C)でもこの90%基準が明記されている。
90%に達するのに何モードくらい必要ですか?
構造によるが、建築なら10〜30モード、機器・衛星なら50〜200モード、配管系は300モード以上必要なこともある。配管は局所モードが多くて有効質量が分散するんだ。ただしモード数を増やしても計算時間は線形にしか増えないから、直接法と比べれば圧倒的に速い。
直接法との比較
| 項目 | モード法(Modal Transient) | 直接法(Direct Transient) |
|---|---|---|
| 支配方程式 | $m$個の独立な1自由度系 | $N \times N$の連立方程式を毎ステップ解く |
| 計算速度(100万DOF) | 非常に高速(モード数に比例) | 遅い($N$に強く依存) |
| 非線形対応 | 不可(接触・塑性・大変形) | 対応可能 |
| 周波数範囲 | 採用モードの最高固有振動数まで | 時間刻みで決まる全帯域 |
| 出力の柔軟性 | 各モードの寄与を個別に評価可能 | 物理量を直接出力 |
| 減衰の扱い | モード減衰(比例減衰前提) | 任意の減衰行列を適用可能 |
直接法と比べてモード法はどれくらい速いんですか? 具体的な数字が知りたいです。
100万自由度のモデルで100モードに縮約すれば、時間積分のコストは$1/10{,}000$以下になる。もちろん固有値解析の前計算コストが加わるが、それでもトータルで直接法の1/50〜1/100程度の時間で済むことが多い。ただし注意——衝撃荷重のように高周波成分が支配的な問題では、必要モード数が膨大になって速度の優位性が薄れる。その場合は直接法か陽解法を検討すべきだ。
NASAが育てたモード重ね合わせ法
モード重ね合わせ法が大規模構造に本格適用されたのは1960年代のNASAのアポロ計画・サターンVロケット開発がきっかけだ。当時のコンピュータで100万自由度の直接法は不可能だったが、構造を数十〜数百モードに縮約することで振動解析が実用化された。NastranのSOL 112(モード法過渡応答)はこの時代に生まれたソルバーアーキテクチャが原型であり、60年以上経った今も航空宇宙業界の標準手法として使われている。
数値解法と実装
Duhamel積分と解析解
各モードの運動方程式を分離できたとして、その1自由度系はどうやって解くんですか?
解析的にはDuhamelの積分(たたみ込み積分)で解ける:
ここで$\omega_{di} = \omega_i\sqrt{1-\zeta_i^2}$は減衰固有角振動数だ。物理的に言えば、「過去のすべての瞬間$\tau$に加わった力$f_i(\tau)$が、指数的に減衰しながら現在$t$まで影響を及ぼす」ことを表している。
積分を解析的に計算できるケースは限られますよね? 地震波のような複雑な波形だと……。
その通り。ステップ荷重や三角波のような単純な波形なら閉じた形で解けるが、実際の荷重波形は数値積分するのが普通だ。
数値時間積分の実装
各モードの時間積分は1自由度系だから、Newmark-$\beta$法や中央差分法をスカラー演算で実行する。$m$モード × $n_t$時間ステップの計算で、各計算はスカラーの加減乗除のみ。Newmark法($\beta=1/4, \gamma=1/2$の台形則)を適用すると:
1自由度なのでこれはスカラー計算——行列の分解は一切不要だ。
各モードを並列に計算できますよね? GPUで加速とかできるんですか?
原理的には完全並列化可能だ。各モードは互いに独立だからね。実際にNastran 2023以降ではモード法のGPUアクセラレーションがサポートされている。ただし実務的には、ボトルネックは時間積分よりも固有値解析と出力の物理座標変換にあることが多い。
残余ベクトル(RESVEC)の役割
モード数を増やせば精度が上がるのはわかりましたが、「残余ベクトル」って何ですか? ソルバーのマニュアルによく出てくるんですけど。
すごく重要な概念だ。モード法では$m$個のモードしか使わないから、荷重ベクトル$\{F\}$のうち採用モードで表現できない成分が切り捨てられる。残余ベクトル(Residual Vector、RESVEC)は、その「切り捨てられた成分」を擬似モードとして追加するものだ。
具体的にはどう計算するんですか?
荷重ベクトルから採用モードの寄与を差し引いた残差を求め、それを剛性行列の逆行列(静的変位)として解く:
この$\{r\}$を正規化して追加モードとする。効果は劇的で、特に荷重点付近の局所応答が大幅に改善される。NastranのPARAM,RESVEC,YES、AbaqusのRESIDUAL MODESオプションが対応する。衝撃解析では必ず有効にすること。
出力の物理座標変換
モード座標$q_i(t)$が得られたあと、実際の変位や応力はどうやって求めるんですか?
物理座標への変換は:
応力は各モードの応力成分にモード座標を掛けて重ね合わせる:$\{\sigma(t)\} = \sum_{i=1}^{m} \{\sigma_i\}\, q_i(t)$。ここで重要な実務テクニックがある——全節点・全時刻ステップで物理座標変換するとデータ量が膨大になるから、着目点(応力評価点やセンサー位置)のみで変換するのが定石だ。NastranではSETカードで出力節点を限定する。
出力データ爆発との戦い
100万自由度 × 1万ステップ × 6応力成分 = 600億個の浮動小数点数。倍精度で480GBだ。モード法で計算は速く終わっても、出力で破綻することがある。ある宇宙機メーカーのエンジニアは「解析3分、後処理30分、出力ファイルのコピーに2時間」と嘆いていた。モード座標$q_i(t)$だけを保存して、後処理時に必要な箇所だけ物理座標変換する運用が合理的だ。
実践ガイド
モード数の決め方
有効質量90%が基準なのはわかりましたが、実務でモード数をどう決めるか、もっと具体的に教えてください。
判断基準は3つある。(1) 有効質量90%以上——これは最低ライン。(2) 荷重の周波数帯域——入力荷重をFFTして、エネルギーの99%が含まれる最高周波数$f_{\max}$を求め、その$1.5\sim2$倍の固有振動数までのモードを採用する。(3) 直接法との比較検証——小規模モデルで直接法と結果を突き合わせ、必要なモード数を事前に見積もる。
例えばロケットの衛星フェアリング分離の衝撃解析なら?
フェアリング分離は2000Hz程度までの帯域が要求される。衛星構造の固有振動数は低次が5〜50Hzくらいだから、2000Hzまでの全モード——場合によっては数千モードが必要になる。このレベルだとモード法の速度優位性が薄れるので、RESVEC併用か直接法への切り替えを検討する。一方、打上げ時の正弦波振動試験(5〜100Hz)なら50〜100モードで十分だ。
時間刻みの設定
時間刻み$\Delta t$の設定基準はありますか?
基本ルールは着目する最高振動数$f_{\max}$の周期の1/20以下にすることだ:
もう1つ重要なのが、荷重波形の変化を捉えるのに十分細かいこと。ステップ荷重の立上り時間が0.001秒なら、$\Delta t = 0.0001$秒程度は必要。モード法はNewmark法(無条件安定)を使うことが多いので安定性の制約はないが、精度の制約がある。
減衰の設定と注意点
モード法での減衰設定ってどうするんですか? レイリー減衰を使うんですよね。
モード法の大きな利点の1つは、各モードに個別の減衰比$\zeta_i$を設定できることだ。レイリー減衰だと$\zeta_i = \frac{\alpha}{2\omega_i} + \frac{\beta\omega_i}{2}$で、低次と高次で減衰比が歪む。モード法なら「1次モードは2%、5次モードは5%」のように物理的に合理的な値を直接指定できる。
減衰比の値はどうやって決めるんですか? 実測データがないことも多いですよね。
実測がなければ経験値を使う。鋼構造:0.5〜2%、コンクリート:3〜5%、ボルト締結構造:2〜4%、ゴム・樹脂:5〜15%。ただしこれはあくまで目安であって、振幅依存性がある。小振幅では低く、大振幅では高くなる。宇宙機器の設計基準(NASA-HDBK-7005)には材料別の推奨値リストがある。
SRS(衝撃応答スペクトル)評価
SRS(Shock Response Spectrum)ってモード法と関係がありますか?
非常に密接な関係がある。SRSは「ある衝撃入力に対して、固有振動数$f_n$の1自由度系が示す最大応答」を$f_n$の関数としてプロットしたものだ。モード法過渡応答の副産物として各モードの最大応答$|q_i|_{\max}$が得られるから、そこからSRSを直接構成できる。宇宙機器の環境試験条件(NASA-STD-7003、ECSS-E-ST-10-03C)はSRSで規定されるので、モード法過渡応答→SRS評価は宇宙業界の標準ワークフローだ。
実務チェックリスト
モード法過渡応答解析の実務で確認すべき項目をまとめよう:
- 固有値解析が正常完了 — 負の固有値、剛体モードの有無を確認
- 有効質量が3方向とも90%以上 — ソルバーの出力テーブルで確認
- RESVEC(残余ベクトル)が有効 — 衝撃荷重では必須
- 時間刻み$\Delta t \leq T_{\min}/20$ — $T_{\min}$は着目する最短周期
- 荷重波形のサンプリングが十分 — 立上り/立下りの急峻な部分を捕捉
- モード減衰$\zeta_i$が全モードに設定済み — 未設定だと$\zeta=0$扱い
- 初期条件(変位・速度)が正しい — 通常はゼロだが、先行ステップの引き継ぎに注意
- 直接法との比較検証 — 小規模モデルで一致を確認
H-IIAロケットの衛星振動解析
JAXA H-IIAロケットの打上げ時、衛星構造には20〜2000Hzの振動環境が作用する。衛星側は1000〜5000自由度に縮退したCraig-Bampton縮約モデルをロケット側の発射台モデルと結合し、上位200モードによるモード重ね合わせ過渡応答で打上げ荷重を評価する。この解析で得られたSRSが衛星の試験レベル(Qualification/Acceptance)の根拠となる。
ソフトウェア比較
Nastran SOL 112
NastranのSOL 112はモード法過渡応答の業界標準だ。入力の骨格はこうなる:
SOL 112 $ モード法過渡応答
CEND
METHOD = 10 $ 固有値解析法の参照
TSTEP = 100 $ 時間ステップ定義
DLOAD = 200 $ 動的荷重定義
PARAM,RESVEC,YES $ 残余ベクトル有効
BEGIN BULK
EIGRL, 10, , , 50 $ 50モードを抽出
TSTEP, 100, 1000, 0.001 $ 1000ステップ、Δt=0.001s
TLOAD1, 200, 300, , 0, 400 $ 時間荷重
TABLED1, 400, ... $ 荷重波形テーブル
ポイントはPARAM,RESVEC,YES。これがないと衝撃解析で結果が直接法と乖離する。
Abaqus *MODAL DYNAMIC
Abaqusでは2つのステップを連続実行する:
** Step 1: 固有値解析
*STEP
*FREQUENCY
50, , $ 50モード抽出
*RESIDUAL MODES $ 残余モード有効
*END STEP
** Step 2: モード法過渡応答
*STEP
*MODAL DYNAMIC
0.001, 1.0 $ Δt=0.001s, 全時間=1.0s
*MODAL DAMPING
1, 50, 0.02 $ モード1〜50に減衰比2%
*CLOAD
... $ 時間荷重
*END STEP
Abaqusの特徴は*MODAL DAMPINGでモードごとに減衰比を指定できること。モード範囲を分けて異なる減衰比を設定するのも容易だ。
Ansys MSUP Transient
Ansys Mechanicalだとどういう設定になりますか?
Ansys MechanicalではGUIで「Modal」解析→「Transient (Mode Superposition)」の順に設定する。APDLレベルではANTYPE,MODAL→ANTYPE,TRANSの流れだ。APDLコマンドの骨格:
/SOLU
ANTYPE,MODAL ! 固有値解析
MODOPT,LANB,50 ! Block Lanczoで50モード
MXPAND,50 ! 50モード展開
SOLVE
FINISH
/SOLU
ANTYPE,TRANS ! 過渡応答
TRNOPT,MSUP ! MSUP(モード重ね合わせ)法
DELTIM,0.001 ! Δt=0.001s
MDAMP,ALL,0.02 ! 全モード減衰比2%
F,NODE,FY,LOAD ! 荷重定義
TIME,1.0
SOLVE
TRNOPT,MSUPが直接法ではなくモード法を選択するキーワードだ。
機能比較表
| 機能 | Nastran SOL 112 | Abaqus *MODAL DYNAMIC | Ansys MSUP |
|---|---|---|---|
| 残余ベクトル | PARAM,RESVEC,YES | *RESIDUAL MODES | RESVEC,ON |
| モード別減衰 | TABDMP1 | *MODAL DAMPING | MDAMP |
| SRS出力 | PARAM,SRS(組込み) | Pythonスクリプト | POST26マクロ |
| GPU加速 | 2023〜(固有値) | 制限的 | 2024R1〜 |
| Craig-Bampton結合 | SEAM+SOL 112 | *SUBSTRUCTURE | CMS |
| 得意分野 | 航空宇宙・防衛 | 汎用・研究 | 自動車・一般産業 |
結局どれを使えばいいんですか?
使い分けの基本はこうだ。航空宇宙のSRS評価・衛星振動→Nastran一択。SRS計算が組み込みで、業界標準のCLA(Coupled Load Analysis)ワークフローと整合する。汎用的なモード法過渡+直接法との切り替えが多い問題→Abaqus。同一モデルでモード法と直接法を容易に切り替えられる。GUIベースの業務フロー・自動車→Ansys Mechanical。Workbenchのワークフロー管理が便利だ。
モード法から直接法への切り替え判断
あるプロジェクトで「モード法でやったのに直接法と2割も結果が違う」という相談を受けたことがある。原因はモデルに非線形ばね要素(GAP要素)が混在していたことだった。GAP要素は接触/離間で剛性が変わるため、固有モードの直交性が崩れる。結局直接法に切り替えて解決したが、「モード法が適用可能か」の判断を最初にやるべきだった。モデルに非線形要素が1つでも含まれていたら、モード法は使えないと思ったほうがいい。
トラブルシューティング
直接法と結果が合わない
モード法と直接法で同じモデル・同じ荷重を解いたのに、応答のピーク値が30%も違います。何が原因でしょう?
直接法との不一致は「モード法あるある」の第1位だ。原因を優先順に並べると:
- モード数の不足 — 有効質量の累積を確認。3方向全てで90%以上あるか? 特にZ方向(重力方向)だけ80%を切っていて見落とすケースが多い
- RESVECが無効 — Nastranのデフォルトは
RESVEC=YESだが、古いバルクデータを流用するとRESVEC=NOになっていることがある - $\Delta t$が大きすぎる — モード法の$\Delta t$が直接法のそれより粗いと、荷重ピークを捕捉できない
- 減衰の設定差 — 直接法のレイリー減衰とモード法のモード減衰で、各モードの実効減衰比が異なっている
- 非線形要素の混入 — GAP要素、CBUSH非線形、MPC依存の要素がモデルに含まれていないか確認
5番の非線形要素の混入が一番怖いですね……。ソルバーはエラーを出さないんですか?
Nastranは非線形要素をSOL 112で使うとWARNINGは出すが、計算自体は止まらない。線形化した状態で計算を続行する。だから結果は出るが、直接法と一致しない。「計算が回る=結果が正しい」とは限らない典型例だ。
衝撃荷重の初期応答が不正確
衝撃荷重を加えた直後の0.001秒間だけ、モード法の応答が直接法と全く合いません。その後は概ね一致するんですが……。
典型的なRESVEC不足の症状だ。衝撃荷重はフーリエ変換すると広帯域——数kHzまで成分が広がる。採用モードの最高固有振動数が500Hzだとすると、500Hz以上の成分が全て欠落する。RESVECはこの高周波成分の静的な寄与を補償するが、動的な寄与(衝撃直後の高周波振動)は補えない。対処法は3つ:
- RESVECを有効にする(まだなら)
- モード数を大幅に増やす(最高固有振動数を荷重帯域の2倍以上に)
- 直接法に切り替える(衝撃直後の高精度が必要なら)
応答が減衰しない
荷重を除荷した後も応答が全く減衰せず、永遠に振動し続けます。
減衰が設定されていない。モード法では減衰を明示的に指定しないとゼロ扱いになるソルバーが多い(Nastranのデフォルトは$\zeta=0$)。設定方法はソルバーごとに異なる:
- Nastran: TABDMP1カードでモード減衰テーブルを定義し、CASEでSDamp=IDを指定
- Abaqus: *MODAL DAMPING行でモード番号範囲と減衰比を指定
- Ansys: MDMPコマンドで個別or全モード一括設定
まず全モード一律2%で試し、実測データがあれば調整するのが実務的だ。
負の固有値が出る
固有値解析で負の固有値が出ました。モード法過渡応答を実行していいんでしょうか?
絶対にそのまま実行してはいけない。負の固有値は剛性行列の正定値性が崩れていることを意味する——つまりモデルのどこかが不安定だ。よくある原因:
- 拘束不足 — 剛体モードが残っている(6つの剛体モードの固有値≒0は正常だが、それ以外の負値は異常)
- 要素のワーピング/歪み — ヤコビアンが負になるほど歪んだ要素がある
- 初期応力の導入ミス — プリストレスの符号が逆、または過大
負の固有値があるとモード座標$q_i(t)$が発散する。まず固有値解析の段階で原因を解消すること。
トラブル解決の鉄則:最小再現モデル
モード法過渡応答でトラブルが起きたら、まず「最小再現モデル」を作ること。100万要素のモデルで1時間回して結果がおかしい、をデバッグするのは非効率的だ。問題を再現する最小構成(10要素の片持ち梁 + 同じ荷重条件)で数秒で解析を回し、同じ症状が出るか確認する。症状が出たら、そのモデルでパラメータを1つずつ変えて原因を特定する。出なかったら、元のモデルから要素を段階的に削って、どの要素を加えると症状が出るかを突き止める。プロのデバッグは「引き算」だ。
関連トピック
なった
詳しく
報告