ブシネスク問題（半無限弾性体の点荷重） — ソルバー別実装詳細

軸対称CAX8Rモデルの場合の要点はこうだ。

• *NODE で $r \geq 0$ の節点を定義。$r = 0$ の節点は自動的に対称軸として扱われる
• *ELEMENT, TYPE=CAX8R で8節点低減積分の軸対称四角形要素を使用
• *BOUNDARY で底面を YSYMM（$u_z = 0$）、側面に無限要素を配置するか $u_r = 0$
• *CLOAD で対称軸上の表面節点に $P/2\pi$（軸対称なので周方向積分分の補正は不要。Abaqusが自動で処理）

注意: 軸対称モデルのCLOADは「単位ラジアンあたりの荷重」ではなく「全周分の荷重」として入力する。これを間違えると$2\pi$倍ずれる。

OUTPUT, HISTORY, VARIABLE=PRESELECT で特定節点の変位を時刻歴出力。応力はEL PRINT で要素積分点値を取得するか、Python ODB APIで座標指定のフィールド出力を抽出する。probeValues メソッドで任意点の補間値を取得できるから、理論値との自動比較に便利だ。

CTRIAX6は6節点三角形で、四角形にする場合はTRAX3やTRAX6を使うが、Abaqusほど使いやすくはない。実務的には3D 1/4対称モデル（CHEXA-20）で解くことが多い。

SOL 101で解き、GPSTRESS(SORT1,REAL)で応力出力を取得する。特定節点での応力を抽出するにはDMAP ALTER か pyNastran の OP2 読み込みが便利だ。

Nastranの応力テンソルの出力順序は $\sigma_{xx}, \sigma_{yy}, \sigma_{zz}, \sigma_{xy}, \sigma_{yz}, \sigma_{xz}$ で、Abaqusは $\sigma_{11}, \sigma_{22}, \sigma_{33}, \sigma_{12}, \sigma_{13}, \sigma_{23}$ だ。特にせん断成分の順序が異なる（Nastranは $xy, yz, xz$、Abaqusは $12, 13, 23$）から、比較時に取り違えないように注意すること。

CalculiXはAbaqus互換の.inpファイルがほぼそのまま使える。C3D20要素で3D 1/4モデルを解く。注意点として、CalculiXのC3D20の積分点配置がAbaqusと微妙に異なることがある。結果は有効数字3〜4桁で一致するが、それ以上の精度が必要なら積分スキームの詳細をソースコード（src/e_c3d.f）で確認すべきだ。

Code_Asterは MODELISATION='AXIS' + QUAD8 で軸対称解析。SALOMEのSMESH モジュールでメッシュを生成し、MED形式で出力する。コマンドファイルのAFFE_MODELEとAFFE_MATERIAUで物理設定を定義する。

ソースコードが公開されているため、要素定式化のアルゴリズムを第三者が直接検証できることだ。実際、CalculiXのC3D20の剛性マトリクス計算はe_c3d.f内のGauss積分で明示されている。商用ソルバーではこの透明性がない。学術論文でOSSベースの検証を行い、その結果を商用ソルバーの独立検証として引用するのは確立された手法だ。