Skip the increment when the linear solver did not converge (for MITC4 test)

Kanato Konorequested to merge
chikasuiro/FrontISTR:bugfix/shell-mitc4-rotation-subdivision into master

線形ソルバーが収束しなかった場合に、その解をそのまま変位増分として適用してしまう問題を修正します。

Close #813 (closed)

背景

tests/lib/static_lib_shell_mitc4 が、OpenMPスレッド並列（2スレッド以上）の実行でのみ、実行ごとに数百倍の実行時間を要することがありました。結果はリファレンスと一致するためテストはPASSし、CIのtestジョブ所要時間の増大としてのみ現れていました。

原因は次の3段の連鎖です。

  1. 全体剛性・内力の組立が要素ループ内の !$omp atomic による加算であるため、加算順序がスレッドのスケジューリングに依存し、実行ごとに丸め誤差が変わる（1スレッドでは決定的）
  2. MITC4シェルの接線剛性は条件数が大きく（DIAG前処理後の推定条件数 約5.6×10⁴）、198自由度の系にCGが300〜10000反復を要するため、1の最下位ビットの差で収束する実行と発散する実行に分かれる
  3. fstr_Newton が線形ソルバーの収束状態を確認せずに解を増分として適用するため、発散した解（相対残差 4.92×10⁸ を観測）由来の巨大な回転増分が ShellUpdateTriadWithIncrement に渡り、nsub = max(1, ceiling(theta_norm/5.0D-2)) という上限のない細分ループが暴走する

299.29秒を要した実行では、全21回の線形ソルバー時間の合計が1.128秒であり、99.6%が線形ソルバー外で消費されていました。実行中のスタックサンプルは全サンプルが fstr_apply_solution_increment → ShellUpdateTriadWithIncrement → ShellRotationVectorToMatrix → sin/cos にありました。

変更内容

1. fstr_Newton に線形ソルバーの収束確認を追加（fistr1/src/analysis/static/fstr_solve_NonLinear.f90）

solve_LINEQ の直後に収束状態を確認し、非収束であれば増分を適用せず、Newton-Raphson反復の失敗として扱います。同一ファイルの fstr_Newton_contactSLag に既存の同型のガードと、メッセージ書式・NRstat_i(knstDRESN) = 4（既存の理由コード MatSolveError）・CutBack_stat の加算まで揃えています。solve_LINEQ はステータスを返さないため、判定には hecmw_mat_get_flag_converged を用いています。

2. hecmw_solve の右辺ゼロ時にステータスフラグを設定（hecmw1/src/solver/main/hecmw_solver.f90）

変更1を入れた時点で、既存テスト tests/analysis/static/amplitude と tests/analysis/static/amplitude_total がserialで失敗しました。調査の結果、これは変更1が既存の潜在不具合を露見させたものでした。

hecmw_solve は右辺ベクトルのノルムが0の場合（hecmw_solve_check_zerorhs）に X = 0 を設定して即座にreturnしますが、flag_converged と flag_diverged を更新しません。このため両フラグには直前の求解の値が残り、自明解が得られている場合でも「非収束」と誤判定されます。X = 0 は右辺が0のときの厳密解であるため、この分岐で flag_converged を1、flag_diverged を0に設定するようにしました。

この潜在不具合は本MRの変更に起因するものではありません。右辺ゼロ時の早期returnは e289281d（2018-10-02）で導入され、flag_converged / flag_diverged 自体は f42fa4bf（2026-05-04）で追加されたため、フラグが更新されない状態は2026-05-04以降ずっと存在していました。amplitude系テストを整備した f736b63b（2026-08-11）のマージ直前のmaster（0b577bfe）を別途ビルドして tests/analysis/static/amplitude を実行したところ、右辺ゼロ経路を示す #### HEC-MW-SOLVER-W-2002: ZERO RHS norm がマージ後と同じ4回出力されており、この経路はそれ以前から通っていたことを確認しています。

なお hecmw1/src/solver/contact/solve_LINEQ_contact_elim.f90 は以前から hecmw_mat_get_flag_diverged を istat として返しており、接触解析ではこの古い値を参照しうる状態にありました（実害の有無は未確認です）。

効果

同一マシン・他プロセスなしの同一条件で、単一ケース K741ULSTRONGAUTOINC を2スレッドで反復実行した結果です。

試行数 平均 最大 線形ソルバー反復上限到達
修正前 60 6.09 s 312.75 s 28%
修正後 100 0.58 s 1.52 s 0%

tests/test.sh で tests/lib/static_lib_shell_mitc4（4ケース合計）を各15回実行した結果です。

構成 修正前 最小〜最大 修正後 最小〜最大
serial (1proc×1thread) 4.18〜4.53 s 4.00〜4.19 s
mpi (2proc×1thread) 5.01〜6.58 s 4.94〜5.48 s
openmp (1proc×2thread) 5.02〜10.84 s 4.95〜5.10 s
hybrid (2proc×2thread) 7.05〜628.40 s 7.00〜7.44 s

スレッド数による所要時間の外れ値が解消しています。

テスト結果

ローカル（macOS 26.5 / arm64、GCC 16.1.0、Open MPI 5.0.9、Release）でのctest結果です。

ラベル 結果
serial 78件中78件PASS
openmp 78件中78件PASS
hybrid 78件中77件PASS

hybridの失敗1件 test_hybrid_analysis/dynamic/autoinc_cutback は、本修正前のバイナリでも2回中2回失敗し、serialでは修正後も通過するため、本修正とは無関係な既存の失敗です。同テストがhybridでのみ失敗する点は、本件と同じOpenMP由来の非決定性が関与している可能性がありますが未確認です。

本MRの範囲外とした事項

  • fstr_Newton_contactALag の同種の欠落: 同一ファイルの3つのNewton-Raphson反復ルーチンのうち、拡張Lagrange法の経路だけが求解ステータスを確認していません（istat を受け取りながら参照していません）。別issueとして切り出します。
  • 組立の !$omp atomic による非決定性そのもの: 本件の引き金ですが、決定化しても1スレッドとは別の丸めに固定されるだけで「発散する側」に固定される可能性が残るため、本MRの修正なしでは対策になりません。
  • 回帰テストの追加: 本不具合は確率的に発現するため、決定的な再現テストを追加できていません。tests/lib/static_lib_shell_mitc4 の所要時間自体が指標になるため、ctestの TIMEOUT プロパティによる検出の余地があります。

レビュー時にご確認いただきたい点

  • flag_converged は求解後に再計算した相対残差が !SOLVER の判定値を下回る場合にのみ1になるため、反復上限到達だけでなく、判定値をわずかに上回った場合や直接法（MUMPS / MKL PARDISO）でも本ガードが発動します。INC_TYPE = AUTO 未指定時は解析停止となるため、CIの with_mumps / with_mkl を含むジョブでの確認をお願いします（ローカルではMUMPS/MKLなしのビルドで検証しています）。
  • 追加したガードのreturnは hecmw_mpc_mat_finalize を経由しないため、MPCを含むモデル（MPC_METHOD=3）ではカットバックのたびに hecMESHmpc / hecMATmpc がリークします。同じ問題は同ルーチンの既存の早期return（iterStatus == kitrDiverged 時）にもありますが、本MRで解放処理を集約すべきかご意見をいただきたいです。

Merge request reports