モンテカルロシミュレーション — CAE不確かさ定量化の基本手法

いい質問だね。ざっくり言うと、入力パラメータの確率分布からランダムにサンプリングして、FEMやCFDの解析を何百回・何千回と繰り返し実行する手法だ。出力のばらつき（分布）を統計的に求めるのが目的だよ。

そう、それがモンテカルロ法の最大の弱点でもある。ただし最大の強みは、入力変数が何次元になっても収束速度が変わらないことだ。感度解析だと変数が10個あれば指数関数的にコストが増えるけど、モンテカルロは10次元でも100次元でも同じペースで収束する。

例えば、自動車部品の材料強度のばらつき（変動係数CoV=5%程度）が疲労寿命にどう影響するかを評価するケース。ヤング率、降伏応力、板厚、荷重振幅など複数のパラメータが同時にばらつくとき、それらの「組み合わせ効果」が出力にどう伝播するかを調べるのがモンテカルロ法の得意領域だ。

モンテカルロ法の核心は大数の法則と中心極限定理だ。$N$ 個のサンプル $\mathbf{X}_1, \mathbf{X}_2, \ldots, \mathbf{X}_N$ を入力分布から独立にサンプリングして、それぞれに対してモデル $Y = f(\mathbf{X})$ を評価する。すると標本平均は：

この標本平均の標準誤差（Standard Error）は次の式で表される：

その通り。これが有名な$1/\sqrt{N}$ 収束則だ。重要なのは：

中心極限定理のおかげで、$N$ が十分大きければ（実務では $N \geq 30$ 程度から）標本平均は正規分布に近似できる。だから95%信頼区間はこう書ける：

相対的な推定精度（Relative Accuracy）を計算しよう。変動係数を $\text{CoV}_Y = \sigma_Y / \mu_Y$ とすると：

例えば出力の変動係数が $\text{CoV}_Y = 0.10$（10%）のとき：

目的によって全然違う。ざっくりまとめるとこんな感じだ：

だから稀事象の推定には純粋なモンテカルロ法ではなく、重点サンプリング（Importance Sampling）やサブセットシミュレーション（Subset Simulation）といった特殊な手法を使う。あるいはサロゲートモデルを挟んで計算コストを下げるんだ。

最もシンプルな方法はそうだ。各入力変数の確率分布から独立に擬似乱数を生成してサンプルを作る。これを純粋乱数サンプリング（Crude Monte Carlo / Simple Random Sampling）と呼ぶ。

鋭いね。特にサンプル数が少ないとき（$N < 100$ 程度）、純粋乱数だと入力空間のカバレッジにムラが出やすい。例えば2変数で100サンプルを振ると、分布域の隅に全然サンプルが行かないことがある。そこで登場するのがLHS（ラテン超方格サンプリング）だ。

LHSのアイデアはシンプルだ。各入力変数の分布を $N$ 個の等確率の層（stratum）に分割して、各層から必ず1点ずつサンプルを取る。そして各変数間の組み合わせをランダムにシャッフルする。

そのイメージで合ってる。数学的には、各変数 $X_j$（$j = 1, \ldots, d$）の累積分布関数 $F_j$ の値域 $[0, 1]$ を $N$ 等分して：

LHSのメリットをまとめると：

「分散低減法」として知られる手法群があるよ。標準誤差 $\text{SE} = \sigma_Y / \sqrt{N}$ の $\sigma_Y$ を実効的に小さくするか、推定量を工夫して同じ $N$ でも精度を上げるんだ。

応答関数が入力に対して単調（例：材料強度が上がれば寿命も伸びる）なら効果が大きい。$\mathbf{U}$ で大きな出力が出たとき、$1-\mathbf{U}$ では小さくなるから、2つの平均を取ると分散が相殺される。ただし非単調な関数だと逆効果になることもあるから注意だ。

基本的な流れはこうだ：

それが一番悩ましいところで、「ゴミを入れればゴミが出る（Garbage In, Garbage Out）」が最もクリティカルに効くのが入力分布の設定だ。理想的にはメーカーの品質データや実験データからフィッティングするんだけど、実務では以下の目安を使うことが多い：

強度は物理的に負の値を取れない。正規分布だと $\mu - 3\sigma$ あたりで負になる可能性があるけど、対数正規分布なら常に正の値しか生成しない。さらに強度データは右裾が長い（高い値の方に引っ張られる）傾向があるから、対数正規の方がフィットしやすいんだ。

2つのアプローチがある。1つ目は事前決定型で、目標精度から必要サンプル数を逆算する：

2つ目は逐次型。解析を進めながら累積平均と累積標準偏差をモニタリングして、安定したら打ち切る。実務ではこんなチェックを行う：

自動車のサスペンションアーム部品を例にしよう。決定論的な解析では「疲労寿命 $N_f = 2 \times 10^6$ サイクルで安全」と出たとする。でも実際の量産品には材料のばらつきがあるよね？

そこでモンテカルロ法だ。こういう手順でやる：

結果、平均寿命は $2.1 \times 10^6$ サイクルだけど、5パーセンタイル値（B5寿命）は $4.5 \times 10^5$ サイクルしかない——つまり20個に1個は設計寿命の1/4以下で壊れる可能性がある、ということが分かる。決定論的解析だけでは見えないリスクだよね。

主要なCAEプラットフォームにはUQ（不確かさ定量化）モジュールが組み込まれている。自分でPythonスクリプトを書く方法もあるし、GUI上で設定する方法もある。

実は最も柔軟で強力なのはオープンソースの組み合わせだったりする：

使える。Dakotaは「ブラックボックス」インターフェースで、入力ファイルのテンプレートにパラメータを埋め込み、ソルバーを外部コマンドとして呼び出し、結果ファイルから出力を読み取る。OpenFOAM、CalculiX、Code_Aster、なんでもOKだ。

そこでサロゲートモデル（代理モデル）を使う。少数の高精度シミュレーション（数十〜数百点）で応答曲面を構築し、モンテカルロサンプリングはその軽い代理モデルに対して行う。

PCEは応答を直交多項式で展開する。例えば入力が正規分布ならエルミート多項式を使う：

直交性のおかげで、平均と分散は係数から直接計算できる：$\mu_Y = c_{\mathbf{0}}$、$\sigma_Y^2 = \sum_{\boldsymbol{\alpha} \neq \mathbf{0}} c_{\boldsymbol{\alpha}}^2$。モンテカルロサンプリングすら不要になるんだ。ただし変数が多い（$d > 10$）と係数の数が爆発するから、スパースPCEや適応的手法を使う。

モンテカルロ法は完全並列（Embarrassingly Parallel）だ。各サンプルは独立に計算できるから、1,000コアあれば1,000サンプルを同時に回せる。これは他のUQ手法（例：適応的PCE）にはない大きな利点だ。

現場でよく見るトラブルをまとめよう：

選択肢は3つある：

その通り。モンテカルロ法は「統計的手法」だから、再現性の確保は科学的信頼性の基盤だ。シード値、サンプリング手法、サンプル数、入力分布のパラメータ——これらを全て記録しておくことが、V&V（検証と妥当性確認）の第一歩なんだよ。