convergence-study

収束性の調査

目標

メッシュまたは時間刻みを細かくしたときに、数値解が期待される収束率で収束することを検証するための、スクリプト駆動型の収束解析を提供します。

要件

• Python 3.10+
• サードパーティ製パッケージは不要 — スクリプトは Python 標準ライブラリ（math）のみを使用します。

収集する入力

入力	説明	例
格子間隔	メッシュサイズの系列（粗いものから細かいものへ）	0.4,0.2,0.1,0.05
時間刻み幅	dt 値の系列	0.04,0.02,0.01
解の値	各細分化レベルでの QoI	1.16,1.04,1.01,1.0025
期待収束次数	数値スキームの形式次数	2.0
安全係数	GCI の安全係数（デフォルトは 1.25）	1.25

スクリプトの出力（JSON フィールド）

スクリプト	主な出力
scripts/h_refinement.py	results.observed_orders、results.mean_order、results.richardson_extrapolated_value、results.convergence_assessment
scripts/dt_refinement.py	h_refinement と同じですが、時間収束を対象とします
scripts/richardson_extrapolation.py	results.extrapolated_value、results.error_estimate、results.observed_order
scripts/gci_calculator.py	results.observed_order、results.gci_fine、results.gci_coarse、results.asymptotic_ratio、results.in_asymptotic_range、results.extrapolated_value、results.notes

ワークフロー

1. 格子／時間刻みの細分化調査を実行し、少なくとも 3 段階で評価します
2. h_refinement.py または dt_refinement.py で 観測収束次数を計算します
3. スキームの期待収束次数と比較します
4. リチャードソン外挿で 離散化誤差を推定します
5. gci_calculator.py を使用して、正式な解の検証用に GCI を報告します
6. 収束結果と異常を 記録します

判断の指針

Do you have 3+ refinement levels?
+-- YES --> Run h_refinement.py or dt_refinement.py
|           +-- Observed order matches expected? --> Solution verified
|           +-- Order too low? --> Check: pre-asymptotic, coding error, insufficient resolution
|           +-- Order too high? --> Check: superconvergence or cancellation effects
+-- NO (only 2 levels) --> Use richardson_extrapolation.py with assumed order
                           (less reliable without order verification)

CLI の例

# Spatial convergence with 4 grid levels
python3 scripts/h_refinement.py --spacings 0.4,0.2,0.1,0.05 --values 1.16,1.04,1.01,1.0025 --expected-order 2.0 --json

# Temporal convergence with 3 timestep levels
python3 scripts/dt_refinement.py --timesteps 0.04,0.02,0.01 --values 2.12,2.03,2.0075 --expected-order 2.0 --json

# Richardson extrapolation with assumed 2nd-order
python3 scripts/richardson_extrapolation.py --spacings 0.02,0.01 --values 1.0032,1.0008 --order 2.0 --json

# GCI for 3-mesh verification
python3 scripts/gci_calculator.py --spacings 0.04,0.02,0.01 --values 1.0128,1.0032,1.0008 --json

エラー処理

エラー	原因	解決方法
spacings and values must have the same length	入力配列の長さが一致しません	同じ長さのリストを指定します
At least 2 refinement levels required	データポイントが少なすぎます	細分化レベルを追加します
Exactly 3 refinement levels required	GCI には 3 段階が必要です	細・中・粗の各レベルを指定します
Oscillatory convergence detected	収束が単調ではありません	メッシュ品質またはスキームを確認します

解釈の指針

状況	意味	対応
観測次数が期待次数と一致する	漸近領域にあることを示す最も強い証拠です	GCI を報告し、外挿します
観測次数 < 期待次数	漸近領域に達していないか、コードに不具合があります	さらに細分化するか、デバッグします
観測次数が負	解が発散しています	実装を確認します
GCI の漸近比が 1.0 に近い	下記の注意事項を参照してください	次数の比較で確認します
GCI の漸近比が 1.0 から大きく離れている	漸近領域にありません	さらに細分化します

漸近比に関する注意事項（一定の細分化比の場合）。 細分化比が等しい場合（r21 == r32。一般的なケース）、漸近比 AR = GCI_coarse / (r^p * GCI_fine) は代数的に f1/f2 に簡約されます。この場合、測定しているのは最も細かい 2 つの QoI 値の相対的な差だけであり、データが仮定したべき乗則に従っているかどうかではありません。そのため、観測次数が期待値から大きく外れていても、AR が 1.0 に近いと誤った安心感を与えることがあります。gci_calculator.py の JSON には、この点を示す notes エントリが出力されます。漸近領域にあるかを実際に判定するには、次の手順を実行します。

1. 観測次数 p をスキームの理論次数／期待次数と比較します。一致していれば、漸近領域にあることを示す有意な証拠になります。
2. より強い検証を行うには、系統的に細分化した格子を 4 つ以上使用し、連続する格子 3 段階の組み合わせで観測次数が一貫していることを確認します（h_refinement.py は組み合わせごとの次数に加えて mean_order も報告します）。

検証チェックリスト

• h_refinement.py／dt_refinement.py が results.mean_order を報告できるよう、体系的に細分化したグリッド／タイムステップを3つ以上使用し、各三つ組の results.observed_orders を記録しました（単一ペアのリチャードソン実行では次数を検証できません）。
• results.mean_order を記録し、results.convergence_assessment が PASS（観測次数がスキームの予想次数の10％以内）と読み取れることを確認しました。FAIL または unknown の場合、結果はまだ検証されていません。
• results.in_asymptotic_range が true であり、GCI または外挿値を引用する前に、notes の項目に漸近前（次数の変動が50％を超える）、負／非正の次数、または誤差差分がゼロであるとの報告がないことを確認しました。
• gci_calculator.py の細分化比 r21／r32 が1.3以上（丸め誤差ノイズの下限）であり、Oscillatory convergence detected エラーが発生していないことを確認しました。
• 報告する離散化不確かさとして results.gci_fine（使用した安全係数も記録：次数を検証済みのグリッドが3つ以上の場合は1.25、次数を仮定したグリッドが2つの場合は3.0）を記録し、最良推定値として results.extrapolated_value を記録しました。
• 細分化比が一定の場合、1.0に近い asymptotic_ratio を漸近範囲の証明として扱いませんでした（これは f1/f2 に縮退します）。観測次数と予想次数の比較によって漸近範囲を確認し、gci_calculator.py の notes にある注意事項を読みました。

よくある落とし穴と誤った自己正当化

ありがちな近道	誤りである理由／すべきこと
「2つのグリッドで値がほぼ一致しているので、収束しています」	2段階だけでは観測次数を推定できません。h_refinement.py／dt_refinement.py を3段階以上で実行してください。2グリッドのリチャードソン実行では次数を 仮定 するため、安全係数には1.25ではなく3.0が必要です。
「漸近比が約1.0なので、漸近範囲に入っています」	細分化比が一定の場合、AR = GCI_coarse/(r^p*GCI_fine) は f1/f2 に縮約され、最も細かい2つの QoI 値の差だけを測定します。観測次数と予想次数を比較し（さらに4つ以上のグリッドで一貫性を確認し）、漸近範囲を検証してください。
「GCI_fine が非常に小さいので、解はグリッドに依存しません」	GCI がほぼゼロでも、QoI が細分化の影響を受けにくい場合や、差分が丸め誤差ノイズの範囲内（比が1.3未満）である場合があります。まず in_asymptotic_range が true であり、細分化比が1.3以上であることを確認してください。
「観測次数が予想より高いので、さらに良い結果です」	次数が形式次数を大幅に上回る場合、通常は超収束または誤差の相殺を示しており、精度がさらに高いことを意味するわけではありません。これを警告として扱い（10％を超えて外れると convergence_assessment は FAIL になります）、より多くのグリッドレベルで検証してください。
「スクリプトが外挿値を出力したので、それを使います」	観測次数が非正の場合、解は発散しており、h_refinement.py は発散を示す注記とともに richardson_extrapolated_value = null を返します。その場合、外挿値や GCI を引用しないでください。
「陰解法／安定なソルバーなので、解析にはどのタイムステップでも問題ありません」	時間的な 安定性 は時間的な 精度 を意味しません。dt_refinement.py でも、スキームの次数を求めるには体系的に縮小したタイムステップが3つ以上必要です。粗すぎる dt の系列では、漸近前の状態が続きます。

セキュリティ

入力検証

• すべての数値パラメーター（spacings、timesteps、values、expected-order、order）が有限の正の数であることを検証します
• カンマ区切りの値リストは長さが一致する必要があり（spacings と values の長さが同じであること）、要素数は10,000件までに制限されます
• GCI 計算ツールでは細分化レベルを正確に3つにする必要があります。リチャードソン外挿には最低2つ必要です
• 安全係数は1.0以上の有限数であることを検証します（Roache では Fs は {1.25, 3.0} を使用します）

ファイルアクセス

• スクリプトは外部ファイルを読み込みません。すべての入力は CLI 引数で指定します
• スクリプトが書き込むのは標準出力（JSON 出力）のみです。エージェントが明示的に Write ツールを使用しない限り、ファイルは作成されません

ツールの制限

• Bash：4つの Python 解析スクリプト（h_refinement.py、dt_refinement.py、richardson_extrapolation.py、gci_calculator.py）を、明示的な引数リストを指定して実行するために使用します
• Read：スクリプトのソースと参考ドキュメントを確認するために使用します

安全対策

• eval()、exec()、動的コード生成は使用しません
• すべてのサブプロセス呼び出しで明示的な引数リストを使用します（shell=True は使用しません）
• スクリプトでは Python 標準ライブラリ（math モジュール）のみを使用します。信頼できないデータの pickle 読み込みやデシリアライズは行いません
• ツールを Bash と Read のみに限定することで、エージェントによるファイルシステムの変更を制限します

参考資料

• references/convergence_theory.md - 形式的な収束次数、両対数解析、漸近領域
• references/gci_guidelines.md - Roache の GCI 手法、ASME V&V 20、安全係数