(57)【要約】【課題】  均一面圧要素に比べて接触面圧が精度良く算出でき、三角形要素に比べてメッシュ分割が容易で、かつ有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単な接触面圧の計算方法を提供する。【解決手段】  接触領域よりやや大きめの四角形領域Ｒを設定し、この領域を直交する２辺の方向に並ぶ複数の四角形要素ｅにメッシュ切りする。任意の四角形要素ｅを細分化する。各四角形要素ｅの内部に形状関数に従う線形な接触面圧分布が生じると仮定する。四角形要素ｅが各節点の接触面法線方向の変位に及ぼす影響を、全四角形要素について重ね合わせた値が、前記各節点の前記接触面法線方向の変位と等しいとする等式を求める。各四角形要素ｅが支持する荷重の和が、接触物体間に作用する外力に等しいとする等式を求める。両等式を満足する各四角形要素の接触面圧を求める。【選択図】        図４
2008262263.tif 000002
【特許請求の範囲】
【請求項１】
接触領域がほぼ平坦である２つの接触物体間の接触面の接触面圧を、コンピュータを用いて計算する方法であって、
前記接触物体の接触面付近の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力過程と、入力された値を用いて前記接触面圧を計算する解析過程と、計算した結果を出力する出力過程とでなり、
前記解析過程として、
前記接触領域よりやや大きめの四角形領域を設定し、この四角形領域を直交する２辺の方向に並ぶ複数の四角形要素にメッシュ切りする過程と、
各四角形要素の内部に形状関数に従う線形な接触面圧分布が生じると仮定し、四角形要素が各節点の接触面法線方向の変位に及ぼす影響を、全四角形要素について重ね合わせた値が、前記各節点の前記接触面法線方向の変位と等しいとする等式を求める過程と、
各四角形要素が支持する荷重の和が、接触物体間に作用する外力に等しいとする等式を求める過程と、
上記両等式を満足する各四角形要素の接触面圧を求める過程と、
を含むことを特徴とする接触面圧の計算方法。
【請求項２】
請求項１において、前記メッシュ切りする過程において、メッシュ切りした各四角形要素のうちの任意の四角形要素を、上記２辺のうちの１辺の方向につき細分化して複数の細分化四角形要素とし、この細分化は、任意の細分化指定位置における前記１辺の方向に沿う両側に、任意の細分化幅Δａを設定して、細分化指定位置からΔａ離れた点から細分化指定位置に近づく毎に二分の１ずつ幅が狭くなる等比級数的な幅を持つ複数の細分化四角形要素を設定する処理とし、前記各等式を求める過程において、細分化された四角形要素については、細分化四角形要素を用いて前記各等式を求める接触面圧の計算方法。
【請求項３】
請求項１または請求項２において、前記四角形要素の内部に生じると仮定する線形な接触面圧分布が、次式（１）で表される分布であり、
【数１】
2008262263.tif 000018
前記四角形要素が各節点の接触面法線方向の変位に及ぼす影響を、全四角形要素について重ね合わせた値が、前記各節点の前記接触面法線方向の変位と等しいとする等式が、次式（９）で表される等式であり、
【数２】
2008262263.tif 000019
前記各四角形要素が支持する荷重の和が、接触物体間に作用する外力に等しいとする等式が、次式（１４）で表される式である、
【数３】
2008262263.tif 000020
接触面圧の計算方法。
【請求項４】
接触領域がほぼ平坦である２つの接触物体間の接触面における表面下応力を、コンピュータを用いて計算する方法において、
前記接触領域よりやや大きめの四角形領域を設定し、この四角形領域を直交する２辺の方向に並ぶ複数の四角形要素にメッシュ切りした場合における各四角形要素の接触面圧が既知の値として与えられている場合の表面下応力の計算方法であって、
各四角形要素の内部に形状関数に従う線形な接触面圧分布が生じると仮定して、前記表面下応力を計算することを特徴とする接触物体の表面下応力の計算方法。
【請求項５】
接触領域がほぼ平坦である２つの接触物体間の接触面における表面下応力を、コンピュータを用いて計算する方法であって、
請求項１ないし請求項３のいずれか１項に記載の接触面圧の計算方法を含み、接触面圧を計算した後、
接触物体間の接触面における表面下応力の計算過程として、
前記外力が前記接触面の法線方向の力である場合は、次式（１８）〜（２３）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々の四角形要素の全域に渡り積分し、かつ全四角形要素の影響の和をとって表面下応力を求め、
【数４】
2008262263.tif 000021
前記外力が前記接触面の接線方向の力である場合は、次式（ｄ）〜（ｉ）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々の四角形要素の全域に渡り積分し、かつ全四角形要素の影響の和をとって表面下応力を求め、
【数５】
2008262263.tif 000022
前記外力が前記接触面の法線方向の力および接線方向の力を含む場合は、前記法線方向の力である場合の表面下応力の計算、および前記接線方向の力である場合の表面下応力の計算の両方を行い、両方の計算結果を重ね合わせる、
接触物体の表面下応力の計算方法。
【請求項６】
接触領域がほぼ平坦である２つの接触物体間の接触面の接触面圧を、コンピュータを用いて計算する装置であって、
前記接触物体の接触面付近の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力処理手段と、入力された値を用いて前記接触面圧を計算する解析手段と、計算した結果を出力する出力処理手段とでなり、
前記解析手段として、
前記接触領域よりやや大きめの四角形領域を設定し、この四角形領域を直交する２辺の方向に並ぶ複数の四角形要素にメッシュ切りするメッシュ分割手段と、
各四角形要素の内部に形状関数に従う線形な接触面圧分布が生じると仮定し、四角形要素が各節点の接触面法線方向の変位に及ぼす影響を、全四角形要素について重ね合わせた値が、前記各節点の前記接触面法線方向の変位と等しいとする等式を求める変位影響算出手段と、
各四角形要素が支持する荷重の和が、接触物体間に作用する外力に等しいとする等式を求める荷重和・外力関係作成手段と、
上記両等式を満足する各四角形要素の接触面圧を求める等式解計算手段と、
を含む接触物体の接触面圧の計算装置。
【発明の詳細な説明】
【技術分野】
【０００１】
この発明は、接触領域がほぼ平坦である２物体間の接触面圧と表面下応力の簡易で高精度な計算方法に関するものであり、例えば、ころ軸受等の転がり軸受における軌道輪と転動体との接触部における接触面圧計算と表面下応力計算に利用される。
【背景技術】
【０００２】
ころと内輪またはころと外輪の接触のように、接触領域がほぼ平坦とみなせる場合の接触面圧計算方法は、半無限体に集中荷重が作用したときの解析解（ブジネスク（Boussinesq）の解やセルッティ（Cerruti ）の解）の重ね合わせにより、数値的に面圧分布を求める手法が知られており、各種文献（例えば非特許文献１）が報告されている。また、この手法をベースとした面圧計算ソフトフェアも市販されている。
【０００３】
しかしほとんど場合、メッシュ内部の面圧は均一と仮定されており、線形な面圧分布を仮定したものは少ない。メッシュ内部に線形面圧分布を仮定した要素を用いる手法としては、栫井の三角形要素（非特許文献２）がある。
【特許文献１】公開昭５９−６９５１１９号
【非特許文献１】ハーネット（Hartnett, M. J.)著, 「転がり要素軸受の接触応力解析」（“The Analysis of Contact Stresses in Rolling Element Bearings”）, ASME J. of Lubr., 101, Jan.,p105, ,(1979).
【非特許文献２】栫井著, 「半無限体近似による三次元転がり接触の専用境界要素法による解析」, 機論, 57, 544, (1991), 147-151
【非特許文献３】川井忠彦著．「機械のための有限要素法入門」（オーム社）第68ページ
【発明の開示】
【発明が解決しようとする課題】
【０００４】
非特許文献１などのように、均一な面圧分布を仮定した場合、メッシュをかなり細かくしないと、接触面圧や表面下応力値が正確に計算できないという問題がある。メッシュを細かくするに従い、計算に時間がかかる。
非特許文献２のように、線形面圧分布を仮定した三角形要素を用いる方法は、三角形であるため、その補間関数の特徴としてメッシュ切りパターンの影響を受け易い（この課題は、例えば非特許文献３に説明されている）。また、ころと内輪のような線接触の場合、接触領域はほぼ四角形に近く、三角形要素であると扱い難いという問題がある。
また、有限要素法や境界要素法は、接触面圧や表面下応力値を精度良く計算できるが、計算時間が多大になるという課題がある。
【０００５】
この発明の目的は、均一面圧要素に比べて接触面圧が精度良く算出でき、また三角形要素に比べてメッシュ切りパターンの影響を小さくできてメッシュ分割が容易であり、かつ有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単な接触面圧の計算方法および計算装置を提供することである。
この発明の他の目的は、均一面圧要素に比べて表面下応力値が精度良く算出でき、また三角形要素に比べてメッシュ切りパターンの影響を小さくできてメッシュ分割が容易であり、かつ有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単な接触物体の表面下応力の計算方法を提供することである。
【課題を解決するための手段】
【０００６】
この発明の接触面圧の計算方法は、接触領域がほぼ平坦である２つの接触物体間の接触面の接触面圧を、コンピュータを用いて計算する方法であって、
前記接触物体の接触面付近の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力過程と、入力された値を用いて前記接触面圧を計算する解析過程と、計算した結果を出力する出力過程とでなり、
前記解析過程として、
前記接触領域よりやや大きめの四角形領域を設定し、この四角形領域を直交する２辺の方向に並ぶ複数の四角形要素にメッシュ切りする過程と、
各四角形要素の内部に形状関数に従う線形な接触面圧分布が生じると仮定し、四角形要素が各節点の接触面法線方向の変位に及ぼす影響を、全四角形要素について重ね合わせた値が、前記各節点の前記接触面法線方向の変位と等しいとする等式を求める過程と、
各四角形要素が支持する荷重の和が、接触物体間に作用する外力に等しいとする等式を求める過程と、
上記両等式を満足する各四角形要素の接触面圧を求める過程と、
を含むことを特徴とする。
上記「節点」とは各要素の結合点を示す単語で、上記要素のような一次要素の場合、節点は図形の頂点になる。なお、二次要素では図形の辺にも節点を設けて、要素を結合させる。
【０００７】
この方法によると、四角形要素の内部の面圧を、形状関数を用いて線形補間（すわなち１次関数で補間）した要素を用いて接触面圧分布を求めるため、均一面圧要素に比べて接触面圧が精度良く算出でき、また三角形要素に比べてメッシュ切りパターンの影響を小さくでき、メッシュ分割が容易である。また、有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単である。
【０００８】
この発明において、前記メッシュ切りする過程において、メッシュ切りした各四角形要素のうちの任意の四角形要素を、上記２辺のうちの１辺の方向につき細分化して複数の細分化四角形要素とし、この細分化は、任意の細分化指定位置における前記１辺の方向に沿う両側に、任意の細分化幅Δａを設定して、細分化指定位置からΔａ離れた点から細分化指定位置に近づく毎に二分の１ずつ幅が狭くなる等比級数的な幅を持つ複数の細分化四角形要素を設定する処理とし、前記各等式を求める過程において、細分化された四角形要素については、細分化四角形要素を用いて前記各等式を求めるようにしても良い。細分化する箇所は何箇所でも良く、例えば、エッジロードが発生する箇所や、面圧分布が急激に変化する箇所に設定する。
このように、任意の四角形要素を細分化した場合、必要部分のみにつき、より正確な面圧分布を求め、また全体のメッシュ切りを細かくする場合と異なり、計算時間の増加が抑えられる。
【０００９】
この発明の接触面圧の計算方法は、具体的には次の方法が採用できる。すなわち、前記四角形要素の内部に生じると仮定する線形な接触面圧分布を、次式（１）で表される分布とする。
【数６】
2008262263.tif 000003
【００１０】
前記四角形要素が各節点の接触面法線方向の変位に及ぼす影響を、全四角形要素について重ね合わせた値が、前記各節点の前記接触面法線方向の変位と等しいとする等式は、次式（９）で表される等式とする。
【数７】
2008262263.tif 000004
【００１１】
前記各四角形要素が支持する荷重の和が、接触物体間に作用する外力に等しいとする等式は、次式（１４）で表される式とする。
【数８】
2008262263.tif 000005
【００１２】
この発明における第１の接触物体の表面下応力の計算方法は、接触領域がほぼ平坦である２つの接触物体間の接触面における表面下応力を、コンピュータを用いて計算する方法において、
前記接触領域よりやや大きめの四角形領域を設定し、この四角形領域を直交する２辺の方向に並ぶ複数の四角形要素にメッシュ切りした場合における各四角形要素の接触面圧が既知の値として与えられている場合の表面下応力の計算方法であって、
各四角形要素の内部に形状関数に従う線形な接触面圧分布が生じると仮定して、前記表面下応力を計算することを特徴とする。
この方法によると、各節点での接触面圧値を用いて面圧分布を線形補間し、表面下応力を計算するため、均一面圧要素に比べて表面下応力が精度良く算出でき、また三角形要素に比べてメッシュ切りパターンの影響を小さくでき、メッシュ分割が容易である。また、有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単である。
【００１３】
この発明における他の接触物体の表面下応力の計算方法は、接触領域がほぼ平坦である２つの接触物体間の接触面における表面下応力を、コンピュータを用いて計算する方法であって、この発明の接触面圧の計算方法を含み、この発明の計算方法で接触面圧を計算した後に、次の各過程を行う。
すなわち、接触物体間の接触面における表面下応力の計算過程として、
前記外力が前記接触面の法線方向の力である場合は、次式（１８）〜（２３）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々の四角形要素の全域に渡り積分し、かつ全四角形要素の影響の和をとって表面下応力を求める。
【数９】
2008262263.tif 000006
【００１４】
前記外力が前記接触面の接線方向の力である場合は、次式（ｄ）〜（ｉ）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々の四角形要素の全域に渡り積分し、かつ全四角形要素の影響の和をとって表面下応力を求める。
【数１０】
2008262263.tif 000007
【００１５】
前記外力が前記接触面の法線方向の力および接線方向の力を含む場合は、前記法線方向の力である場合の表面下応力の計算、および前記接線方向の力である場合の表面下応力の計算の両方を行い、両方の計算結果を重ね合わせる。
【００１６】
この発明の接触面圧の計算装置は、接触領域がほぼ平坦である２つの接触物体間の接触面の接触面圧を、コンピュータ１を用いて計算する装置であって、
前記接触物体の接触面付近の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力処理手段４と、入力された値を用いて前記接触面圧を計算する解析手段７と、計算した結果を出力する出力処理手段５とでなる。
前記解析手段７には、
前記接触領域よりやや大きめの四角形領域を設定し、この四角形領域を直交する２辺の方向に並ぶ複数の四角形要素にメッシュ切りするメッシュ分割手段８と、
各四角形要素の内部に形状関数に従う線形な接触面圧分布が生じると仮定し、四角形要素が各節点の接触面法線方向の変位に及ぼす影響を、全四角形要素について重ね合わせた値が、前記各節点の前記接触面法線方向の変位と等しいとする等式を求める変位影響算出手段９と、
各四角形要素が支持する荷重の和が、接触物体間に作用する外力に等しいとする等式を求める荷重和・外力関係作成手段１０と、
上記両等式を満足する各四角形要素の接触面圧を求める等式解計算手段１１と、
を含む。
【００１７】
この構成の接触面圧の計算装置によると、この発明の接触面圧の計算方法が実施でき、四角形要素の内部の面圧を、形状関数を用いて線形補間した要素を用いて接触面圧分布を求めるため、均一面圧要素に比べて接触面圧が精度良く算出でき、また三角形要素に比べてメッシュ切りパターンの影響を小さくでき、メッシュ分割が容易である。また、有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単である。
【００１８】
なお、この発明の接触面圧の計算方法、計算装置、および表面下応力の計算方法は、例えば、前記２つの接触物体が、ころ軸受の軌道輪ところであっても良い。その場合、前記直交する２辺の方向を、ころの軸方向、およびころの転動方向とする。前記細分化を行う場合は、ころの軸方向に細分化する。
【発明の効果】
【００１９】
この発明の接触面圧の計算方法および計算装置は、四角形要素の内部の面圧を、形状関数を用いて線形補間した要素を用いて接触面圧分布を求めるため、均一面圧要素に比べて接触面圧が精度良く算出でき、また三角形要素に比べてメッシュ切りパターンの影響を小さくでき、メッシュ分割が容易である。また、有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単である。
この発明の接触物体の表面下応力の計算方法は、面圧分布を線形補間し、表面下応力を計算するため、均一面圧要素に比べて表面下が精度良く算出でき、また三角形要素に比べてメッシュ切りパターンの影響を小さくでき、メッシュ分割が容易である。また、有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単である。
【発明を実施するための最良の形態】
【００２０】
この発明の一実施形態を図１ないし図７と共に説明する。図１は、この発明方法の各過程の流れ図を、図２はこの発明方法を実施する装置の概念構成のブロック図を、図３は２つの接触物体Ｍ１，Ｍ２の例をそれぞれ示す。
【００２１】
この接触面圧および表面下応力計算方法は、接触領域がほぼ平坦である２つの接触物体Ｍ１，Ｍ２（図３）の接触面における接触面圧および表面下応力を、コンピュータを用いて計算する方法であって、計算に必要な各種のデータを入力する入力過程（Ｓ１）と、入力された値を用いて接触面圧および表面下応力を計算する解析過程（Ｓ２）と、計算した結果を出力する出力過程（Ｓ３）とでなる。解析過程（Ｓ２）として、順次行われる次のメッシュ分割過程（Ｔ１）、変位影響算出過程（Ｔ２）、荷重和・外力関係作成過程（Ｔ３）、等式解計算過程（Ｔ４）、および表面下応力計算過程（Ｔ５）を含む。メッシュ分割過程（Ｔ１）は、細分化過程（Ｔ１ａ）を含む。
【００２２】
計算装置となるコンピュータ１（図２）は、中央処理装置（ＣＰＵ）およびメモリ等の記憶手段（いずれも図示せず）を有し、このコンピュータ１のハードウェア（オペレーションシステムを含む）、およびこのコンピュータ１に実行させる接触面圧・表面下応力計算プログラム（図示せず）により、図２に概念構成で示す各手段が構成される。
【００２３】
入力処理手段４は、入力装置２から、計算に必要な各種のデータを入力して入力情報記憶手段６に記憶させる手段である。この入力装置２から入力して入力処理手段４により入力情報記憶手段６に記憶させる過程が、図１（Ａ）の入力過程（Ｓ１）である。入力するデータは、接触物体Ｍ１，Ｍ２の接触面付近の形状、弾性係数（ここでは、ヤング率Ｅと、ポアソン比ν）、および作用する外力Ｆ等である。ここで言う「形状」は、大きさを含む概念である。入力処理手段４は、上記の形状等の情報の他に、後述の解析手段７でメッシュ分割する枡目の各軸方向の大きさや四角形要素の個数（ｎa ，ｎb ）、細分化指定位置Ａや細分化幅ΔＡ等を、入力装置２から入力して入力情報記憶手段６に記憶させる機能を備えるものとしても良い。入力処理手段４は、入力装置２から入力された情報に対して、解析手段７で演算するための前処理を行う機能を有するものとしても良い。
入力装置２は、キーボード等であっても、また通信手段や記憶素子等であっても良い。入力情報記憶手段６は、コンピュータ１の上記記憶手段における所定の記憶領域である。出力処理手段５は、解析手段７で計算した結果を、液晶ディスプレイ，プリンタ，あるいは通信機器等からなる出力装置３に出力する手段である。
【００２４】
解析手段７は、入力処理手段４で入力されて入力情報記憶手段６に記憶されたデータを用いて前記面圧および表面下応力を計算する手段であり、メッシュ分割手段８、変位影響算出手段９、荷重和・外力関係作成手段１０、等式解計算手段１１、および表面下応力計算手段１２を備える。メッシュ分割手段８は、細分化部８ａを有する。解析手段７に備えられた上記各手段８〜１２、および細分化部８ａは、それぞれ図１（Ｂ）の流れ図におけるメッシュ分割過程（Ｔ１）、変位影響算出過程（Ｔ２）、荷重和・外力関係作成過程（Ｔ３）、等式解計算過程（Ｔ４）、表面下応力計算過程（Ｔ５）、およびメッシュ分割過程（Ｔ１）における細分化過程（Ｔ１ａ）を実行する手段であり、その実行に必要な計算式や設定値等を記憶している。上記各手段８〜１２、および細分化部８ａの具体的な内容は、次に述べる上記各過程（Ｔ１）〜（Ｔ５），（Ｔ１ａ）の説明の通りである。
【００２５】
図１（Ａ）の解析過程（Ｓ２）を、同図（Ｂ）および図３〜図７と共に説明する。
メッシュ分割過程（Ｔ１）では、２つの接触物体Ｍ１，Ｍ２（図３）の接触領域よりやや大きめの四角形領域Ｒを設定し（図４（Ａ））、この四角形領域を直交する２辺の方向に並ぶ複数の四角形要素ｅ（ｅij：ｉ＝１〜ｎa ，ｊ＝１〜ｎb ）にメッシュ切りする。  この場合に、２つの接触物体Ｍ１，Ｍ２が、図３のように、ころ軸受におけるころと内輪（または外輪）等の軌道輪である場合、ころ軸方向をｘ方向、ころ転動方向をｙ方向にとる。ついで、図４のように、この前記四角形領域Ｒのｘ方向、ｙ方向を、それぞれ幅２×ａ，２×ｂの四角形要素ｅで、ｎa ，ｎb 個に均等分割する。
【００２６】
さらに、細分化過程（Ｔ１ａ）として、ころ軸方向ｘの任意の位置（図ではＡ点）で、四角形要素ｅを軸方向に細分化する。細分化の方法は、細分化する指定位置Ａの両側に任意の細分化幅Δａを初期値として設定し、Δａ離れた点から指定位置Ａに近づく毎に二分の一ずつ幅が狭くなる等比級数的方法をとる。最小の細分化幅をどの程度とするかは、適宜設定しておく。
細分化指定位置Ａと細分化幅Δａは、軸方向ｘの任意の位置に何個でも設定でき、一般にはエッジロードが発生する箇所や、面圧分布が急激に変化する箇所に設定する。ころ軸受の場合、例えば、接触物体Ｍ１であるころの端部を細分化指定位置Ａとする。
なお、細分化された四角形要素ｅについては、後述の各等式（９），（１４）等を求める過程において、細分化四角形要素ｅ′を用いて求める。したがって、ｘ軸方向の四角形要素ｅの個数ｎa は、細分化に従い、元の個数よりも増えることにある。
【００２７】
変位影響算出過程（Ｔ２）は、各四角形要素ｅの内部に形状関数に従う線形な接触面圧分布が生じると仮定し、四角形要素ｅが各節点Ｔの接触面法線方向ｚの変位に及ぼす影響を、全四角形要素ｅについて重ね合わせた値が、前記各節点Ｔの前記接触面法線方向の変位と等しいとする等式（９）を求める過程である。
すなわち、四角形要素ｅの内部に、以下のような線形な面圧分布を仮定する。
【数１１】
2008262263.tif 000008
【００２８】
四角形要素ｅijが節点Ｔkl（ｘK,ｙL ）のｚ方向（法線方向：接触面に垂直方向）変位に及ぼす影響は次式で表される。
【数１２】
2008262263.tif 000009
これを全要素について重ね合わせることにより、節点Ｔkl（ｘK,ｙL ）のｚ方向変位（Ｄij）が、次式（９）で示されるように得られる。
【００２９】
【数１３】
2008262263.tif 000010
【００３０】
なお、上記（１０′）式は、一般的にしめせば、次の式（１０）のとおりである。
【数１４】
2008262263.tif 000011
【００３１】
荷重和・外力関係作成過程（Ｔ３）では、各四角形要素ｅが支持する荷重の和が、接触物体間に作用する外力に等しいとする等式（１４）を求める。すなわち、外力Ｆとの釣合から次式が成り立つ。
【数１５】
2008262263.tif 000012
【００３２】
等式解計算過程（Ｔ４）は、上記の２つの等式（９），（１４）を満足するように、未知数である各四角形要素ｅの接触面圧Ｐijと、二物体間の接近量αを求める。例えば、接触面圧Ｐijの解法には逐次反復法、接近量αの解法にはニュートン法を用い、図６の流れ図で示される方法で解けば、解が得られる。
【００３３】
図６の流れ図は、等式（９），（１４）を、コンピュータにおける等式解計算手段１１によって解くアルゴリズムの一例を示す。
まず、接近量αの初期値を設定する（ステップＨ１）。そのα値に応じて各四角形要素ｅijの接触面圧Ｐijを設定する（ステップＨ２）。この場合に、ｉ＝１〜ｎa の範囲、およびｊ＝１〜ｎb の範囲で全ての四角形要素ｅijの接触面圧Ｐijを設定する。その設定した各四角形要素ｅijの接触面圧Ｐijで（９）式が成り立つか否かを判定し、成り立たない場合は、四角形要素ｅijの接触面圧の設定過程（ステップＨ２）に戻る（Ｈ３）。ステップＨ２では、各四角形要素ｅijの接触面圧Ｐijを再度設定し直し、その設定した四角形要素ｅijの接触面圧Ｐijで（９）式が成り立つか否かを判定し、成り立たない場合は四角形要素ｅijの接触面圧Ｐijの設定過程（Ｈ２）に戻る。このような処理を（９）式が成り立つまで繰替し、（９）式が成り立つと、ステップＨ４に進む。
ステップＨ４では、（１４）式が成り立つか否かを判定し、成り立たない場合は、ステップＨ１に戻ってα値を、再度設定し直す。この後、前記と同様にしてステップＨ２〜Ｈ４，Ｈ１の処理を繰り返し、（１４）式が成り立つと、処理を終了する。
このようにして求めた、式（９）、（１４）が成り立つときの、ステップＨ２で設定した各四角形要素ｅijの接触面圧Ｐijが、等式解計算手段１１の出力となる接触面圧の値である。
なお、四角形要素ｅijの接触面圧の最初の値は何でもよいが、できるだけ実際の面圧に近い値を設定した方が早く収束する。
【００３４】
表面下応力計算過程（Ｔ５）では、表面下応力を以下のように求める。ここでは剪断応力τzyについて示す。
一つの四角形要素ｅに（１）式のような線形接触面圧分布を仮定したとき、表面下の任意の点Ｑkl（ｘk , ｙl , ｚ）に働く剪断応力τzyijklは次式のようになる。
【数１６】
2008262263.tif 000013
【００３５】
点Ｑkl（ｘk , ｙl , ｚ）に働く剪断応力τzyijklは、すべての四角形要素ｅijについて重ね合わせて次式で与えられる。
【数１７】
2008262263.tif 000014
【００３６】
上記内容は、次式で与えられる単位集中荷重作用時の表面下の応力のうち、τyzを用いているが、他の成分についても同様に求めることができる。
【数１８】
2008262263.tif 000015
【００３７】
また、上記内容は法線方向荷重作用時のブジネスク（Boussinesq）の解を用いたものであるが、接線方向荷重の作用時は以下の示すセルッティ（Cerruti ）の解（下記）を用いて、上記手法をそのまま適用すれば求めることができる。
法線方向と接線方向の両方の荷重が作用したときは両者を重ね合わせる。
【００３８】
【数１９】
2008262263.tif 000016
【００３９】
つぎに、この実施形態の接触面圧および表面下応力の計算方法，計算装置ににより実際に計算した結果につき、従来例と比較して説明する。面圧計算結果は一定面圧を仮定した従来の方法とあまり差はないので省略し、剪断応力についての計算例を掲載する。
図８のような二次元の三角形面圧分布（面に直行する方向には無限に面圧分布がある）が与えられたときの表面下剪断応力τyzの理論解〔ジョンソン，ケー，エル（Johnson, K. L.）著，「接触メカニズム」（,"Contact Mechanics",） Cambridge Univ. Press, 1985. 〕を次式に示す。
【数２０】
2008262263.tif 000017
【００４０】
Ｐ０＝300 、ｂ＝0.226mm としたときの表面下0.01mm（z=0.01）の位置での剪断応力計算結果を図９に示す。
なお、三角形面圧を表現するために一定面圧要素で行ったＹ方向分割を図１０に、この手法の分割を図１１に示す。線形面圧要素では面圧を２分割で表現している。
図９からわかるように、一定面圧要素では精度良く計算できていないのに対して、線形面圧要素ではＹ方向の分割数が少ないにも拘らず理論解と良く一致している。
【００４１】
以上のようにこの接触面圧および表面下応力の計算方法，計算装置によると、四角形要素の内部の面圧を、形状関数を用いて線形補間した要素を用いて接触面圧分布、表面下応力を求めるため、均一面圧要素に比べて接触面圧が精度良く算出できる。また、三角形要素に比べてメッシュ切りパターンの影響を小さくでき、メッシュ分割が容易である。さらに、有限要素法や境界要素法と比べて計算時間が短く、使用法が簡単である。
【図面の簡単な説明】
【００４２】
【図１】（Ａ）はこの発明の一実施形態に係る接触面圧および表面下応力の計算方法を示す流れ図、（Ｂ）はその解析過程の流れ図である。
【図２】この発明の一実施形態に係る接触面圧および表面下応力の計算装置の概念構成を示すブロック図である。
【図３】同計算方法の計算対象とする接触物体の一例を示す破断斜視図である。
【図４】（Ａ）は同実施形態における計算方法のメッシュ切り過程の説明図、（Ｂ）はその細分化方法の説明図である。
【図５】同計算方法における要素内局所座標の説明図である。
【図６】同計算方法における等式解計算過程の流れ図である。
【図７】同計算方法における表面化応力の説明図である。
【図８】従来の計算方法における三角形面圧分布の説明図である。
【図９】この計算方法における線形面圧要素を用いた場合と、理論解と、一定面圧分布とをにつき求めた表面下せん断応力の説明図である。
【図１０】従来の一定面圧要素での面圧分布の説明図である。
【図１１】同実施形態における線形面圧要素のｙ方向位置と面圧の関係を示す説明図である。
【符号の説明】
【００４３】
１…コンピュータ
４…入力処理手段
５…出力処理手段
６…入力情報記憶手段
７…解析手段
８…メッシュ分割手段
８ａ…細分化手段
９…変位影響算出手段
１０…荷重和・外力関係作成手段
１１…等式解計算手段
１２…表面下応力計算手段
Ｍ１，Ｍ２…接触物体
ｅ…四角形要素
ｅ′…細分化四角形要素
【図１】
2008262263.tif 000023
【図２】
2008262263.tif 000024
【図３】
2008262263.tif 000025
【図４】
2008262263.tif 000026
【図５】
2008262263.tif 000027
【図６】
2008262263.tif 000028
【図７】
2008262263.tif 000029
【図８】
2008262263.tif 000030
【図９】
2008262263.tif 000031
【図１０】
2008262263.tif 000032
【図１１】
2008262263.tif 000033