MITC4シェル: 線形ソルバー非収束時に有限回転の細分ループが暴走することで、OpenMP実行で実行時間が数百倍に増大する？

概要

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

原因は、線形ソルバーが収束しなかった場合にその解をそのまま変位増分として適用し、有限回転の更新ルーチンが持つ「回転増分を0.05rad刻みに細分するループ」の反復回数が制限なく増大することにある。OpenMPは、この状態を引き起こす丸め誤差の非決定性をもたらす引き金となっている。

現象

CIでの観測

testステージの所要時間（同一pipeline、秒）:

ジョブ 所要時間
ubuntu2204/openmp/openmp 144.8
ubuntu2404/openmp/openmp 719.7
ubuntu2404/hybrid/openmp 678.0
ubuntu2204/hybrid/openmp 2178.0
ubuntu2204/serial/serial 136.3
ubuntu2404/serial/serial 137.3

同じラベルでもジョブによって10倍以上ばらつく。

ローカルでの再現

FrontISTR master 8320e93c（Version 5.9）を macOS 26.5 / arm64 / GCC 16.1.0 / Open MPI 5.0.9、-DCMAKE_BUILD_TYPE=Release でビルドし、tests/test.sh 経由で lib/static_lib_shell_mitc4（4ケース合計）を反復実行した結果:

構成 プロセス×スレッド 試行数 最小 最大
serial 1×1 15 4.18 s 4.53 s
mpi 2×1 15 5.01 s 6.58 s
openmp 1×2 15 5.02 s 10.84 s
hybrid 2×2 13 7.05 s 628.40 s

スレッド数が2以上の構成（openmp, hybrid）でのみ外れ値が出る。MPIの有無は関係しない。

単一ケース K741ULSTRONGAUTOINC を2スレッドで84回反復した場合、平均4.68秒に対し最大299.42秒、37%の実行で線形ソルバーが反復上限に到達した。1スレッドでは84回すべてが同一の反復回数（CG総反復13481回）となり、実行時間も一定である。

原因

以下の3段の連鎖による。

(1) OpenMPが実行ごとの丸め誤差の差を生む

全体剛性・内力の組立が要素ループの !$omp atomic による加算で行われているため、加算順序がスレッドのスケジューリングに依存する。

  • hecmw1/src/solver/matrix/hecmw_mat_ass.f90 （D/AL/AU への加算、97行・118行・129行ほか）
  • fistr1/src/analysis/static/fstr_Update.f90:270 （QFORCE への加算）
  • fistr1/src/analysis/static/fstr_CreateMatrix.f90:183

結果として、同一バイナリ・同一入力でも実行ごとに組み上がる行列・残差ベクトルが最下位ビットで変化する。1スレッドでは加算順序が一意に定まるため完全に決定的となる。

なお線形ソルバー側にスレッド依存の数値差はない。内積 hecmw_innerProduct_R（hecmw1/src/solver/communication/hecmw_solver_misc.F90:46）はOpenMP化されておらず逐次であり、行列ベクトル積（hecmw_solver_las_66.F90）とDIAG前処理（hecmw_precond_DIAG_66.f90）はいずれも行ごとに独立で、総和順序がスレッド数に依存しない。

(2) 接線剛性が悪条件で、CGが発散することがある

対象モデルは33節点・20要素（198自由度）と極小だが、MITC4シェルの接線剛性は条件数が大きい（ESTCOND によるDIAG前処理後の推定条件数は約5.6×10⁴）。CGは198次元の系に対し300〜10000反復を要しており、次元を大きく超えている（丸め誤差により直交性が失われている）。

このため、(1)の最下位ビットの差によって、収束する実行と発散する実行に分かれる。発散した場合、相対残差が4.92×10⁸まで増大したまま反復上限10000で打ち切られる。

(3) 非収束の解がそのまま適用され、細分ループが暴走する

fistr1/src/analysis/static/fstr_solve_NonLinear.f90:136-143 は、solve_LINEQ の後に線形ソルバーの収束フラグ（hecmw_mat_get_flag_converged / hecmw_mat_get_flag_diverged）を確認せず、解 hecMAT%X をそのまま fstr_apply_solution_increment に渡す。

fistr1/src/common/fstr_FiniteRotationKinematics.f90:247:

nsub = max(1, ceiling(theta_norm/5.0D-2))
theta_step = theta_inc/dble(nsub)
do isub = 1, nsub
  ...
  call ShellRotationVectorToMatrix(theta_phys, rot_inc)
  triad_trial = matmul(rot_inc, triad_base)
  call ShellOrthonormalizeTriad(triad_trial, triad_base)
  ...
end do

回転増分を0.05rad刻みに細分するループで、反復回数 nsub に上限がない。発散した線形解由来の巨大な theta_inc が渡されると、1節点の更新だけで数百秒を要する。

実測（299.29秒を要した実行のログとスタックサンプル）:

項目 値
総実行時間 299.29 s
全21回の線形ソルバー時間の合計 1.128 s
線形ソルバー外の時間 298.2 s（99.6%）
サブステップ数 4（正常時と同じ）

実行中のスタックサンプルは全サンプルが fstr_apply_solution_increment → ShellUpdateTriadWithIncrement → ShellRotationVectorToMatrix → sin/cos にあった。

その後はAUTOINCのカットバックにより最終的には収束するため、結果はリファレンスと一致しテストはPASSする。

なお nsub は既定の4バイト整数であり、theta_norm が約1.07×10⁸を超えると ceiling の結果が32bit整数の範囲を超える。

再現方法

cd tests
./test.sh -f <build>/fistr1/fistr1 -e <build>/hecmw1/tools/hecmw_part1 \
          -r <build>/hecmw1/tools/rmerge -d lib/static_lib_shell_mitc4 -p 1 -t 2

1回で再現するとは限らないため、10〜60回程度繰り返すこと。K741ULSTRONGAUTOINC.msh の所要時間が通常1〜2秒のところ、数十〜数百秒になる実行が現れる。-t 1 では発生しない。

影響範囲

  • 有限回転運動学を使うMITC4シェル要素（!ELEMENT, TYPE=741、4節点、弾性材料、TL/UL）の非線形解析
  • OpenMPスレッド数2以上での実行（openmp・hybrid の両ラベル）
  • テストは失敗しないため、CIの所要時間の増大としてのみ現れる
  • 実運用のモデルでも、線形ソルバーが収束しない場面で同様の停滞が起こりうる

補足: 検証して否定した仮説

  • 前処理とOpenMPの相性: 当初、PRECOND = 3 がSSORで、スレッド数2以上でRCM+マルチカラー並べ替えに切り替わることを疑った。しかし hecmw1/src/solver/precond/hecmw_precond.f90:37-41 のとおり case(3) はDIAG（ブロック対角）であり、マルチカラーSSORの経路には入らない（実行ログのバナーも ### 6x6 BLOCK CG, DIAG, 1）。DIAGは行ごとに独立であり、スレッド数に依存する数値差を生じない。
  • ソルバー設定の調整: PRECOND = 1（SSOR）に変更すると2スレッドで190〜212秒/回、1スレッドでも10分超と、かえって悪化する。反復上限を10000から500に下げても停滞は発生する（12回中1回）。設定の問題ではない。

環境

  • FrontISTR master 8320e93c5b8394d4854530e28665d58a0a0747a5（Version 5.9）
  • macOS 26.5 / arm64、GCC 16.1.0、Open MPI 5.0.9、CMake 4.4.2、-DCMAKE_BUILD_TYPE=Release
  • CI側の観測は ubuntu2204 / ubuntu2404 イメージ