ヘルツ接触理論(球−平面)
理論と物理
概要
先生、Hertz接触はFEAの接触解析の検証で一番使われる問題と聞きましたが、なぜですか?
2つの弾性体の接触で、接触半径、最大面圧、接触変位の全てに閉形式の解析解が存在するからだ。1882年にHertzが導出した。球と平面の接触が最も基本的で、非線形接触解析のCode Verificationの入り口になる。
非線形なのに解析解があるんですか?
Hertz理論の巧妙さはそこにある。接触面積が荷重とともに変化する幾何学的非線形を含むが、弾性半空間のBoussinesq解の重ね合わせと幾何学的適合条件を連立することで閉形式が得られる。ただし前提条件がある。接触面が平面にくらべて十分小さいこと、弾性限度内であること、摩擦がないことだ。
支配方程式
具体的な式を教えてください。
半径 $R$ の弾性球が弾性平面に荷重 $P$ で押し付けられる場合。等価弾性係数 $E^ = [(1-\nu_1^2)/E_1 + (1-\nu_2^2)/E_2]^{-1}$、等価半径 $R^ = R$(平面は $R_2 = \infty$)として
接触半径: $a = \left(\frac{3PR^}{4E^}\right)^{1/3}$
最大面圧: $p_0 = \frac{2aE^}{\pi R^} = \frac{1}{\pi}\left(\frac{6PE^{2}}{R^{2}}\right)^{1/3}$
近づき量: $\delta = \frac{a^2}{R^} = \left(\frac{9P^2}{16R^E^{*2}}\right)^{1/3}$
面圧分布はどうなりますか?
接触面内の面圧は半楕円分布だ。
中心で最大 $p_0$、接触縁 $r = a$ でゼロ。この分布はBoussinesq解の重ね合わせで再現でき、接触面下の応力場も厳密に計算できる。最大von Mises応力は表面直下 $z \approx 0.48a$ で生じ、転がり疲労の起点になる。
ベンチマーク数値例
具体的な数値での検証例をお願いします。
鋼球($R = 10$ mm、$E = 200$ GPa、$\nu = 0.3$)vs 鋼平面(同材料)。$P = 100$ N。
$E^* = 200/(2 \times (1-0.09)) = 109.9$ GPa
$a = (3 \times 100 \times 0.01 / (4 \times 109.9 \times 10^9))^{1/3} = 0.0727$ mm
$p_0 = 3P/(2\pi a^2) = 9013$ MPa
$\delta = a^2/R = 0.000529$ mm
| 物理量 | 理論値 | Abaqus (C3D20R) | 誤差 |
|---|---|---|---|
| 接触半径 $a$ [mm] | 0.0727 | 0.0731 | 0.55% |
| 最大面圧 $p_0$ [MPa] | 9013 | 8940 | 0.81% |
| 近づき量 $\delta$ [mm] | 0.000529 | 0.000527 | 0.38% |
1%以内の精度ですね。メッシュの細かさはどの程度ですか?
接触面内に20要素以上を配置した結果だ。接触半径 $a$ が微小だから、接触部の要素サイズは $a/10 \approx 7$ μm にする必要がある。遠方は粗くできるからバイアスメッシュが必須だ。
各項の物理的意味
- 保存量の時間変化項:対象とする物理量の時間的変化率を表す。定常問題では零となる。【イメージ】浴槽にお湯を張るとき、水位が時間と共に上がる——この「時間あたりの変化速度」が時間変化項。バルブを閉じて水位が一定になった状態が「定常」であり、時間変化項はゼロ。
- フラックス項(流束項):物理量の空間的な輸送・拡散を記述する。対流と拡散の2種類に大別される。【イメージ】対流は「川の流れがボートを運ぶ」ように流れに乗って物が運ばれること。拡散は「インクが静止した水中で自然に広がる」ように濃度差で物が移動すること。この2つの輸送メカニズムの競合が多くの物理現象を支配する。
- ソース項(生成・消滅項):物理量の局所的な生成または消滅を表す外力・反応項。【イメージ】部屋の中でヒーターをつけると、その場所に熱エネルギーが「生成」される。化学反応で燃料が消費されると質量が「消滅」する。外部から系に注入される物理量を表す項。
仮定条件と適用限界
- 連続体仮定が成立する空間スケールであること
- 材料・流体の構成則(応力-歪み関係、ニュートン流体則等)が適用範囲内であること
- 境界条件が物理的に妥当かつ数学的に適切に定義されていること
次元解析と単位系
| 変数 | SI単位 | 注意点・換算メモ |
|---|---|---|
| 代表長さ $L$ | m | CADモデルの単位系と一致させること |
| 代表時間 $t$ | s | 過渡解析の時間刻みはCFL条件・物理的時定数を考慮 |
数値解法と実装
接触アルゴリズムの選択
FEAの接触アルゴリズムにはどんな種類がありますか?
主に3種類だ。
| 手法 | 原理 | 精度 | 計算コスト |
|---|---|---|---|
| ペナルティ法 | 侵入量に比例した反力 | ペナルティ定数依存 | 低 |
| ラグランジュ乗数法 | 侵入量=0を厳密に課す | 高精度 | 高(自由度増加) |
| 拡張ラグランジュ法 | ペナルティ+ラグランジュの混合 | 高精度 | 中 |
Hertz検証にはラグランジュ乗数法か拡張ラグランジュ法を推奨する。ペナルティ法だとペナルティ係数の値で面圧が変わるから、検証には適さない。
各ソルバーでのデフォルト設定はどうなっていますか?
Abaqusはデフォルトで拡張ラグランジュ法(SURFACE INTERACTION + CONTACT PAIR)。NastranはSOL 101のGLUE/SLIDE接触でペナルティベース。Ansysはデフォルトで拡張ラグランジュ。接触問題のV&Vでは使用したアルゴリズムを明記することが必須だ。
メッシュ設計
接触問題のメッシュ設計で最も重要なことは何ですか?
接触面の要素サイズが接触半径 $a$ の1/10以下であること。粗いメッシュだと接触面が階段状になり、面圧分布が不正確になる。
推奨構成:
- 接触面近傍: $h_{elem} = a/10$。HEX20の構造メッシュ
- 接触面から離れた領域: バイアスで粗くする。比率2.0程度
- 球の1/4対称モデル(2面の対称条件)
- 平面側は半径 $10a$ 以上のモデルサイズ
接触面にはどちらの面をmaster/slaveにすべきですか?
Abaqusの慣例では凸面(球)をslave、平面をmasterにする。slaveの方がメッシュを細かくすべきだ。NastranやAnsysでも同様の考え方で、剛性が高い面をmasterにする。両面のメッシュ密度が大きく異なる場合は面圧の精度に影響するから、接触面では両側のメッシュ密度を揃えるのが望ましい。
収束のポイント
接触解析が収束しないことがよくありますが、対処法は?
Hertz問題のような滑らかな接触は通常問題なく収束する。収束しない場合:
1. 初期接触のギャップを確認。ギャップが大きすぎると初期ステップで接触が検出されない
2. 荷重を複数ステップに分割。NLGEOM=ON で幾何非線形を有効にする
3. 接触の安定化(Abaqusの*CONTACT CONTROLS, STABILIZE)を一時的に有効にする
4. スムーシングオプション(Abaqusの*SURFACE SMOOTHING)で接触面の不連続を軽減
NLGEOMが必要なのはなぜですか?
Hertz接触では球と平面の近づきに伴って接触面積が変化する。これは幾何学的非線形だからNLGEOM=ONが必須。ONにしないとソルバーは初期のギャップ状態に基づいた線形近似をしてしまい、正しい接触面積が得られない。
低次要素
計算コストが低く実装が簡単だが、精度は限定的。粗いメッシュでは大きな誤差が生じる可能性がある。
高次要素
同一メッシュでより高い精度を達成。計算コストは増加するが、必要な要素数は少なくなる場合が多い。
ニュートン・ラフソン法
非線形問題の標準的手法。収束半径内で2次収束。$||R|| < \epsilon$ で収束判定。
時間積分
実践ガイド
検証手順
Hertz接触のV&V検証を系統的に行うにはどうすればいいですか?
段階的に進める。
1. Boussinesq問題の検証(前段): 集中荷重下の半空間の応力場が理論値に一致することを確認
2. Hertz接触(同材・球-平面): 最も基本的な設定で接触半径、面圧、近づき量を検証
3. Hertz接触(異材・球-球): 等価弾性係数と等価半径の計算が正しいか確認
4. 摩擦付きHertz接触: Mindlin理論との比較。摩擦力分布の検証
5. 弾塑性接触: Hardnessとの関係(Meyer則)。理論解はないのでクロスチェック
各段階で何を評価指標にすべきですか?
ステップ1: 応力の理論値との一致。ステップ2-3: 接触半径(5%以内)、最大面圧(5%以内)、荷重-変位関係の$P \propto \delta^{3/2}$の指数の一致。ステップ4: すべり開始荷重。ステップ5: 残留変形と面圧分布のソルバー間一致。
転がり疲労への応用
Hertz接触が転がり軸受の設計にどう関係するんですか?
軸受の接触応力分布はHertz理論で算出し、接触面下の交番せん断応力が転がり疲労寿命を支配する。ISO 281の動定格荷重はHertz理論に基づくLundberg-Palmgren理論で計算される。
FEAでHertz接触を解き、最大せん断応力 $\tau_{max}$ の位置と値を検証することで、軸受設計ツールの妥当性を確認できる。$\tau_{max} = 0.31 p_0$($z = 0.48a$)が理論値だ。
楕円接触(線接触に近い場合)はどう扱いますか?
Hertzの一般理論は楕円接触にも対応しており、楕円積分を含む式で記述される。実用的にはHamrock-Dowsonの近似式が使える。ローラー軸受やギヤの歯面接触がこのケースに相当する。FEAでは3Dモデルで直接解くが、理論値との比較で精度を担保する。
結果の評価方法
接触面積の正確な評価方法を教えてください。
FEAでは接触状態(open/closed/sliding)のフラグで接触領域を判定する。Abaqusでは COPEN(接触開口量)がゼロ以下の節点が接触状態にある。この節点群の外縁が接触半径に対応する。
より正確には、面圧が正の領域の面積を積分する。CPRESS > 0 の要素面積の合計が接触面積 $\pi a^2$ と一致するか確認する。メッシュが粗いと接触縁の階段状離散化で誤差が出るから、細かいメッシュが必要だ。
解析フローのたとえ
解析フローは「科学実験」に似ている。仮説(解析モデル)を立て、実験(計算実行)し、結果を検証し、仮説を修正する——このPDCAサイクルが品質の高い解析を生む。
初心者が陥りやすい落とし穴
最もよくある失敗は「結果の検証を怠る」こと。美しいコンター図が得られても、それが物理的に正しいとは限らない。必ず理論解、実験データ、またはベンチマーク問題との比較を行うこと。
境界条件の考え方
境界条件は「実験の治具」に相当する。治具の設計が不適切であれば実験結果が無意味になるように、CAEでも境界条件が現実を正しく表現しているかが最も重要。
ソフトウェア比較
Abaqusでの実装
AbaqusでHertz接触を設定する具体的な方法を教えてください。
CONTACT PAIR で球面をslave surface、平面をmaster surfaceに指定する。摩擦なしならSURFACE INTERACTION + *SURFACE BEHAVIOR, PRESSURE-OVERCLOSURE=HARD を設定。
NLGEOM=YES のSTATICステップで荷重をRampで印加。自動増分(CONTROLS, PARAMETERS=FIELD で接触の安定化制御)を使い、最終荷重で接触が安定するまで反復。
出力: CPRESS(面圧)、COPEN(ギャップ)、U(変位)をフィールド出力。History出力でRF(反力)とU(変位)を記録し、$P-\delta$曲線を作成。
General Contact と Contact Pair のどちらを使うべきですか?
V&Vでは Contact Pair を推奨する。General Contact はアルゴリズムが自動選択されて透明性が低い。Contact Pair なら Surface-to-Surface discretization + Augmented Lagrange を明示的に指定でき、レポートに設定を明記しやすい。
Nastranでの実装
Nastranでの接触設定は?
SOL 400(非線形)+ BCTABLE/BCONPRG で接触定義。BCBODY1でmaster、BCBODY2でslaveを指定。BCTPARM で接触アルゴリズムのパラメータを制御。
Nastranの接触はAbaqusに比べて設定項目が多く複雑だが、ペナルティ法とラグランジュ法の切り替えがBCTPARMで明示的に制御できる。
AbaqusとNastranで結果は一致しますか?
同等のメッシュ密度・接触アルゴリズムで、接触半径と最大面圧は2〜3%以内で一致する。荷重-変位曲線の$P \propto \delta^{3/2}$の指数は両者とも理論値の1.5に非常に近い値(1.48〜1.52)を示す。
オープンソースでの検証
CalculiXやCode_Asterでも接触解析はできますか?
CalculiXは*CONTACT PAIRをAbaqus互換で実装しており、HEX20要素のface-to-face接触が可能。Hertz問題では商用ソルバーと遜色ない結果が得られる。
Code_Asterは DEFI_CONTACT + FORMULATION='DISCRETE' または 'CONTINUE' で接触定義。CONTINUEが拡張ラグランジュに相当し、Hertz問題に適している。
OSSの接触解析の信頼性はどう評価しますか?
Hertz問題が最良のリトマス試験だ。OSSの結果が理論値に5%以内で一致し、メッシュ収束が確認できれば、その接触アルゴリズムの基本的な信頼性は保証される。その上で商用ソルバーとのクロスチェックを行えば、独立した二重検証になる。
選定で最も重要な3つの問い
- 「何を解くか」:ヘルツ接触理論(球−平面)に必要な物理モデル・要素タイプが対応しているか。例えば、流体ではLES対応の有無、構造では接触・大変形の対応能力が差になる。
- 「誰が使うか」:初心者チームならGUIが充実したツール、経験者ならスクリプト駆動の柔軟なツールが適する。自動車のAT車(GUI)とMT車(スクリプト)の違いに似ている。
- 「どこまで拡張するか」:将来の解析規模拡大(HPC対応)、他部門への展開、他ツールとの連携を見据えた選択が長期的なコスト削減につながる。
先端技術
弾塑性接触
荷重が大きくて塑性変形が起きる場合はどうなりますか?
最大von Mises応力が降伏応力に達すると塑性変形が始まる。降伏開始荷重は $P_y = (\pi \sigma_y)^3 R^{2} / (6E^{2}) \times C_y^3$ で、$C_y \approx 1.6$ は理論的に求められている(Tabor, 1951)。
塑性変形が進むと面圧分布がHertzの半楕円分布から一様分布に近づき、接触半径は弾性予測より大きくなる。この領域では解析解がないから、FEAのクロスバリデーションか実験との比較がValidationの手段になる。
FEAで弾塑性接触を解くときの注意点は?
材料モデルの選択が結果を大きく左右する。等方硬化と移動硬化で除荷後の残留変形が変わる。J2塑性(von Mises)が標準だが、接触面直下の3軸応力状態でのモデルの妥当性を考慮する必要がある。メッシュは弾性ケースの2〜3倍の密度が必要で、塑性域の要素サイズは接触半径の1/20以下にする。
動的接触(衝突)
球が平面に衝突する場合はどうですか?
Hertz理論を準静的に適用すると、衝突時間$T$と最大変形$\delta_{max}$が予測できる。
$m$ は球の質量、$v$ は衝突速度。陽解法(LS-DYNA、Abaqus/Explicit)の検証に使える。反発係数 $e = 1$(完全弾性衝突)の再現が検証指標だ。
CFL条件はどう設定しますか?
接触部の要素サイズが小さいとCFL条件 $\Delta t \leq h_{min}/c_L$($c_L$は縦波速度)で時間刻みが極めて小さくなる。質量スケーリング(人工的に密度を増加)で$\Delta t$を大きくする方法があるが、衝突解析では慣性が重要だから適用には注意が必要。
粗さを考慮した接触
実際の表面は粗さがありますが、理論との乖離はどうなりますか?
Greenwood-Williamson(1966)モデルが代表的だ。表面粗さを多数のアスペリティ(微小突起)のHertz接触として統計的に扱う。粗さの影響で見かけの接触面積が理論値より大幅に小さくなり、実接触面は名目接触面積の1〜10%程度になる。
FEAではフラクタル表面の直接モデル化や、粗さの統計的パラメータを入力する表面相互作用モデル(Abaqusの*SURFACE BEHAVIOR with roughnessなど)が利用可能だ。ただし計算コストが膨大になるため、マルチスケール手法やROMが研究されている。
トラブルシューティング
収束しない場合
接触解析が収束しないときのデバッグ手順を教えてください。
以下の順で確認する。
1. 接触の初期状態: 球と平面が初期にわずかに重なっている(overclosure)か、離れすぎていないか確認。初期ギャップが大きいと接触検出ができない
2. 荷重ステップの分割: 最終荷重の1/10から始めて段階的に増やす。AbaqusではStep内のInitial/Maximum incrementを小さく設定
3. 接触安定化: *CONTACT CONTROLS, STABILIZE で微小な粘性ダンピングを追加し収束を助ける。最終結果で安定化の影響がないことを確認
4. 要素タイプ: HEX8は接触面でチャタリング(接触/非接触の振動)が起きやすい。HEX20に切り替える
チャタリングとは何ですか?
接触縁の節点が反復計算中にopen/closedを繰り返す現象だ。これが収束を妨げる最大の要因。surface-to-surface離散化(node-to-surfaceではなく)を使うと改善される。Abaqusでは*CONTACT PAIR, INTERACTION=... の TYPE=SURFACE TO SURFACE を指定する。
面圧が理論値と合わない場合
面圧分布がHertzの半楕円にならないのですが。
いくつかの原因が考えられる。
- メッシュ不足: 接触面内に最低10要素が必要。5要素以下では半楕円分布がギザギザに離散化される
- ペナルティ法の影響: ペナルティ係数が小さいと侵入量が大きくなり面圧が過小になる。ラグランジュ法に切り替える
- 球面の幾何精度: CADの球面がファセット化されていると接触が不均一になる。二次要素(HEX20)の中間節点を球面上に正確に配置する
中間節点の位置って重要なんですか?
非常に重要だ。HEX20の中間節点が辺の正確な中点にある(直線辺)のと、球面の正確な位置にある(曲線辺)のとで、接触面の幾何精度が全く異なる。Abaqusの*SURFACE, TYPE=ELEMENT で curved elements を有効にするか、メッシュ生成時に中間節点を正確な曲面上に配置することで改善する。
荷重-変位関係の検証
$P \propto \delta^{3/2}$の検証方法を教えてください。
反力 $P$ と近づき量 $\delta$ をHistory出力で記録し、両対数プロットで傾きが1.5になるか確認する。最小二乗法でフィットし、指数 $n$ を算出する。理論値は$n = 1.5$で、$n = 1.48 \sim 1.52$ が許容範囲。
この指数のチェックは個別の値($a$、$p_0$)の比較より堅牢だ。メッシュ粗さの影響を受けにくく、接触アルゴリズム全体の妥当性を示す指標になる。
「解析が合わない」と思ったら
- まず深呼吸——焦って設定をランダムに変えると、問題がさらに複雑になる
- 最小再現ケースを作る——ヘルツ接触理論(球−平面)の問題を最も単純な形で再現する。「引き算のデバッグ」が最も効率的
- 1つだけ変えて再実行——複数の変更を同時に行うと、何が効いたか分からなくなる。科学実験と同じ「対照実験」の原則
- 物理に立ち返る——計算結果が「重力に逆らって物が浮く」ような非物理的な結果なら、入力データの根本的な間違いを疑う
関連トピック
なった
詳しく
報告