(57)【特許請求の範囲】
【請求項１】
き裂先端を原点とし、き裂の上下面を各々θ＝±πとする極座標系（ｒ，θ）において、
き裂の周辺に配置した特異要素を用いて有限要素解析を行い、特異要素上の節点のｘ、ｙ方向各々の変位Ｕ（ｒ，θ）、Ｖ（ｒ，θ）を求め、
前記節点変位を用いて、以下の式（１）により所要の応力拡大係数を求め、
【数１】
0003564538.tif 000019
前記節点変位を用いて、以下の式（２）により誤差指標を求める一連の処理を特異要素数の異なる条件で行い、
【数２】
0003564538.tif 000020
実際的なθ方向の分割を行ったときの特異要素数の条件で、前記式（１）により求めた所要の応力拡大係数に含まれる誤差を前記式(２)により求めた誤差指標により評価する、
ことを特徴とする応力拡大係数の誤差評価方法。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
この発明は、応力拡大係数の誤差評価方法に係り、特に、機械構造物、土木構造物などのインフラストラクチャに欠陥が発見された場合の構造健全性評価に適用可能な応力拡大係数誤差評価指標に関するものである。
【０００２】
【従来の技術】
任意の荷重を受ける構造物中のき裂の健全性評価を行うために、有限要素解析を行い、その結果をもとに応力拡大係数（Ｋ値）を評価することが広く行われている。これは、Ｋ値の解が存在しない場合に特に有効な手法となるが、その実用にあたっては、得られたＫ値解の精度判定が重要な課題となる。
【０００３】
有限要素解析にて、き裂を含む構造を扱う場合、き裂先端における応力の特異性を表現・評価するために多くの手法が考案されている。この中で有力なものの一つがＢａｒｓｏｕｍとＨｅｎｓｈｅｌｌ、Ｓｈａｗが独立に提案した特異要素を用いて特異性を表し、Ｔｒａｃｅｙの式によりＫ値を評価する手法（変位法：Ｄｉｓｐｌａｃｅｍｅｎｔ  ＣｏｒｒｅｌａｔｉｏｎＴｅｃｈｎｉｑｕｅ，以下ＤＣＴ）である。この手法の特徴は、比較的粗い要素分割で実用的に十分な精度のＫ値を得ることが期待できる点にあり、過去において多くの研究者がこの特異要素を用いたＫ値評価において、Ｋ値精度を確保するための特異要素寸法選定に関する研究を行ってきた。しかしながら、現在においては要素寸法のみならず、荷重条件もＫ値に影響を及ぼすことが指摘され、あらゆる条件を満足する特異要素の最適寸法は存在しないとの見解に至っている。
【０００４】
一般に、有限要素解析によるＫ値の解は、要素数を増加させることによりその精度は向上する。しかし、工学問題である以上（特に、特異要素を用いるメリットを最大限生かすためにも）少ない要素数で十分な精度の解が得られることが望ましく、その実現のためには、一度の解析を通じてＫ値の解と共にその誤差の程度を見積もることができ、それを踏まえた補正や得られた解の実用的観点からの適否を判断することができれば都合がよい。
【０００５】
一般に、弾性問題における有限要素解の誤差評価においてはＺｉｅｎｋｉｅｗｉｃｚ、Ｚｈｕらが提案したエネルギノルムを用いての誤差指標の考え方を適用することが考えられる。しかしながら、この誤差指標は、Ｋ値の次元を持たないため、この指標が小さくなる場合にＫ値の誤差も小さくなることは期待できても、そのＫ値の誤差の程度を知ることはできない。
【０００６】
その結果、数回の解析によりその収束をみるという性格のものであったため、解析に習熟・労力を要し、結果的に解析に要するコストは大きなものとならざるをえなかった。
【０００７】
【発明が解決しようとする課題】
上述したように、一般に、有限要素解析によるＫ値の解は、要素数を増加させることによりその精度は向上する。しかしながら、工学問題である以上少ない要素数で十分な精度の解が得られることが望ましく、その実現のためには、一度の解析を通じてＫ値の解と共にその誤差の程度を見積もることができ、それを踏まえた補正や得られた解の実用的観点からの適否を判断できることが望まれている。
【０００８】
そこで、この発明は、上述した問題点に鑑みなされたものであって、その目的は、ただ一度の解析により求めた応力拡大係数中の誤差の程度を評価し、結果的に、精度の高い応力拡大係数を短時間で算出することができる応力拡大係数の誤差評価方法を提供することにある。
【０００９】
【課題を解決するための手段】
上記課題を解決し目的を達成するために、
請求項１に記載の応力拡大係数の誤差評価方法は、
き裂先端を原点とし、き裂の上下面を各々θ＝±πとする極座標系（ｒ，θ）において、
き裂の周辺に配置した特異要素を用いて有限要素解析を行い、特異要素上の節点のｘ、ｙ方向各々の変位Ｕ（ｒ，θ）、Ｖ（ｒ，θ）を求め、
前記節点変位を用いて、以下の式（１）により所要の応力拡大係数を求め、
【数３】
0003564538.tif 000002
前記節点変位を用いて、以下の式（２）により誤差指標を求める一連の処理を特異要素数の異なる条件で行い、
【数４】
0003564538.tif 000003
実際的なθ方向の分割を行ったときの特異要素数の条件で、前記式（１）により求めた所要の応力拡大係数に含まれる誤差を前記式(２)により求めた誤差指標により評価する、
ことを特徴とする。
【００１０】
この発明の応力拡大係数の誤差評価方法によれば、き裂先端変位場の関数形が既知であり、特異要素についてはこの変位の関数形を一部表現しうること、また応力拡大係数をき裂先端要素の変位のみから評価しうることに着目し、応力拡大係数の次元を有する誤差指標（式（２））を提案するものである。これにより、ただ一度の解析を通じて所要の応力拡大係数ＫＤＣＴと共にその誤差ΔＫＤＣＥを見積もることができ、それを踏まえた補正や得られた解の実用的観点からの適否を判断することができる。したがって、短時間で応力拡大係数を求めるとともにその誤差の程度を評価することが可能であり、精度の高い応力拡大係数を算出することができる。
【００１１】
この発明に係る応力拡大係数の誤差評価方法は、大型汎用構造解析ソフトウエアに組み込むことにより有用である。
【００１２】
例えば、この発明は、
き裂先端を原点とし、き裂の上下面を各々θ＝±πとする極座標系（ｒ，θ）において、
き裂の周辺に配置した特異要素を用いて有限要素解析を行い、特異要素上の節点のｘ、ｙ方向各々の変位Ｕ（ｒ，θ）、Ｖ（ｒ，θ）を求め、
前記節点変位を用いて、以下の式（１）により所要の応力拡大係数を求め、
【数５】
0003564538.tif 000004
前記節点変位を用いて、以下の式（２）により前記式（１）によって求められた所要の応力拡大係数に含まれる誤差を評価する、
【数６】
0003564538.tif 000005
ことを特徴とする応力拡大係数誤差評価プログラムを記憶した記憶媒体にも適用可能である。
【００１３】
また、この発明は、
応力拡大係数誤差評価プログラムを記憶した記憶媒体と、
応力拡大係数の誤差を評価するための所定条件の入力を受け付ける入力手段と、
前記入力手段を介して入力された所定条件に基づいて、前記記憶媒体に記憶した前記応力拡大係数誤差評価プログラムを実行する制御手段と、を備えたシステムであって、
前記制御手段は、前記応力拡大係数誤差評価プログラムに基づいて、
き裂先端を原点とし、き裂の上下面を各々θ＝±πとする極座標系（ｒ，θ）において、
き裂の周辺に配置した特異要素を用いて有限要素解析を行い、特異要素上の節点のｘ、ｙ方向各々の変位Ｕ（ｒ，θ）、Ｖ（ｒ，θ）を求め、
前記節点変位を用いて、以下の式（１）により所要の応力拡大係数を求め、
【数７】
0003564538.tif 000006
前記節点変位を用いて、以下の式（２）により前記式（１）によって求められた所要の応力拡大係数に含まれる誤差を評価する、
【数８】
0003564538.tif 000007
ことを特徴とする応力拡大係数誤差評価プログラムを実行するシステムにも適用可能である。
【００１４】
【発明の実施の形態】
以下、この発明の応力拡大係数の誤差評価方法の一実施の形態について図面を参照して説明する。
【００１５】
まず、この発明の応力拡大係数（Ｋ値）の誤差評価方法に適用される誤差指標ΔＫＤＣＥ、すなわち、ＤＣＥ  Ｉｎｄｅｘ（Ｄｉｓｐｌａｃｅｍｅｎｔ  Ｃｏｒｒｅｌａｔｉｏｎ  Ｅｒｒｏｒ  Ｉｎｄｅｘ）の導出過程について説明する。
【００１６】
まず、図１に示すように、き裂２０の先端を原点、き裂２０の上下面を各々θ＝±πとする極座標系（ｒ，θ）を考える（ここでは、三角形特異要素数８個の場合を示し、図中の黒丸は要素上の節点を示す）。このとき、ｘ、ｙ方向各々の変位ｕ、ｖの一般解は、
【数９】
0003564538.tif 000008
として与えられる。ここで、Ｇはせん断弾性係数、関数ｆＩｕｎ、ｆＩｖｎ、ｆＩＩｕｎ、ｆＩＩｖｎは次式により与えられ、添字のＩ、ＩＩは各々モードＩ、ＩＩに、そしてｕ、ｖは各々ｘ、ｙ方向の変位に対応する。モードＩとは、き裂２０に対してｙ方向に沿った引っ張り方向の変位に対応し、モードＩＩとは、き裂２０に対してｘ方向に沿ったせん断方向の変位に対応する。κは平面応力・平面ひずみ両状態に対し，き裂先端の変位を共通の式で表すときに用いる係数に対応し、νをポアッソン比として平面ひずみ場にて（３−４ν）、平面応力場に対して（３−ν）／（１＋ν）となる量である。
【００１７】
【数１０】
0003564538.tif 000009
【００１８】
ここでｆＩｖｎ、ｆＩＩｕｎが奇関数、ｆＩｕｎ、ｆＩＩｖｎが偶関数であることより、ｘ軸をはさんで上下対称点の相対変位を考えると、それらは単一のモードで表わされる（モード分離）。
【００１９】
【数１１】
0003564538.tif 000010
【００２０】
一方、特異要素を用いる有限要素解析では、特異要素の一辺上の変位（Ｕ（ｒ，θ）、Ｖ（ｒ，θ））を、解析の結果得られたその辺θ＝θ上のｑｕａｒｔｅｒ  ｐｏｉｎｔ（ｒ＝Ｌ／４）、ｅｎｄ  ｐｏｉｎｔ（ｒ＝Ｌ）の節点変位を用いることにより次式として表わされる。
【００２１】
【数１２】
0003564538.tif 000011
【００２２】
特異要素についてもＵ＊（ｒ，θ）≡Ｕ（ｒ，θ）−Ｕ（ｒ，−θ）、Ｖ＊（ｒ，θ）≡Ｖ（ｒ，θ）−Ｖ（ｒ，−θ）を導入することにより、式（３）にならったモード分離を行い次式を得ることができる。
【００２３】
【数１３】
0003564538.tif 000012
【００２４】
この特異要素を用いた計算結果により、Ｋ値を評価する手法は種々あるが、よく用いられるものの一つにＴｒａｃｅｙの評価式がある。このＫ値評価式では、き裂面上（θ＝π）で式（５）と変位の一般解式（３）にて級数の第２項までをとる場合と対応づけることにより、Ｋ値すなわちＫＤＣＴを評価している。具体的には、ｆＩｖ１（π）＝ｆＩＩｕ１（π）＝κ＋１、ＡＩ１＝ＫＩＤＣＴ／（２π）１／２、ＡＩＩ１＝−ＫＩＩＤＣＴ／（２π）１／２、Ｇ’≡（２π／Ｌ）１／２Ｇ／（１＋κ）、とおくことにより、次式が得られる。
【００２５】
【数１４】
0003564538.tif 000013
【００２６】
この式（６）は、具体的には、次式（６）’に示すように、図１に示す特異要素のき裂２０’の上下面における４点ｂ、ｃ、ｄ、ｅの変位で表すことができる。このとき、点ｂ及び点ｄは、上述したｑｕａｒｔｅｒ  ｐｏｉｎｔに相当し、点ｃ及び点ｅは、上述したｅｎｄ  ｐｏｉｎｔに相当し、それぞれの座標は、ｂ：（Ｌ／４，π）、ｃ：（Ｌ，π）、ｄ：（Ｌ／４，−π）、ｅ：（Ｌ，−π）である。
【００２７】
【数１５】
0003564538.tif 000014
【００２８】
このとき、ＫＤＣＴではθ＝πに対し、Ｋ値を評価しているが、特異要素は、変位の解析解が示すべきθ特性が保証されているわけではないので、他の特異要素の辺上変位からＫ値を評価する場合のＫ値は、一般に式（６）で求めたＫＤＣＴと異なる。
【００２９】
さて、き裂先端の十分に小さい領域を選んだ結果、Ｏ（ｒ３／２）以上の高次項の影響が十分小さいとの条件が満足される場合には、式（３）にて級数の第２項まで考えることにより、この領域の変位を十分な精度で表わすことができる。ここで、き裂面上に着目し、ｆＩｖ２（π）＝０、ｆＩＩｕ２（π）＝０であることより、式（３）は真のＫ値、すなわちＫＩ、ＫＩＩを用いて次のように表わされる。
【００３０】
【数１６】
0003564538.tif 000015
【００３１】
一方、特異要素のき裂面上のＵ＊（ｒ，π）、Ｖ＊（ｒ，π）は、式（６）にて定義したＫＩＤＣＴ、ＫＩＩＤＣＴを用いて次式で表わされる。
【００３２】
【数１７】
0003564538.tif 000016
【００３３】
ここで、有限要素解析（ＦＥＡ）による各節点の変位は、変位関数に含まれる未定係数（特異要素の場合ｒ１／２、ｒの項の未定係数）を、エネルギの意味で最適になるように決定して与えられるものであり、真の変位の解が特異要素の変位関数に含まれていない限り、式（８）のｒ／Ｌの係数は零とはならない。
【００３４】
しかしながら、ここで用いる要素は、適合要素であることから、要素を限りなく小さくしていくとき（ｒ方向だけでなくθ方向の分割も小さくする）、有限要素解は、真の解に近づくので、式（８）の第２項は０に、ＫＤＣＴは、真の解に収束する。
【００３５】
したがって、式（８）の第２項は、Ｋ値の誤差の目安となる指標に対応し、式（８）のｒ／Ｌの係数に（−Ｇ’／２）を乗じてＫ値の次元を持たせた次のΔＫＤＣＥを誤差指標とし、ＤＣＥ  Ｉｎｄｅｘと称することにする。
【００３６】
【数１８】
0003564538.tif 000017
【００３７】
この誤差指標ＤＣＥ  Ｉｎｄｅｘの特長は、特異要素の変位から直接評価が可能であること、要素を無限に小さくしていくと零に収束すること、さらにＫ値そのものの次元を持つという点にある。その結果、従来提案されている他の誤差指標にありがちな数回の解析により、その収束をみるという性格のものではないため、ただ一度の解析によりＫ値誤差を推定することができる。これは、解析に要する労力、結果的にコストを大幅に低減することが可能となる。
【００３８】
次に、この誤差指標を用いた応力拡大係数の誤差評価方法の一例について説明する。以下に示す例では、特に、上述した誤差指標ΔＫＤＣＥを組み込んだ応力拡大係数誤差評価プログラムを記憶した記憶媒体、及び、この応力拡大係数誤差評価プログラムを実行するシステムについて説明する。
【００３９】
図２に示すように、このシステム１は、制御手段として機能するＣＰＵ１０と、入力手段として機能する入力装置１２と、計算結果などを表示する表示装置１４と、記憶媒体として機能する記憶装置１６を備えて構成されている。この記憶装置１６としては、フロッピー（登録商標）ディスクやハードディスクなどの磁気ディスク、コンパクトディスクなどの光ディスク、及び、これらを駆動する駆動装置などを含んで構成されている。
【００４０】
記憶装置１６は、応力拡大係数誤差評価プログラムなどの他に種々の情報を記憶している。この応力拡大係数誤差評価プログラムは、特異要素を用いた有限要素解析を行うプログラムや、応力拡大係数ＫＤＣＴ及びこの応力拡大係数ＫＤＣＴに含まれる誤差ΔＫＤＣＥを算出するプログラムなどを含んでいる。
【００４１】
入力装置１２は、応力拡大係数ＫＤＣＴの誤差ΔＫＤＣＥを評価するために必要な所定条件の入力を受け付ける。ＣＰＵ１０は、入力装置１２を介して入力された所定条件に基づいて、記憶装置１６に記憶された応力拡大係数誤差評価プログラムを実行する。表示装置１４は、ＣＰＵ１０により応力拡大係数誤差評価プログラムにしたがって算出された応力拡大係数ＫＤＣＴ、及び算出された応力拡大係数ＫＤＣＴに含まれる誤差ΔＫＤＣＥといった計算結果などを表示する。
【００４２】
次に、上述したシステム１において、応力拡大係数誤差評価プログラムに基づいて、応力拡大係数ＫＤＣＴ及びこの応力拡大係数ＫＤＣＴに含まれる誤差ΔＫＤＣＥを算出する方法について図３に示したフローチャートを参照して説明する。
【００４３】
すなわち、ＣＰＵ１０は、記憶装置１６に記憶される応力拡大係数誤差評価プログラムに基づいて、必要な所定条件の入力を受け付ける（ＳＴ１）。このとき、入力される所定条件としては、き裂長さａ、特異要素数ｍ、特異要素寸法Ｌ、ヤング率Ｅやポアッソン比νなどの物性値、作用する応力σなどの境界条件、領域分割した接点の座標などである。
【００４４】
続いて、ＣＰＵ１０は、き裂周辺に配置した特異要素を用いて有限要素解析を行い、特異要素上の節点のｘ、ｙ方向各々の変位Ｕ（ｒ，θ）、Ｖ（ｒ，θ）を求める（ＳＴ２）。特に、ここでは、図１に示したように、き裂２０の周辺の４点について、有限要素解析を行って変位を求める。
【００４５】
続いて、ＣＰＵ１０は、ステップＳＴ２において求めた節点変位を用いて、上述した式（６）（実質的に式（６）’と同一）により、所要の応力拡大係数ＫＤＣＴを求める（ＳＴ３）。
【００４６】
続いて、ＣＰＵ１０は、ステップＳＴ２において求めた節点変位を用いて、上述した式（９）により、ステップＳＴ３において求めた所要の応力拡大係数ＫＤＣＴに含まれる誤差ΔＫＤＣＥを見積もる（ＳＴ４）。
【００４７】
続いて、ＣＰＵ１０は、最終的に正しいと思われる応力拡大係数Ｋ値をＫＤＣＴ−ΔＫＤＣＥとして推定する（ＳＴ５）。
【００４８】
次に、解析解が存在する二次元弾性問題に対し、特異要素を用いた有限要素解析を行い、式（６）により求めた応力拡大係数ＫＤＣＴと、解析解Ｋｒｅｆの差Ｋｅｒｒｏｒ＝（ＫＤＣＴ−Ｋｒｅｆ）を式（９）により求めた誤差ΔＫＤＣＥと、を比較した例について説明する。
【００４９】
いずれも、ヤング率Ｅ＝２０６ＧＰａ、ポアッソン比ν＝０．３とした。κは平面応力条件として評価した。また、特異要素数ｍは８から始め１６、２４、３０と増加させ、特異要素寸法比Ｌ／ａについてはＬ／ａを１／３、１／６、１／１２、１／２４と変化させた（具体的には特異要素のみを再分割）。
【００５０】
図４は、一様引張σ＝９．８ＭＰａを受ける幅Ｗ＝１０ｍｍ、長さＨ＝２Ｗ＝２０ｍｍの片側き裂はりに無次元き裂長さξ＝ａ／Ｗ＝０．１、０．３のき裂が存在する場合を示すものであり、図５は、図４の場合についてその各々に対してＫＩＤＣＴを求め、これと次式ＫＩｒｅｆより求まるＫＩｅｒｒｏｒをΔＫＩＤＣＥと比較した結果（ｍ＝１６）を示す。
【００５１】
【数１９】
0003564538.tif 000018
【００５２】
図中、特に説明はしていないが、各ａ／Ｗ毎にＬ／ａに対応する４個のデータがあり、Ｌ／ａとΔＫＩＤＣＥの大小関係が対応している。図中、ＫＩｅｒｒｏｒとΔＫＩＤＣＥの差は、各点と原点を通る傾き１の直線との垂直方向距離として表され、その最大値はａ／Ｗ＝０．１のＬ／ａ＝１／３のときの０．０１５６、すなわちＫＩｒｅｆを基準にしてその１．５６％である。
【００５３】
Ｌ／ａを小さくしていくと共にΔＫＩＤＣＥは減少し、それに従ってＫＩｅｒｒｏｒも減少していくことがわかる。図示の範囲内では、ΔＫＩＤＣＥは、ＫＩｅｒｒｏｒと同程度の値をとるとして良いようであり、上述した応力拡大係数誤差評価プログラムに基づいた計算結果の妥当性を裏付けている。
【００５４】
なお、上述した比較例のほかに、種々の条件の異なるケースにおいても同様の検証を行った結果、この発明に係る応力拡大係数の誤差評価方法は、極めて精度が高いことが確認された。
【００５５】
上述したように、この発明によれば、特異要素を用いる有限要素解析結果をもとに応力拡大係数（Ｋ値）を評価するための誤差評価指標を開発し、これをＤＣＥ  Ｉｎｄｅｘと名付けた。このＤＣＥ  Ｉｎｄｅｘは、き裂先端変位場の関数形が既知であり、特異要素についてはこの変位の関数形を一部表現しうること、またＫ値をき裂先端変位のみから評価しうることに着目し、Ｋ値の次元を有する誤差指標として開発したものである。
【００５６】
ＤＣＥ  Ｉｎｄｅｘは、Ｋ値に含まれる誤差そのものではないが、上述したように解析解が分かっている問題への適用結果によると、ＤＣＥ  Ｉｎｄｅｘは、特に多くの問題で重要となるモードＩのＫ値評価において、実際的なθ方向の分割を用いた有限要素解析により、Ｋ値誤差そのものに近い値をとるものとなり、実際的な方法としてそれによるＫ値補正に有用であることが確認された。
【００５７】
したがって、短時間で応力拡大係数ＫＤＣＴを求めるとともにその応力拡大係数に含まれる誤差ΔＫＤＣＥの程度を評価することが可能であり、精度の高い応力拡大係数を算出することができる。
【００５８】
【発明の効果】
以上説明したように、この発明によれば、ただ一度の有限要素解析を通じて所要の応力拡大係数ＫＤＣＴと共にその誤差をΔＫＤＣＥにより見積もることができ、それを踏まえた補正や得られたＫ値解の実用的観点からの適否を判断することができる。この結果、従来提案されている他の誤差指標にありがちな数回の解析により、その収束をみるという性格のものではないため、ただ一度の解析によりＫ値誤差を推定することができる。したがって、ただ一度の解析により求めた応力拡大係数ＫＤＣＴをΔＫＤＣＥにより補正することが可能であり、結果的に、精度の高い応力拡大係数（Ｋ値）を短時間で算出することが可能な応力拡大係数の誤差評価方法を提供することができる。
【図面の簡単な説明】
【図１】図１は、き裂まわりに８個の三角形特異要素を配置した例を示すもので、節点の位置を示すと同時に変数の定義を示す。
【図２】図２は、この発明の応力拡大係数の誤差評価方法が適用されるシステムの構成を概略的に示すブロック図である。
【図３】図３は、図２に示したシステムに適用される応力拡大係数の誤差評価方法を説明するためのフローチャートである。
【図４】図４は、この発明の適用例を示すものであり、対象とした構造及び荷重条件を示す図である。
【図５】図５は、実際のＫ値誤差ＫｅｒｒｏｒとＫ値誤差指標ΔＫＤＣＥを比較した図である。
【符号の説明】
１…システム
１０…ＣＰＵ
１２…入力装置
１４…表示装置
１６…記憶装置
２０…き裂
【図１】
0003564538.tif 000021
【図２】
0003564538.tif 000022
【図３】
0003564538.tif 000023
【図４】
0003564538.tif 000024
【図５】
0003564538.tif 000025