JP5026225B2 - 共形接触下の接触面圧および表面下応力計算方法・計算装置 - Google Patents
共形接触下の接触面圧および表面下応力計算方法・計算装置 Download PDFInfo
- Publication number
- JP5026225B2 JP5026225B2 JP2007277066A JP2007277066A JP5026225B2 JP 5026225 B2 JP5026225 B2 JP 5026225B2 JP 2007277066 A JP2007277066 A JP 2007277066A JP 2007277066 A JP2007277066 A JP 2007277066A JP 5026225 B2 JP5026225 B2 JP 5026225B2
- Authority
- JP
- Japan
- Prior art keywords
- contact
- surface pressure
- calculation
- displacement
- calculating
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000004364 calculation method Methods 0.000 title claims description 138
- 238000000034 method Methods 0.000 claims description 131
- 238000006073 displacement reaction Methods 0.000 claims description 78
- 238000004458 analytical method Methods 0.000 claims description 44
- 238000005096 rolling process Methods 0.000 claims description 19
- 238000013459 approach Methods 0.000 claims description 17
- 238000012937 correction Methods 0.000 claims description 16
- 230000014509 gene expression Effects 0.000 claims description 14
- 238000012545 processing Methods 0.000 claims description 14
- 238000005259 measurement Methods 0.000 claims description 13
- 230000000694 effects Effects 0.000 claims description 6
- CIWBSHSKHKDKBQ-JLAZNSOCSA-N Ascorbic acid Chemical compound OC[C@H](O)[C@H]1OC(=O)C(O)=C1O CIWBSHSKHKDKBQ-JLAZNSOCSA-N 0.000 description 24
- 230000010354 integration Effects 0.000 description 7
- 230000006870 function Effects 0.000 description 4
- 238000012887 quadratic function Methods 0.000 description 4
- 238000004891 communication Methods 0.000 description 2
- 238000007796 conventional method Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000012886 linear function Methods 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 1
- 239000004973 liquid crystal related substance Substances 0.000 description 1
- 239000011159 matrix material Substances 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
Images
Description
この発明は、接触領域が平坦と見なせない接触物体間の接触面圧および表面下応力の計算を、有限要素法や境界要素法を使用することなく行う計算方法および装置に関する。
接触面圧の計算問題は下記の二つの分野に分けられる。
一つは非共形接触問題（Counter-formal（またはNon-Conforming）Problem ）であり、接触領域が接触部の物体半径に比べて小さい場合の問題である。
他の一つは、共形接触問題（Conformal （またはConforming）Problem ）であり、接触領域が接触部の物体半径に比べて大きい場合の問題である。
一つは非共形接触問題（Counter-formal（またはNon-Conforming）Problem ）であり、接触領域が接触部の物体半径に比べて小さい場合の問題である。
他の一つは、共形接触問題（Conformal （またはConforming）Problem ）であり、接触領域が接触部の物体半径に比べて大きい場合の問題である。
非共形接触問題（Counter-formal Problem）では、接触領域はほぼ平面上にあると見なせるため、比較的解き易い。
特に接触面が二次関数で表される曲面の場合は、ヘルツ（Hertz ）接触理論などが広く用いられており、数表を用いて電卓レベルで計算が可能である。二次関数で表すことができない曲面の場合は、半無限体に集中荷重が作用する場合の解析解（ブジネスク（Boussinesq）の解の重ね合わせにより、数値的に面圧分布を求める手法が知られている。また、これらと類似の解法で各種文献（例えば非特許文献１や非特許文献２）が出ており、この手法をベースとした面圧計算ソフトウェアも市販されている。
特開昭５９−６９５１９号公報
ハーネット（Hartnett, M. J. ）著 ,「転がり要素軸受の接触応力解析」（"The Analysis of Contact Stresses in Rolling Element Bearings"）, ASMEJ. of Lubr., 101, Jan.p105,(1979),
栫井著, 「三次元弾性接触問題の境界要素解析」, トライボロジスト, 36,6,(1991), 429-433
特に接触面が二次関数で表される曲面の場合は、ヘルツ（Hertz ）接触理論などが広く用いられており、数表を用いて電卓レベルで計算が可能である。二次関数で表すことができない曲面の場合は、半無限体に集中荷重が作用する場合の解析解（ブジネスク（Boussinesq）の解の重ね合わせにより、数値的に面圧分布を求める手法が知られている。また、これらと類似の解法で各種文献（例えば非特許文献１や非特許文献２）が出ており、この手法をベースとした面圧計算ソフトウェアも市販されている。
しかし、共形接触問題（Conformal Problem ）については、上記のような解析手法は適用できないため、有限要素法や境界要素法を用いて解く以外に計算方法がない。
有限要素法や境界要素法を使用して計算する場合、計算時間が多大になるという課題がある。
有限要素法や境界要素法を使用して計算する場合、計算時間が多大になるという課題がある。
このような課題を解消する方法，装置として、本出願人は、次に示す方法，装置を先に提案した（特願２００７−０２６５７７号）。
この提案例は、半無限体に集中荷重が作用したときの変位や表面下応力の解析解を共形接触問題に適用して面圧および表面下応力を計算する手法であり、法線方向荷重作用時はブジネスク（Boussinesq）の解を用いて重ね合わせる。この際に、接触２物体の曲面近傍に基準となる曲面を設置して、その基準曲面をメッシュ分割し、２物体間の距離を計測する点をメッシュ上に設置し、各セルに分布荷重が作用したときの２物体間計測点への変位量の影響を算出する。その位置での接触判定によって接触領域を決め、かつ各セルの面圧の合計が外力と釣合うように２物体間の接近量を決めて、接触面圧を算出する。
この提案例は、半無限体に集中荷重が作用したときの変位や表面下応力の解析解を共形接触問題に適用して面圧および表面下応力を計算する手法であり、法線方向荷重作用時はブジネスク（Boussinesq）の解を用いて重ね合わせる。この際に、接触２物体の曲面近傍に基準となる曲面を設置して、その基準曲面をメッシュ分割し、２物体間の距離を計測する点をメッシュ上に設置し、各セルに分布荷重が作用したときの２物体間計測点への変位量の影響を算出する。その位置での接触判定によって接触領域を決め、かつ各セルの面圧の合計が外力と釣合うように２物体間の接近量を決めて、接触面圧を算出する。
上記提案例によると、有限要素法や境界要素法を使用する方法に比べて、比較的簡単に短時間で、かつ精度良く、共形接触における接触面圧の計算、さらに表面下の応力の計算をすることができる。
しかし、上記提案例では、３方向の変位を考慮していたため、まだ計算方法が複雑で、かつ数値積分を行うために計算時間がかかる。そのため、より高速な計算を可能にすることが望まれる。
しかし、上記提案例では、３方向の変位を考慮していたため、まだ計算方法が複雑で、かつ数値積分を行うために計算時間がかかる。そのため、より高速な計算を可能にすることが望まれる。
この発明の目的は、有限要素法や境界要素法を使用することなく、非共形接触問題で行われている数値計算方法を拡張して共形接触問題にも適用できるようにし、かつ上記提案例を改善して、１方向の変位の考慮で計算できて、簡単に短時間で、かつ精度良く、共形接触における接触面圧の計算、さらに表面下の応力の計算をできるようにした方法および装置を提供することを目的とする。
この発明の共形接触下の接触面圧および表面下応力計算方法は、半無限体に集中荷重が作用したときの変位や表面下応力の解析解を共形接触問題に適用して、面圧および表面下応力を計算する手法であって、解析解として、法線方向荷重作用時はブジネスク（Boussinesq）の解を用いる。
具体的には、この発明の共形接触下の接触面圧計算方法は、外力が作用したときに互いに共形接触（接触領域がほぼ平面上にあると見なせない接触）を行う２つの接触物体Ｍ１，Ｍ２の接触面における面圧を、コンピュータを用いて計算する方法であって、
前記接触物体の接触面の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力過程（Ｓ１）と、入力された値を用いて前記面圧を計算する解析過程（Ｓ２）と、計算した結果を出力する出力過程（Ｓ３）とでなり、
前記解析過程（Ｓ２）として、次の各過程（Ｔ１，Ｔ２，Ｔ４，Ｔ５）を含む。
前記接触物体の接触面の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力過程（Ｓ１）と、入力された値を用いて前記面圧を計算する解析過程（Ｓ２）と、計算した結果を出力する出力過程（Ｓ３）とでなり、
前記解析過程（Ｓ２）として、次の各過程（Ｔ１，Ｔ２，Ｔ４，Ｔ５）を含む。
すなわち、前記２つの接触物体Ｍ１，Ｍ２の互いに接する面の間に接触面形状に近い形状の基準曲面Ｓを設定してこの基準曲面ＳをセルＣが縦横に並ぶメッシュに分割し、接触物体Ｍ１，Ｍ２間の距離を計算する点である計算点Ａをセル上に定めるメッシュ分割過程（Ｔ１）と、
各セルＣに単位分布荷重が作用したときの前記接触物体Ｍ１，Ｍ２の表面の前記各計算点Ａに対応する点に生じる変位である影響量を算出する変位影響算出過程（Ｔ２）と、
各セルＣにつき面圧Ｐと前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体間の接近量と初期隙間の差に釣り合い、かつ各セルＣの面圧Ｐの合計が外力Ｆと釣り合うという釣り合い条件式から、各セルの面圧Ｐを算出する接触面圧計算過程（Ｔ３）とを含む。なお、前記「接近量」は、前記単位分布荷重が与えられることにより互いに近づく方向に移動する距離であり、前記「初期隙間」は、前記単位分布荷重が与えられる前の距離である。
各セルＣに単位分布荷重が作用したときの前記接触物体Ｍ１，Ｍ２の表面の前記各計算点Ａに対応する点に生じる変位である影響量を算出する変位影響算出過程（Ｔ２）と、
各セルＣにつき面圧Ｐと前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体間の接近量と初期隙間の差に釣り合い、かつ各セルＣの面圧Ｐの合計が外力Ｆと釣り合うという釣り合い条件式から、各セルの面圧Ｐを算出する接触面圧計算過程（Ｔ３）とを含む。なお、前記「接近量」は、前記単位分布荷重が与えられることにより互いに近づく方向に移動する距離であり、前記「初期隙間」は、前記単位分布荷重が与えられる前の距離である。
前記変位影響算出過程（Ｔ２）および接触面圧計算過程（Ｔ３）では、計測点間の距離の算出に上記基準曲面上の距離を用いる。また、基準曲面上の位置の違いによる法線方向の修正を行って面圧を算出する。すなわち、影響係数の算出に基準曲面の法線方向の修正を行って面圧を算出する。
なお、上記基準曲面Ｓは、２つの接触物体Ｍ１，Ｍ２の互いに接する面の接触面形状に近い形状であることが好ましい。前記変位影響算出過程（Ｔ２）において、各セルＣに作用する単位分布荷重は、例えば、均一な面圧とするが、１次関数または２次関数に従って変化する面圧を仮定しても良い。接触物体Ｍ１，Ｍ２間の距離は、離散的な計測点Ａで評価するものとするが、計測点Ａは、例えばメッシュの中心にとる。
前記変位影響算出過程（Ｔ２）では、単位集中荷重作用時の変位解析解（法線方向の荷重に対しては例えばブジネスク（Boussinesq）の解を１つの要素全域にわたり積分して、要素１つの面圧分布によって生じる２つの接触物体Ｍ１，Ｍ２の面上の点（上記計測点Ａに対応）の変位量を求める。上記の「要素」は、セルのことである。
前記変位影響算出過程（Ｔ２）では、単位集中荷重作用時の変位解析解（法線方向の荷重に対しては例えばブジネスク（Boussinesq）の解を１つの要素全域にわたり積分して、要素１つの面圧分布によって生じる２つの接触物体Ｍ１，Ｍ２の面上の点（上記計測点Ａに対応）の変位量を求める。上記の「要素」は、セルのことである。
この発明方法は、前記２つの物体Ｍ１，Ｍ２が、転がり軸受の内輪とボールとである場合に適用することができる。その場合、前記メッシュ分割過程（Ｔ１）では、前記基準曲面Ｓを円筒面とする。
前記変位影響量算出過程（Ｔ２）では、前記円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面Ｓ上の法線方向に離れる座標位置をｚ、これら３軸方向の変位を、それぞれｕ，ｖ，ｗとして、
前記外力が法線方向の力である場合は、次式（１）で表される単位集中荷重作用時の変位解析解ｗを、全てのセルについて重ね合わせることで前記影響量を計算する。
前記変位影響量算出過程（Ｔ２）では、前記円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面Ｓ上の法線方向に離れる座標位置をｚ、これら３軸方向の変位を、それぞれｕ，ｖ，ｗとして、
前記外力が法線方向の力である場合は、次式（１）で表される単位集中荷重作用時の変位解析解ｗを、全てのセルについて重ね合わせることで前記影響量を計算する。
ただし、（１）式により与えられる変位ｗは位置（θ，ｙ）の要素が位置（θ′，ｙ′）に及ぼす変位。
すなわち、この発明方法において、前記変位影響量算出過程（Ｔ２）では、単位集中荷重作用時の変位解析解（Boussinesqの解）を１つの要素全域にわたり積分して、要素１つの面圧分布によって生じる２つの接触物体Ｍ１，Ｍ２の表面上の点（上記計測点Ａに対応）の変位量を求める。
位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’）に及ぼす単位面圧負荷時の変位は（１）式で与える。このとき、基準曲面Ｓの法線方向が位置により変わるためCOS(θ−θ’) がかかってくる。
ここで重要なのは、(1) 要素間のθ方向距離として基準曲面上の距離（Ｒ０θ−Ｒ０θ’）を用いる点、(2) 基準曲面Ｓの法線方向の修正を、上記COS(θ−θ’) で行うこと（すなわち、影響係数の算出に基準曲面Ｓの法線方向の修正を行って面圧を算出すること）の２点である。
上記（１）式は、解析的に積分可能であり、数値積分は不要となる。
位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’）に及ぼす単位面圧負荷時の変位は（１）式で与える。このとき、基準曲面Ｓの法線方向が位置により変わるためCOS(θ−θ’) がかかってくる。
ここで重要なのは、(1) 要素間のθ方向距離として基準曲面上の距離（Ｒ０θ−Ｒ０θ’）を用いる点、(2) 基準曲面Ｓの法線方向の修正を、上記COS(θ−θ’) で行うこと（すなわち、影響係数の算出に基準曲面Ｓの法線方向の修正を行って面圧を算出すること）の２点である。
上記（１）式は、解析的に積分可能であり、数値積分は不要となる。
前記接触面圧計算過程（Ｔ３）は、釣り合い条件式から各セルの面圧Ｐを算出する過程であり、例えば次のように計算する。
すなわち、次の（２）式が釣り合うように各セルの面圧分布を決める。
すなわち、次の（２）式が釣り合うように各セルの面圧分布を決める。
また、次の（３）式が釣り合うように接近量をαを求める。
なお、（２）式は、各セルにつき面圧と前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体間の接近量と初期隙間の差に釣り合うという釣り合い条件式であり、また（１）式で計算される変位ｗを全要素について重ね合わせる式である。（３）式は、各セルの面圧の合計が外力と釣り合うという釣り合い条件式である。
この発明方法によると、半無限体に集中荷重が作用したときの変位の解析解を共形接触問題に適用し、また接触面圧計算時の計測点間の距離の算出に上記基準曲面上の距離を用い、かつ影響係数の算出に基準曲面の法線方向の修正を行って面圧を算出するようにしたため、比較的簡単に、したがって短時間で、共形接触における接触面圧の計算を行うことができる。
また、有限要素法を用いて接触面圧を計算した場合と同程度に精度の良い計算結果が得られることが、この発明方法と有限要素法との計算結果を比較することで、確認できた。
この発明は、このように、従来では共形接触には半無限体が適用できないと考えられていた共形接触問題に、ブジネスクの半無限体近似の解を適用し、良好な面圧の解析結果を得ることができる。
なお、この発明方法は、非共形接触の場合にも適用可能であるが、非共形接触の場合はより簡単な計算で面圧を求めることができるため、共形接触の場合に効果的となる。
また、有限要素法を用いて接触面圧を計算した場合と同程度に精度の良い計算結果が得られることが、この発明方法と有限要素法との計算結果を比較することで、確認できた。
この発明は、このように、従来では共形接触には半無限体が適用できないと考えられていた共形接触問題に、ブジネスクの半無限体近似の解を適用し、良好な面圧の解析結果を得ることができる。
なお、この発明方法は、非共形接触の場合にも適用可能であるが、非共形接触の場合はより簡単な計算で面圧を求めることができるため、共形接触の場合に効果的となる。
この発明方法において、表面下応力を計算する場合は、上記のいずれかの方法で各セルＣの面圧Ｐを計算した後に、表面下応力を計算する。
この表面下応力計算過程において、
セルに対する荷重が法線方向の力である場合は、次式（４）〜（９）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々のセルＣの全域にわたり積分して、かつ全セルＣの影響の和をとって表面下応力を求める。
この表面下応力計算過程において、
セルに対する荷重が法線方向の力である場合は、次式（４）〜（９）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々のセルＣの全域にわたり積分して、かつ全セルＣの影響の和をとって表面下応力を求める。
セルに対する荷重が法線方向の力である場合は、次のように計算する。前記２つの接触物体Ｍ１，Ｍ２が転がり軸受のボールと内輪とである場合、前記メッシュ分割過程（Ｔ１）では、前記のように前記基準曲面Ｓを円筒面とする。前記変位影響量算出過程（Ｔ２）では、前記円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をθ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚとして、
位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’，ｚ）に及ぼす単位面圧負荷時の剪断応力τzyを、次式（１０）で求める。
位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’，ｚ）に及ぼす単位面圧負荷時の剪断応力τzyを、次式（１０）で求める。
このように表面下応力を計算するについても、計算が比較的簡単で短時間で行え、また精度の良い表面下応力の計算が行える。
この発明方法において、前記解析過程の途中で、前記メッシュに分割する大きさを修正してその修正結果により、メッシュ分割後の各過程を再度行うメッシュ修正過程（Ｔ６）を含むようにしても良い。
接触物体の接触面の形状によっては、接触領域の縁がギザギザになるなどして、離散化の誤差が大きくなる場合があるが、解析途中でメッシュを修正することで、離散化の誤差を低減させることができる。
接触物体の接触面の形状によっては、接触領域の縁がギザギザになるなどして、離散化の誤差が大きくなる場合があるが、解析途中でメッシュを修正することで、離散化の誤差を低減させることができる。
なお、この発明方法は、前記２つの接触物体Ｍ１，Ｍ２が、転がり軸受における内輪または外輪である軌道輪と、ボール等の転動体とである場合に、効果的に適用できる。転がり軸受の軌道輪とボール等の転動体との接触においても、接触領域がほぼ平面上にあると見なせず、共形接触の問題として面圧等を求めることが必要であるため、この発明方法を効果的に適用することができる。
この発明の共形接触下の接触面圧計算装置は、外力が作用したときに互いに共形接触を行う２つの接触物体の接触面における面圧を、計算する装置であって、
前記接触物体Ｍ１，Ｍ２の接触面の形状、弾性係数、ポアソン比、および作用する外力Ｆを入力して記憶手段６に記憶させる入力処理手段４と、入力された値を用いて前記面圧を計算する解析手段７と、計算した結果を出力する出力処理手段５とでなり、
前記解析手段７として、次の各手段８，９，１１，１２を備える。
前記接触物体Ｍ１，Ｍ２の接触面の形状、弾性係数、ポアソン比、および作用する外力Ｆを入力して記憶手段６に記憶させる入力処理手段４と、入力された値を用いて前記面圧を計算する解析手段７と、計算した結果を出力する出力処理手段５とでなり、
前記解析手段７として、次の各手段８，９，１１，１２を備える。
すなわち、前記２つの接触物体Ｍ１，Ｍ２の互いに接する面の間に、基準曲面Ｓを設定してこの基準曲面ＳをセルＣが縦横に並ぶメッシュに分割し、接触物体Ｍ１，Ｍ２間の距離を計算する点である計算点Ａをセル上に定めるメッシュ分割手段８と、
各セルＣに均一な単位分布荷重が作用したときの前記接触物体Ｍ１，Ｍ２の表面の前記各計算点Ａに対応する点の変位への影響量を算出する変位影響算出手段９と、
各セルＣにつき面圧Ｐと前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体間の接近量と初期隙間の差に釣り合、かつ各セルＣの面圧Ｐの合計が外力と釣り合うという釣り合い条件式から、各セルの面圧Ｐを算出する接触面圧計算手段１０とを含む。
なお、前記「接近量」は、前記単位分布荷重が与えられることにより互いに近づく方向に移動する距離であり、前記「初期隙間」は、前記単位分布荷重が与えられる前の距離である。
前記変位影響算出手段９および接触面圧計算手段１０は、計測点間の距離の算出に上記基準曲面Ｓ上の距離を用い、かつこの基準曲面Ｓ上の位置の違いによる法線方向の修正を行って面圧を算出する。すなわち、影響係数の算出に基準曲面Ｓの法線方向の修正を行って面圧を算出する。
前記メッシュ修正手段１２は、前記接触面圧計算手段１０による前記各セルＣの面圧Ｐの算出よりも後に、前記メッシュ分割手段８で前記メッシュに分割する大きさを修正してその修正結果により、前記変位影響算出手段９および前記接触面圧計算手段１０による算出を再度行わせる。
各セルＣに均一な単位分布荷重が作用したときの前記接触物体Ｍ１，Ｍ２の表面の前記各計算点Ａに対応する点の変位への影響量を算出する変位影響算出手段９と、
各セルＣにつき面圧Ｐと前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体間の接近量と初期隙間の差に釣り合、かつ各セルＣの面圧Ｐの合計が外力と釣り合うという釣り合い条件式から、各セルの面圧Ｐを算出する接触面圧計算手段１０とを含む。
なお、前記「接近量」は、前記単位分布荷重が与えられることにより互いに近づく方向に移動する距離であり、前記「初期隙間」は、前記単位分布荷重が与えられる前の距離である。
前記変位影響算出手段９および接触面圧計算手段１０は、計測点間の距離の算出に上記基準曲面Ｓ上の距離を用い、かつこの基準曲面Ｓ上の位置の違いによる法線方向の修正を行って面圧を算出する。すなわち、影響係数の算出に基準曲面Ｓの法線方向の修正を行って面圧を算出する。
前記メッシュ修正手段１２は、前記接触面圧計算手段１０による前記各セルＣの面圧Ｐの算出よりも後に、前記メッシュ分割手段８で前記メッシュに分割する大きさを修正してその修正結果により、前記変位影響算出手段９および前記接触面圧計算手段１０による算出を再度行わせる。
この発明の表面下応力計算装置において、前記２つの接触物体Ｍ１，Ｍ２が転がり軸受の内輪とボールとである場合に、次の構成としても良い。
前記メッシュ分割手段８は、前記基準曲面Ｓを円筒面とする。
前記変位影響量算出手段９は、前記円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚ、これら３軸方向の変位を、それぞれｕ，ｖ，ｗとして、
前記外力が法線方向の力である場合は、次式（１）で表される単位集中荷重作用時の変位解析解ｗを、全てのセルについて重ね合わせることで前記影響量を計算する。
前記メッシュ分割手段８は、前記基準曲面Ｓを円筒面とする。
前記変位影響量算出手段９は、前記円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚ、これら３軸方向の変位を、それぞれｕ，ｖ，ｗとして、
前記外力が法線方向の力である場合は、次式（１）で表される単位集中荷重作用時の変位解析解ｗを、全てのセルについて重ね合わせることで前記影響量を計算する。
ただし、（１）式により与えられる変位ｗは位置（θ，ｙ）の要素が位置（θ′，ｙ′）に及ぼす変位。
この場合に、前記接触面圧計算手段１０は、前記（１）式で計算されるｚ軸方向の変位ｗを全セルについて重ね合わせる式である次式（２）を、前記釣り合い条件式の一つとして用いて各セルの面圧分布を求める。
また、前記釣り合い条件式の他の一つとして、次式（３）を用いて前記接近量αを求める。
この構成の面圧計算装置によると、この発明方法を実施して、共形接触を行う接触面における面圧を、比較的簡単に短時間で計算でき、かつ精度良く計算することができる。
この発明の面圧計算装置において、さらに表面下応力の計算機能を付加する場合は、次の表面下応力計算手段１３を設けることで、表面下応力計算装置を構成する。
前記表面下応力計算手段１３は、セルに対する荷重が法線方向の力である場合は、次式（４）〜（９）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々のセルの全域にわたり積分して、かつ全セルの影響の和をとって表面下応力を求める。
前記表面下応力計算手段１３は、セルに対する荷重が法線方向の力である場合は、次式（４）〜（９）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々のセルの全域にわたり積分して、かつ全セルの影響の和をとって表面下応力を求める。
この構成の場合、この発明の表面下応力計算方法を実施して、共形接触を行う接触面における表面下応力を、比較的簡単に短時間で計算でき、かつ精度良く計算することができる。
この発明の共形接触下の接触面圧計算方法および計算装置によると、半無限体に集中荷重が作用したときの変位の解析解を共形接触問題に適用し、また接触面圧計算時の計測点間の距離の算出に上記基準曲面上の距離を用い、かつ影響係数の算出に基準曲面の法線方向の修正を行って面圧を算出するようにしたため、比較的簡単に、したがって短時間で、共形接触における接触面圧の計算が行え、しかも精度の良い計算結果を得ることができる。特に、解析解として、法線方向荷重作用時はブジネスクの解を重ね合わせるようにした場合は、有限要素法を用いて計算した場合と同程度に精度の良い計算結果を得ることができる。特に、前記解析過程の途中で、前記メッシュに分割する大きさを修正してその修正結果により、メッシュ分割後の各過程を再度行うメッシュ修正過程を含むため、接触物体の接触面の形状によっては、接触領域の縁がギザギザになるなどして、離散化の誤差が大きくなる場合があるが、解析途中でメッシュを修正することで、離散化の誤差を低減させることができる。
この発明の共形接触下の接触面圧，表面下応力の計算方法および計算装置によると、接触面圧に加えて、さらに表面下応力についても、比較的簡単に短時間で、また精度良く、共形接触における表面下応力の計算を行うことができる。特に、メッシュ修正手段を設けたため、解析途中でメッシュを修正することにより、離散化の誤差を低減させることができる。
この発明の共形接触下の接触面圧，表面下応力の計算方法および計算装置によると、接触面圧に加えて、さらに表面下応力についても、比較的簡単に短時間で、また精度良く、共形接触における表面下応力の計算を行うことができる。特に、メッシュ修正手段を設けたため、解析途中でメッシュを修正することにより、離散化の誤差を低減させることができる。
この発明の一実施形態を図１ないし図８と共に説明する。図１は、この発明方法の各過程の流れ図を、図２はこの発明方法を実施する装置の概念構成のブロック図を、図３は接触面圧計算過程の具体例を、図４は接触モデルを、図５は基準曲面の例をそれぞれ示す。
この共形接触下の接触面圧および表面下応力計算方法は、外力が作用したときに互いに共形接触を行う２つの接触物体Ｍ１，Ｍ２（図４）の接触面における面圧等を、コンピュータを用いて計算する方法である。２つの接触物体Ｍ１，Ｍ２は、例えば、転がり軸受におけるボールと内輪とである。
この計算方法は、前記接触物体Ｍ１，Ｍ２の接触面の形状、弾性係数（ここでは、ヤング率Ｅと、ポアソン比ν）、および作用する外力Ｆを入力して記憶手段に記憶させる入力過程（Ｓ１）と、入力された値を用いて前記面圧を計算する解析過程（Ｓ２）と、計算した結果を出力する出力過程（Ｓ３）とでなる。解析過程（Ｓ２）は、順次行われる次のメッシュ分割過程（Ｔ１）、変位影響算出過程（Ｔ２）、接触面圧計算過程（Ｔ３）、および表面下応力計算過程（Ｔ４）を含む。
この計算方法は、前記接触物体Ｍ１，Ｍ２の接触面の形状、弾性係数（ここでは、ヤング率Ｅと、ポアソン比ν）、および作用する外力Ｆを入力して記憶手段に記憶させる入力過程（Ｓ１）と、入力された値を用いて前記面圧を計算する解析過程（Ｓ２）と、計算した結果を出力する出力過程（Ｓ３）とでなる。解析過程（Ｓ２）は、順次行われる次のメッシュ分割過程（Ｔ１）、変位影響算出過程（Ｔ２）、接触面圧計算過程（Ｔ３）、および表面下応力計算過程（Ｔ４）を含む。
計算装置となるコンピュータ１（図２）は、中央処理装置（ＣＰＵ）およびメモリ等の記憶手段（いずれも図示せず）を有し、このコンピュータ１のハードウェア（オペレーションシステムを含む）、およびこのコンピュータ１に実行させる接触面圧・表面下応力計算プログラム（図示せず）により、図２に概念構成で示す各手段が構成される。
入力処理手段４は、入力装置２から、接触物体Ｍ１，Ｍ２の接触面の形状、弾性係数、ポアソン比、および作用する外力Ｆを入力して入力情報記憶手段６に記憶させる手段である。ここで言う「形状」は、大きさを含む概念である。例えば、前記接触物体Ｍ１，Ｍ２の接触面の形状のデータとしては、ボールからなる接触物体Ｍ１の半径Ｒ1 と、内輪からなる接触物体Ｍ２の軌道面の溝半径Ｒ2 軸受のPCD 、接触角とを入力する。入力処理手段４は、上記の形状等の情報の他に、後述の解析手段７でメッシュ分割する枡目の各軸方向の大きさやセル個数（ＮＡ，ＮＢ）等を、入力装置２から入力して入力情報記憶手段６に記憶させる機能を備えるものとしても良い。また、入力処理手段４は、入力装置２から入力された情報に対して、解析手段７で演算するための前処理を行う機能を有するものとしても良い。入力装置２は、キーボード等であっても、また通信手段や記憶素子等であっても良い。入力情報記憶手段６は、コンピュータ１の上記記憶手段における所定の記憶領域である。出力処理手段５は、解析手段７で計算した結果を、液晶ディスプレイ，プリンタ，あるいは通信機器等からなる出力装置３に出力する手段である。
入力処理手段４は、入力装置２から、接触物体Ｍ１，Ｍ２の接触面の形状、弾性係数、ポアソン比、および作用する外力Ｆを入力して入力情報記憶手段６に記憶させる手段である。ここで言う「形状」は、大きさを含む概念である。例えば、前記接触物体Ｍ１，Ｍ２の接触面の形状のデータとしては、ボールからなる接触物体Ｍ１の半径Ｒ1 と、内輪からなる接触物体Ｍ２の軌道面の溝半径Ｒ2 軸受のPCD 、接触角とを入力する。入力処理手段４は、上記の形状等の情報の他に、後述の解析手段７でメッシュ分割する枡目の各軸方向の大きさやセル個数（ＮＡ，ＮＢ）等を、入力装置２から入力して入力情報記憶手段６に記憶させる機能を備えるものとしても良い。また、入力処理手段４は、入力装置２から入力された情報に対して、解析手段７で演算するための前処理を行う機能を有するものとしても良い。入力装置２は、キーボード等であっても、また通信手段や記憶素子等であっても良い。入力情報記憶手段６は、コンピュータ１の上記記憶手段における所定の記憶領域である。出力処理手段５は、解析手段７で計算した結果を、液晶ディスプレイ，プリンタ，あるいは通信機器等からなる出力装置３に出力する手段である。
解析手段７は、入力処理手段４で入力されて入力情報記憶手段６に記憶されたデータを用いて前記面圧および表面下応力を計算する手段であり、メッシュ分割手段８、変位影響算出手段９、接触面圧計算手段１０、および表面下応力計算手段１１を備える。解析手段７に備えられた上記各手段８〜１２は、それぞれ図１（Ｂ）の流れ図におけるメッシュ分割過程（Ｔ１）、変位影響算出過程（Ｔ２）、接触面圧計算過程（Ｔ３）、および表面下応力計算過程（Ｔ４）を実行する手段であり、その実行に必要な計算式や設定値等を記憶している。上記各手段８〜１１の具体的な内容は、次に述べる上記各過程（Ｔ１）〜（Ｔ４）の説明の通りである。
解析手段７は、上記の他に、後に説明するメッシュ修正手段１２を有している。
解析手段７は、上記の他に、後に説明するメッシュ修正手段１２を有している。
図１（Ａ）の解析過程（Ｓ２）を、同図（Ｂ）および図４〜図６と共に説明する。
この解析仮定（Ｓ２）は、法線方向加重が作用したときのハーネット（Hartnett) の手法を拡張して共形接触問題（Conformal Problem ）を解く手法である。
この解析仮定（Ｓ２）は、法線方向加重が作用したときのハーネット（Hartnett) の手法を拡張して共形接触問題（Conformal Problem ）を解く手法である。
まず、メッシュ分割過程（Ｔ１）として、図４のように、接触物体Ｍ１，Ｍ２の間に基準曲面Ｓを設定する。この基準曲面Ｓを、図５，図６のようにセルＣが縦横に並ぶメッシュに分割し、接触物体Ｍ１，Ｍ２間の距離を計算する点である計算点ＡをセルＣ上に定める。
また、上記基準曲面Ｓは、２つの接触物体Ｍ１，Ｍ２の接触面形状に近い形状の曲面とする。
接触物体Ｍ１，Ｍ２は、例えば図１１に示す転がり軸受におけるボール２３と内輪２１とである。同図の転がり軸受は、内輪２１と外輪２２の断面円弧面状の軌道面間に複数のボール２３を介在させ、ボール２３を保持する保持器２４を設けたアンギュラ玉軸受、または深溝玉軸受等の玉軸受である。
接触物体Ｍ１，Ｍ２が、上記のように転がり軸受におけるボールと内輪とである場合、基準曲面Ｓは円筒面とみなす。同図では、基準曲面Ｓは、内輪である接触物体Ｍ２の軌道面とボールである接触物体Ｍ１の外球面との間に介在する円筒面であって、その円筒中心軸は内輪である接触物体Ｍ２の軸心と平行としている。
また、上記基準曲面Ｓは、２つの接触物体Ｍ１，Ｍ２の接触面形状に近い形状の曲面とする。
接触物体Ｍ１，Ｍ２は、例えば図１１に示す転がり軸受におけるボール２３と内輪２１とである。同図の転がり軸受は、内輪２１と外輪２２の断面円弧面状の軌道面間に複数のボール２３を介在させ、ボール２３を保持する保持器２４を設けたアンギュラ玉軸受、または深溝玉軸受等の玉軸受である。
接触物体Ｍ１，Ｍ２が、上記のように転がり軸受におけるボールと内輪とである場合、基準曲面Ｓは円筒面とみなす。同図では、基準曲面Ｓは、内輪である接触物体Ｍ２の軌道面とボールである接触物体Ｍ１の外球面との間に介在する円筒面であって、その円筒中心軸は内輪である接触物体Ｍ２の軸心と平行としている。
なお、基準曲面Ｓについて、接触面形状の正確な形状は本来わからないが、ボールと内輪の接触では接触面形状は円筒面の一部に近いとみなせる（解析上もこの基準面を採用している）。また、凸球面と凹球面の接触では接触面形状はほぼ球面とみなせる。従来の「接触面が平坦と見なせる」場合は平面である。複雑な曲面同士の接触では上記のような簡単な曲面で表現できない場合が出てくる。基準面形状をどのようにするかは解析精度に影響してくるため、できるだけ接触面形状に近い方が良いが、複雑にすると計算が難しくなるというデメリットがある。従って「基準曲面Ｓ」は解析する人がケースバイケースで適宜決定していく。
計算点Ａは、接触物体Ｍ１，Ｍ２間の距離を離散的な点で評価するものとするためにセルＣの位置を代表させる点である。計算点Ａは、セルＣ内の適宜の位置に設定すれば良いが、ここではセルＣの中心にとる。
メッシュ分割する枡目の大きさ、つまりセルＣの大きさは、適宜設定すればよいが、入力処理手段４からの入力等によって自由に設定可能としても良い。また、基準曲面Ｓの広さは、接触物体Ｍ１，Ｍ２の接触領域を全て含むように、接触領域より広めに適宜設定する。ここでは、セルＣがＸ軸方向にＮＡ個、Ｙ軸方向にＮＢ個並ぶ基準曲面Ｓとする。
メッシュ分割する枡目の大きさ、つまりセルＣの大きさは、適宜設定すればよいが、入力処理手段４からの入力等によって自由に設定可能としても良い。また、基準曲面Ｓの広さは、接触物体Ｍ１，Ｍ２の接触領域を全て含むように、接触領域より広めに適宜設定する。ここでは、セルＣがＸ軸方向にＮＡ個、Ｙ軸方向にＮＢ個並ぶ基準曲面Ｓとする。
変位影響算出過程（Ｔ２）では、メッシュ分割過程（Ｔ１）の後、各セルＣに均一な分布荷重Ｐが作用したときの接触物体Ｍ１，Ｍ２の表面の前記各計算点Ａに対応する点の変位への影響量を算出する。この場合に、ここでは、個々のセルＣには全体に均一な面圧Ｐが作用すると仮定するが、１次関数または２次関数等の所定の関数に従いセルＣ内の位置によって変わる面圧が作用すると仮定しても良い。
このとき、互いに直交する３軸方向の座標位置をｘ，ｙ，ｚ、これら３軸方向の変位を、それぞれｕ，ｖ，ｗとして計算する。なお、ｘ，ｙ軸方向の変位ｕ，ｖは、後に、単位集中荷重作用時の変位解析解（Boussinesqの解）を説明するために記載したが、この発明では、これらｘ，ｙ軸方向の変位ｕ，ｖを用いず、ｚ軸方向の変位ｗのみを用いて計算する。
２つの接触物体Ｍ１，Ｍ２が、上記のように転がり軸受のボールと内輪であって、基準面Ｓを円筒面とする場合、この円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚ、これら３軸方向の変位を、それぞれｕ，ｖ，ｗとする。なお、基準面Ｓ上のｘ座標位置は、角度θにより表現される極座標Ｒo θで表される。
２つの接触物体Ｍ１，Ｍ２が、上記のように転がり軸受のボールと内輪であって、基準面Ｓを円筒面とする場合、この円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚ、これら３軸方向の変位を、それぞれｕ，ｖ，ｗとする。なお、基準面Ｓ上のｘ座標位置は、角度θにより表現される極座標Ｒo θで表される。
前記変位影響算出過程（Ｔ２）では、単位集中荷重作用時の変位解析解（Boussinesqの解）を１つの要素全域にわたり積分して、要素１つの面圧分布によって生じる２接触物体Ｍ１，Ｍ２の表面上の点（上記計測点Ａに対応）の変位量を求める。
すなわち、位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’）に及ぼす単位面圧負荷時の変位を、次の（１）式で与える。このとき、基準曲面の法線方向が位置により変わるためCOS(θ−θ’) がかかってくる。なお、ここでは、一つのセルＣを一つの要素として計算する。
ここで重要なのは、(1) 要素間のθ方向距離として基準面上の距離（Ｒ０θ−Ｒ０θ’）を用いる点、および(2) 基準面の法線方向の修正を、COS(θ−θ’) の乗算により行うこと（すなわち、影響係数の算出に基準曲面Ｓの法線方向の修正を行って面圧を算出すること）の２点である。
すなわち、位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’）に及ぼす単位面圧負荷時の変位を、次の（１）式で与える。このとき、基準曲面の法線方向が位置により変わるためCOS(θ−θ’) がかかってくる。なお、ここでは、一つのセルＣを一つの要素として計算する。
ここで重要なのは、(1) 要素間のθ方向距離として基準面上の距離（Ｒ０θ−Ｒ０θ’）を用いる点、および(2) 基準面の法線方向の修正を、COS(θ−θ’) の乗算により行うこと（すなわち、影響係数の算出に基準曲面Ｓの法線方向の修正を行って面圧を算出すること）の２点である。
なお、前記の式（１）は、単位集中荷重作用時の変位解析解（Boussinesqの解）から、後述のように導かれる。
前記接触面圧計算過程（Ｔ３）では、各セルＣにつき面圧と前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体Ｍ１，Ｍ２間の接近量と初期隙間の差に釣り合い、かつ各セルＣの面圧の合計が外力と釣り合ういう釣り合い条件式から、各セルの面圧を算出する。
接触面圧計算過程（Ｔ３）は、具体的には、次の（２）式が釣り合うように面圧分布を決め、また（３）式が釣り合うように、接近量αを求める。（２）式は、各セルＣijにつき面圧Ｐklと前記単位分布荷重による影響量ｆijklの積の合計が接触領域内部において接触物体Ｍ１，Ｍ２間の接近量と初期隙間の差Ｄijに釣り合うという条件式であり、（１）式で求められる変位ｗを全要素について重ね合わせる式である。（３）式は各セルＣijの面圧Ｐijの合計が外力Ｆと釣り合ういう条件式である。
なお、Ｃijは、ｘ軸方向にｉ番目でｙ軸方向にｊ番目にあるセルＣを示す。Ｑklは、ｘ軸方向にｋ番目でｙ軸方向にｌ番目の点Ａを示す（なお、ｉ：１〜ＮＡ、ｊ：１〜ＮＢ、ｋ：１〜ＮＡ、ｌ：１〜ＮＢ）。
なお、Ｃijは、ｘ軸方向にｉ番目でｙ軸方向にｊ番目にあるセルＣを示す。Ｑklは、ｘ軸方向にｋ番目でｙ軸方向にｌ番目の点Ａを示す（なお、ｉ：１〜ＮＡ、ｊ：１〜ＮＢ、ｋ：１〜ＮＡ、ｌ：１〜ＮＢ）。
上記の（２）式および（３）式が釣り合うように各セルＣの面圧値を決めれば、離散化された全体の面圧分布が求まる。
なお、上記のように（２）式を解く過程で、２つの接触物体Ｍ１，Ｍ２間の接近量と初期隙間の差Ｄijが正の値となるセルＣijについては、そのセルＣijの面圧は零としており、この接近量と初期隙間の差Ｄijが零または負の値となるセルＣijの並び範囲が、接触領域である。
（２）式は、解析的に積分可能であり、数値積分は不要となる。
なお、上記のように（２）式を解く過程で、２つの接触物体Ｍ１，Ｍ２間の接近量と初期隙間の差Ｄijが正の値となるセルＣijについては、そのセルＣijの面圧は零としており、この接近量と初期隙間の差Ｄijが零または負の値となるセルＣijの並び範囲が、接触領域である。
（２）式は、解析的に積分可能であり、数値積分は不要となる。
図３に、（２）式、（３）式を、コンピュータにおける接触面圧計算手段１０によって解くアルゴリズムの一例を示す。
まず、α値（：２つの接触物体Ｍ１，Ｍ２の接近量）の初期値を設定する（ステップＨ１）。
そのα値に応じて各セルＣijの接触面圧を設定する（ステップＨ２）。この場合に、ｉ＝１〜ＮＡの範囲、およびｊ＝１〜ＮＢの範囲で全てのセルＣijの接触面圧Ｐijを設定する。その設定した各セルＣijの接触面圧Ｐijで（２）式が成り立つか否かを判定し、成り立たない場合は、セルＣijの接触面圧の設定過程（ステップＨ２）に戻る（Ｈ３）。ステップＨ２では、各セルＣijの接触面圧Ｐijを再度設定し直し、その設定したセルＣijの接触面圧Ｐijで（２）式が成り立つか否かを判定し、成り立たない場合はセルＣijの接触面圧Ｐijを設定過程（Ｈ２）に戻る。このような処理を（２）式が成り立つまで繰替し、（２）式が成り立つと、ステップＨ４に進む。
ステップＨ４では、（３）式が成り立つか否かを判定し、成り立たない場合は、ステップＨ１に戻ってα値を、再度設定し直す。
この後、前記と同様にしてステップＨ２〜Ｈ４，Ｈ１の処理を繰り返し、（３）式が成り立つと、処理を終了する。
このようにして求めた、式（２）、（３）が成り立つときの、ステップＨ２で設定した各セルＣijの接触面圧Ｐijが、接触面圧計算手段１０の出力となる接触面圧の値である。 なお、セルＣijの接触面圧の最初の値は何でもよいが、できるだけ実際の面圧に近い値を設定した方が早く収束する。
まず、α値（：２つの接触物体Ｍ１，Ｍ２の接近量）の初期値を設定する（ステップＨ１）。
そのα値に応じて各セルＣijの接触面圧を設定する（ステップＨ２）。この場合に、ｉ＝１〜ＮＡの範囲、およびｊ＝１〜ＮＢの範囲で全てのセルＣijの接触面圧Ｐijを設定する。その設定した各セルＣijの接触面圧Ｐijで（２）式が成り立つか否かを判定し、成り立たない場合は、セルＣijの接触面圧の設定過程（ステップＨ２）に戻る（Ｈ３）。ステップＨ２では、各セルＣijの接触面圧Ｐijを再度設定し直し、その設定したセルＣijの接触面圧Ｐijで（２）式が成り立つか否かを判定し、成り立たない場合はセルＣijの接触面圧Ｐijを設定過程（Ｈ２）に戻る。このような処理を（２）式が成り立つまで繰替し、（２）式が成り立つと、ステップＨ４に進む。
ステップＨ４では、（３）式が成り立つか否かを判定し、成り立たない場合は、ステップＨ１に戻ってα値を、再度設定し直す。
この後、前記と同様にしてステップＨ２〜Ｈ４，Ｈ１の処理を繰り返し、（３）式が成り立つと、処理を終了する。
このようにして求めた、式（２）、（３）が成り立つときの、ステップＨ２で設定した各セルＣijの接触面圧Ｐijが、接触面圧計算手段１０の出力となる接触面圧の値である。 なお、セルＣijの接触面圧の最初の値は何でもよいが、できるだけ実際の面圧に近い値を設定した方が早く収束する。
この接触面圧計算方法によると、このように、半無限体に集中荷重が作用したときの変位の解析解を共形接触問題に適用し、また接触面圧計算時の計測点間の距離の算出に上記基準曲面Ｃ上の距離を用い、かつ影響係数の算出に基準曲面Ｃの法線方向の修正を行って面圧を算出するようにしたため、比較的簡単に、したがって短時間で、共形接触における接触面圧の計算を行うことができる。
また、有限要素法を用いて接触面圧を計算した場合と同程度に精度の良い計算結果が得られることが、この実施形態の方法と有限要素法との計算結果を比較することで、確認できた。
この実施形態は、このように、従来では共形接触には半無限体が適用できないと考えられていた共形接触問題に、ブジネスクの半無限体近似の解を適用し、良好な面圧の解析結果を得ることができる。
なお、この発明方法は、非共形接触の場合にも適用可能であるが、非共形接触の場合はより簡単な計算で面圧を求めることができるため、共形接触の場合に効果的となる。
また、有限要素法を用いて接触面圧を計算した場合と同程度に精度の良い計算結果が得られることが、この実施形態の方法と有限要素法との計算結果を比較することで、確認できた。
この実施形態は、このように、従来では共形接触には半無限体が適用できないと考えられていた共形接触問題に、ブジネスクの半無限体近似の解を適用し、良好な面圧の解析結果を得ることができる。
なお、この発明方法は、非共形接触の場合にも適用可能であるが、非共形接触の場合はより簡単な計算で面圧を求めることができるため、共形接触の場合に効果的となる。
表面下応力計算過程（Ｔ４）では、上記のように求めた面圧分布に対して、次式（４）〜（９）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々のセルＣの全域にわたり積分して、かつ全セルＣの影響の和をとって表面下応力が求まる。
セルＣに対する荷重が法線方向の力である場合は、次のように計算する。前記２つの接触物体Ｍ１，Ｍ２が転がり軸受のボールと内輪とである場合、前記メッシュ分割過程（Ｔ１）では、前記のように前記基準曲面Ｓを円筒面とする。前記変位影響量算出過程（Ｔ２）では、前記円筒面からなる基準曲面Ｓ上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をθ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚとして、
位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’，ｚ）に及ぼす単位面圧負荷時の剪断応力τzyを、次式（１０）で求める。
位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’，ｚ）に及ぼす単位面圧負荷時の剪断応力τzyを、次式（１０）で求める。
この場合も、積分は解析的に求まるため、数値積分は不要である。
このように表面下応力を計算するについても、計算が比較的簡単で短時間で行え、また精度の良い表面下応力の計算が行える。
なお、上記実施形態において、接触物体Ｍ１，Ｍ２の形状によっては、接触領域の縁がギザギザにあることがあるが、そのような場合は、解析途中で、例えば接触面圧計算過程（Ｔ３）の後、メッシュ分割の修正、つまりセルＣの各軸ｘ，ｙ方向の大きさを修正してその修正結果により、メッシュ分割後の各過程を再度行わせるメッシュ修正過程（Ｔ５）を含むようにしても良い。
図２のメッシュ修正手段１２は、このようなメッシュ修正過程（Ｔ５）を実行する手段である。
図２のメッシュ修正手段１２は、このようなメッシュ修正過程（Ｔ５）を実行する手段である。
接触物体Ｍ１，Ｍ２が、転がり軸受におけるボールと内輪である場合における、この実施形態、前述の提案例（特願２００７−０２６５７７号）、および有限要素法（ＦＥＭ）の計算結果を、図９に示す。３者とも良く一致しており、どれも使用可能レベルである。一方、従来手法（接触面が平坦という前提の解析）では、図１０のようであり、正しくない結果となる。
前述の提案例の手法では、３方向の変位を考慮していたため、計算時間が長くかかっていたが（ＦＥＭよりは早い）が、この実施形態は、これを１方向とすることができた。また、全順の提案例は数値積分で行っていたが、この実施形態では解析積分に置き換えることができ、計算の簡素化と高速化が達成できた。
計算の所要時間の例を示すと、前述の提案例の手法では１９７．０８秒であったのに対して、実施形態の手法では１１２．４５秒となり、従来手法の５７％であった。
計算の所要時間の例を示すと、前述の提案例の手法では１９７．０８秒であったのに対して、実施形態の手法では１１２．４５秒となり、従来手法の５７％であった。
つぎに、上記の式（１）とブジネスク（Boussinesq）の解との関係を説明する。
作用荷重が法線方向である場合の単位集中荷重作用時の変位解析解であるブジネスク（Boussinesq）の解は、次式（１１）〜（１３）によって現される。
なお、図７はブジネスクの解の座標系を示す。同図のような半無限体平面の原点Ｏに単位荷重が作用するとき、任意の点であるＱ点のｘ，ｙ，ｚ方向の変位を、それぞれｕ，ｖ，ｗとする。
なお、図７はブジネスクの解の座標系を示す。同図のような半無限体平面の原点Ｏに単位荷重が作用するとき、任意の点であるＱ点のｘ，ｙ，ｚ方向の変位を、それぞれｕ，ｖ，ｗとする。
この式から、次のようにして、上記の式（１）が導かれる。
これを要素eが点Qに及ぼす変位として相対距離を用いて積分したものが(1)式である。本手法では(1)式のｘには基準面上の距離を用いるので、ｘの代わりにR0θを用いる。なお、積分に当たっては以下の変数変換
ｘ＝Ｒ0θ ⇒ ｄｘ＝Ｒ0ｄθ
を用いて（１）式が得られる。
ｘ＝Ｒ0θ ⇒ ｄｘ＝Ｒ0ｄθ
を用いて（１）式が得られる。
１…コンピュータ
４…入力処理手段
５…出力処理手段
６…入力情報記憶手段
７…解析手段
８…メッシュ分割手段
９…変位影響算出手段
１０…接触面圧計算手段
１１…表面下応力計算手段
Ｍ１，Ｍ２…接触物体
Ｓ…基準曲面
Ｃ…セル
４…入力処理手段
５…出力処理手段
６…入力情報記憶手段
７…解析手段
８…メッシュ分割手段
９…変位影響算出手段
１０…接触面圧計算手段
１１…表面下応力計算手段
Ｍ１，Ｍ２…接触物体
Ｓ…基準曲面
Ｃ…セル
Claims (10)
- 外力が作用したときに互いに共形接触を行う２つの接触物体の接触面における面圧を、コンピュータを用いて計算する方法であって、
前記接触物体の接触面の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力過程と、入力された値を用いて前記面圧を計算する解析過程と、計算した結果を出力する出力過程とでなり、
前記解析過程として、
前記２つの接触物体の互いに接する面の間に、基準曲面を設定してこの基準曲面をセルが縦横に並ぶメッシュに分割し、接触物体間の距離を計算する点である計算点をセル上に定めるメッシュ分割過程と、
各セルに単位分布荷重が作用したときの前記２つの接触物体の表面の前記各計算点に対応する点に生じる変位である影響量を算出する変位影響算出過程と、
各セルにつき面圧と前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体間の、前記単位分布荷重が与えられることにより互いに近づく方向に移動する距離である接近量と、前記単位分布荷重が与えられる前の距離である初期隙間の差に釣り合い、かつ各セルの面圧の合計が外力と釣り合うという釣り合い条件式から、各セルの面圧を算出する接触面圧計算過程と、
この解析過程の途中で、前記メッシュに分割する大きさを修正してその修正結果により、メッシュ分割後の前記各過程を再度行うメッシュ修正過程とを含み、
前記変位影響算出過程および接触面圧計算過程では、計測点間の距離の算出に上記基準曲面上の距離を用い、かつこの基準曲面上の位置の違いによる法線方向の修正を行って面圧を算出する、
ことを特徴とする共形接触下の接触面圧計算方法。 - 請求項１において、前記２つの接触物体が転がり軸受の内輪とボールとであって、
前記メッシュ分割過程では、前記基準曲面を円筒面とし、
前記変位影響量算出過程では、前記円筒面からなる基準曲面上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚ、前記ｚ軸方向の変位をｗとして、
前記外力が法線方向の力である場合は、次式（１）で表される単位集中荷重作用時の変位解析解により前記影響量を計算する、
ただし、（１）式により与えられる変位ｗは位置（θ，ｙ）の要素が位置（θ′，ｙ′）に及ぼす変位、
共形接触下の接触面圧計算方法。 - 請求項２において、前記接触面圧計算過程では、前記（１）式で計算されるｚ軸方向の変位をｗを全セルについて重ね合わせる式である次式（２）を、前記釣り合い条件式の一つとして用いて各セルの面圧分布を求め、
- 請求項１ないし請求項３のいずれか１項において、前記メッシュ修正過程は、前記接触面圧計算過程の後に行う共形接触下の接触面圧計算方法。
- 外力が作用したときに互いに共形接触を行う２つの接触物体の接触面における面圧および表面下応力を、コンピュータを用いて計算する方法であって、
請求項１ないし請求項４のいずれか１項に記載の接触面圧の計算方法と、この計算方法で各セルの面圧を計算した後に、表面下応力を計算する過程を含み、
この表面下応力計算過程として、
前記外力が法線方向の力である場合は、次式（４）〜（９）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々のセルの全域にわたり積分して、かつ全セルの影響の和をとって表面下応力を求める、
- 請求項５において、前記２つの接触物体が転がり軸受の内輪とボールとであって、
前記メッシュ分割過程では、前記基準曲面を円筒面とし、
前記変位影響量算出過程では、前記円筒面からなる基準曲面上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をθ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚとして、
位置（θ，ｙ）の要素ｅが位置（θ’，ｙ’，ｚ）に及ぼす単位面圧負荷時の剪断応力τzyを、次式（１０）で求める、
- 外力が作用したときに互いに共形接触を行う２つの接触物体の接触面における面圧を計算する装置であって、
前記接触物体の接触面の形状、弾性係数、ポアソン比、および作用する外力を入力して記憶手段に記憶させる入力処理手段と、入力された値を用いて前記面圧を計算する解析手段と、計算した結果を出力する出力処理手段とでなり、
前記解析手段として、
前記２つの接触物体の互いに接する面の間に、基準曲面を設定してこの基準曲面をセルが縦横に並ぶメッシュに分割し、接触物体間の距離を計算する点である計算点をセル上に定めるメッシュ分割手段と、
各セルに単位分布荷重が作用したときの前記２つの接触物体の表面の前記各計算点に生じる変位である影響量を算出する変位影響算出手段と、
各セルにつき面圧と前記単位分布荷重による影響量の積の合計が接触領域内部において接触物体間の、前記単位分布荷重が与えられることにより互いに近づく方向に移動する距離である接近量と、前記単位分布荷重が与えられる前の距離である初期隙間の差に釣り合い、かつ各セルの面圧の合計が外力と釣り合うという釣り合い条件式から、各セルの面圧を算出する接触面圧計算手段と、
この接触面圧計算手段による前記各セルの面圧の算出よりも後に、前記メッシュ分割手段で前記メッシュに分割する大きさを修正してその修正結果により、前記変位影響算出手段および前記接触面圧計算手段による算出を再度行わせるメッシュ修正手段とを含み、
前記変位影響算出手段および接触面圧計算手段は、計測点間の距離の算出に上記基準曲面上の距離を用い、かつこの基準曲面上の位置の違いによる法線方向の修正を行って面圧を算出する、
ことを特徴とする共形接触下の接触面圧計算装置。 - 請求項７において、前記２つの接触物体が転がり軸受の内輪とボールとであって、
前記メッシュ分割手段は、前記基準曲面を円筒面とし、
前記変位影響量算出手段は、前記円筒面からなる基準曲面上の、円筒中心軸に垂直な平面で切断した断面の円弧に沿う方向の座標位置をｘ、円筒中心軸に沿う方向の座標位置をｙ、基準曲面上の法線方向に離れる座標位置をｚ、前記ｚ方向の変位をｗとして、
前記外力が法線方向の力である場合は、次式（１）で表される単位集中荷重作用時の変位解析解により前記影響量を計算する、
ただし、（１）式により与えられる変位ｗは位置（θ，ｙ）の要素が位置（θ′，ｙ′）に及ぼす変位、
共形接触下の接触面圧計算装置。 - 請求項８において、前記接触面圧計算手段は、前記（１）式で計算されるｚ軸方向の変位ｗを全セルについて重ね合わせる式である次式（２）を、前記釣り合い条件式の一つとして用いて各セルの面圧分布を求め、
- 外力が作用したときに互いに共形接触を行う２つの接触物体の接触面における面圧および表面下応力を計算する装置であって、
請求項７ないし請求項９のいずれか１項に記載の接触面圧の計算装置と、表面下応力計算手段とを含み、
この表面下応力計算手段は、
前記外力が法線方向の力である場合は、次式（４）〜（９）で与えられる単位集中荷重作用時の表面下の応力を用いて、個々のセルの全域にわたり積分して、かつ全セルの影響の和をとって表面下応力を求める、
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2007277066A JP5026225B2 (ja) | 2007-10-25 | 2007-10-25 | 共形接触下の接触面圧および表面下応力計算方法・計算装置 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2007277066A JP5026225B2 (ja) | 2007-10-25 | 2007-10-25 | 共形接触下の接触面圧および表面下応力計算方法・計算装置 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2009104492A JP2009104492A (ja) | 2009-05-14 |
| JP5026225B2 true JP5026225B2 (ja) | 2012-09-12 |
Family
ID=40706084
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2007277066A Active JP5026225B2 (ja) | 2007-10-25 | 2007-10-25 | 共形接触下の接触面圧および表面下応力計算方法・計算装置 |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP5026225B2 (ja) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105893760A (zh) * | 2016-04-01 | 2016-08-24 | 北京航空航天大学 | 基于广义磨损量的刀具磨损评价和预测方法及系统 |
| DE10123512B4 (de) | 2000-05-18 | 2023-01-26 | Thk Co., Ltd. | Kugellager und Verfahren zu seiner Herstellung |
Families Citing this family (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102278966B (zh) * | 2010-06-13 | 2015-07-15 | 罗伯特·博世有限公司 | 磨损量确定方法、滑动接触件制造方法及滑动接触件 |
| JP6085748B2 (ja) * | 2012-05-24 | 2017-03-01 | 長谷 陽夫 | ２物体間の接触面圧計算方法および２物体間の接触面圧計算コンピュータプログラム |
| CN113607317B (zh) * | 2021-08-04 | 2022-08-16 | 大连理工大学 | 一种滚道接触应力的间接测量方法及系统 |
-
2007
- 2007-10-25 JP JP2007277066A patent/JP5026225B2/ja active Active
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| DE10123512B4 (de) | 2000-05-18 | 2023-01-26 | Thk Co., Ltd. | Kugellager und Verfahren zu seiner Herstellung |
| CN105893760A (zh) * | 2016-04-01 | 2016-08-24 | 北京航空航天大学 | 基于广义磨损量的刀具磨损评价和预测方法及系统 |
| CN105893760B (zh) * | 2016-04-01 | 2018-11-16 | 北京航空航天大学 | 基于广义磨损量的刀具磨损评价和预测方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| JP2009104492A (ja) | 2009-05-14 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Le et al. | Surface corrections for peridynamic models in elasticity and fracture | |
| Mirzaei et al. | Direct meshless local Petrov–Galerkin (DMLPG) method: a generalized MLS approximation | |
| Chin et al. | Modeling crack discontinuities without element‐partitioning in the extended finite element method | |
| Kania et al. | A catalogue capacity of slewing bearings | |
| JP5026225B2 (ja) | 共形接触下の接触面圧および表面下応力計算方法・計算装置 | |
| Niu et al. | Topology optimization of continuum structures for the uniformity of contact pressures | |
| Chang | Dissipative, noniterative integration algorithms with unconditional stability for mildly nonlinear structural dynamic problems | |
| JP5046671B2 (ja) | 共形接触下の接触面圧および表面下応力計算方法 | |
| Zhang et al. | On the locking free isogeometric formulations for 3-D curved Timoshenko beams | |
| Koop et al. | Passive dynamic biped walking—Part II: stability analysis of the passive dynamic gait | |
| Boffy et al. | Multigrid solution of the 3D stress field in strongly heterogeneous materials | |
| Dunaj et al. | Static condensation in modeling roller guides with preload | |
| Acosta et al. | Composite material point method (CMPM) to improve stress recovery for quasi-static problems | |
| Kim et al. | Finite element analysis for contact problems | |
| Bureick et al. | Fast converging elitist genetic algorithm for knot adjustment in B-spline curve approximation | |
| Higa et al. | Wheelset-rail mechanical model for a steady-state dynamic condition and prediction of rolling contact fatigue locci | |
| Haug | An index 0 Differential-Algebraic equation formulation for multibody dynamics: Holonomic constraints | |
| Cao et al. | Key research on computer aided tolerancing | |
| Blanco-Lorenzo et al. | Approximating the influence coefficients of non-planar elastic solids for conformal contact analysis | |
| CN109766637B (zh) | 基于Krigng代理模型的桥式起重机结构可靠性优化方法 | |
| Hartmann et al. | On accuracy and efficiency of constrained reinitialization | |
| Malikov et al. | Analytical and numerical calculation of the force and power requirements for air bending of structured sheet metals | |
| Sokołowski et al. | Topological derivatives for optimization of plane elasticity contact problems | |
| Wong et al. | Study of the discrete shear gap technique in timoshenko beam elements | |
| Doronin et al. | Problems on comparing analytical and numerical estimations of stressed-deformed state of structure elements |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20100928 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20120221 |
|
| A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20120420 |
|
| TRDD | Decision of grant or rejection written | ||
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20120619 |
|
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20120620 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20150629 Year of fee payment: 3 |
|
| R150 | Certificate of patent or registration of utility model |
Ref document number: 5026225 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |