(57)【要約】【課題】上界近似値及び下界近似値による挟み打ちにより危険設計又は安全設計に要する精度の高い近似解を得ることのできる近似解の収束手段の提供【解決手段】片持平面板の荷重中点の変位の真の解は上界近似値７１及び下界近似値の間にあって、総自由度の増加に追従して収束する上界近似値７１及び下界近似値の間で収束する。ここで、上界近似値７１は所定の応力関数に基いて数値解析装置（上界近似値取得手段）によって算出し、取得したものであり、下界近似値７２は新エネルギー原理に基く式から導出される変位関数に基いて算出し、取得したものである。また、総自由度数は一要素当たりの自由度数×メッシュ数である。【選択図】　図７
2004046380.tif 000002
【特許請求の範囲】
【請求項１】
負荷の対象領域の総自由度を増加するにつれて単調減少に変位と応力の値をとりながら収束する下界近似値を得る下界値取得手段と、前記対象領域の総自由度を増加するにつれて単調増加に変位と応力の値をとりながら収束する上界近似値を得る上界値取得手段と、を備え、
前記下界値取得手段によって取得された近似値と前記上界値取得手段によって取得された近似値との差が所定の閾値より小となる下界近似値又は上界近似値を前記対象領域全体の変位及び応力値として取得する、
ことを特徴とする、近似解の収束手段。
【請求項２】
負荷の対象領域を小領域に分割し、分割された各小領域について前記負荷により生じた変位を計算し、前記計算により取得した前記各小領域の変位を隣接する各小領域の内部境界又は外部境界との相関関係により座標変換し、前記座標変換された前記各小領域の変位を基に前記対象領域全体の変位を算出し、応力を得る有限要素法において、
前記対象領域の総自由度を増加するにつれて単調減少に変位と応力の値をとりながら収束する下界近似値を得る下界値取得手段と、前記対象領域の総自由度を増加するにつれて単調増加に変位と応力の値をとりながら収束する上界近似値を得る上界値取得手段と、を備え、
前記下界値取得手段によって取得された近似値と前記上界値取得手段によって取得された近似値との差が所定の閾値より小となる下界近似値又は上界近似値を前記対象領域全体の変位及び応力値として取得する、
ことを特徴とする、有限要素法における近似解の収束手段。
【請求項３】
前記下界近似値は安全設計に用いられることを特徴とする請求項１又は２記載の近似解の収束手段。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、コンピュータ等の数値解析装置を用いて固体力学境界値及び固有値等の近似解を求める際の近似解の収束方法に関する。
【０００２】
【従来の技術】
工学応用分野ではコンピュータ等の数値解析装置を用いて固体力学境界値及び固有値等の近似解を求め、設計等に用いている。
【０００３】
例えば、従来の有限要素法（以下、節点有限要素法と記す）ではメッシュ分割を細かくした場合、つまり、要素分割数を増加するにつれて図９（ａ）に示すように単調増加に変位と応力の値をとりながら収束する近似解（上界近似解と記す）を用いて設計を行っている。このために数値計算時に、計算結果である近似解の検証、つまり、算出値と閾値との比較を行い、算出値が閾値以下となったときに算出値を解とするようにしている。
【０００４】
【発明が解決しようとする課題】
しかしながら、有限要素法では下界解の導出ができないため、上述したように危険側の設計、つまり、上界近似値を用いた設計を行うことはできるが、図９（ｂ）に示すように要素分割数を増加するにつれて単調減少に変位と応力の値をとりながら収束する近似解（下界近似解と記す）を用いた設計ができないので、構造物の剛性を一般に過大に評価することとなり、設計上非安全側（危険設計）の答えを常に与えることとなるといった問題点がある。また、上界近似値の収束性が保証されていないため、実験によって確かめるか、解析解との対比が必要となるといった問題点があった。
【０００５】
本発明は上記問題点を解決するために創案されたものであり、上界近似値及び下界近似値による挟み打ちにより危険設計又は安全設計に要する精度の高い近似解を得ることのできる近似解の収束手段の提供を目的とする。
【０００６】
【課題を解決するための手段】
上記課題を解決するために、第１の発明の近似解の収束手段は、負荷の対象領域の総自由度を増加するにつれて単調減少に変位と応力の値をとりながら収束する下界近似値を得る下界値取得手段と、対象領域の総自由度を増加するにつれて単調増加に変位と応力の値をとりながら収束する上界近似値を得る上界値取得手段と、を備え、下界値取得手段によって取得された近似値と上界値取得手段によって取得された近似値との差が所定の閾値より小となる下界近似値又は上界近似値を対象領域全体の変位及び応力値として取得する、ことを特徴とする。
【０００７】
また、第２の発明の近似解の収束手段は、負荷の対象領域を小領域に分割し、分割された各小領域について負荷により生じた変位を計算し、計算により取得した各小領域の変位を隣接する各小領域の内部境界又は外部境界との相関関係により座標変換し、座標変換された各小領域の変位を基に対象領域全体の変位を算出し、応力を得る有限要素法において、対象領域の総自由度を増加するにつれて単調減少に変位と応力の値をとりながら収束する下界近似値を得る下界値取得手段と、対象領域の総自由度を増加するにつれて単調増加に変位と応力の値をとりながら収束する上界近似値を得る上界値取得手段と、を備え、下界値取得手段によって取得された近似値と上界値取得手段によって取得された近似値との差が所定の閾値より小となる下界近似値又は上界近似値を対象領域全体の変位及び応力値として取得する、ことを特徴とする。
【０００８】
また、第３の発明は上記第１又は第２の発明の近似解の収束手段において、下界近似値は安全設計に用いられることを特徴とする。
【０００９】
【発明の実施の形態】
実施例の説明に先立ち、まず、無節点有限要素法の概要を説明する。無節点有限要素法は、本願出願人により同日出願の発明（名称「無節点有限要素法」）であり、要素間の境界辺上の変位及び応力を節点を用いずに求める方法である。そのため、無節点有限要素法では次式１のエネルギー原理式で表される無節点要素の概念を導入している。
【００１０】
【数１】
2004046380.tif 000003
ここに、ｕｉ：変位成分、ｔｉ：境界力成分、σｉｊ：応力成分、
バーｕｉ：規定変位　　バーｔｉ：表面荷重成分又は隣接要素から外力として作用する境界力
バーｐｉ：物体力成分
である。
【００１１】
より具体化するため、平面応力について図１の三角形要素を例として無節点有限要素法による定式化について説明する。
図１は無節点有限要素法における無節点要素（小領域）の一実施例を示す図であり、符号１０は無節点要素としての三角形要素、符号１１〜１３は隣接要素（無節点要素）の例である。また、符号ｓ１、ｓ２、ｓ３は無節点要素１０の辺、符号ｓ３、ｓ４、ｓ５は隣接要素１１の辺、符号ｓ６、ｓ７、ｓ１は隣接要素１２の辺、符号ｓ８、ｓ９、ｓ２は隣接要素１３の辺であり、無節点要素１０と隣接要素１１〜１３が共有する辺は内部境界、辺ｓ４〜ｓ９は外部境界である。この場合、図１の三角形要素１０の内部の点Ｐの変位ｕを次式（式２）のように有限べき級数に展開する。
【００１２】
【数２】
2004046380.tif 000004
ここで、ａ１、ａ２、・・・は未知数である。
【００１３】
式２を座標と未知数に分離し、マトリックス−ベクトル形式で表すと次式３のようになる。
【００１４】
【数３】
2004046380.tif 000005
上記式３からひずみεが与えられる。また、ひずみ−応力関係式（σ＝Ｄ（ε−ε０）＋σ０　；Ｄは材料定数を含んだ弾性マトリックス）から応力σが与えられる。また、三角形要素１０の点Ｐの応力σと境界ｓ１〜ｓ３の境界力ｔｉは次式（式４、式５）のように表すことができる。
【００１５】
【数４】
2004046380.tif 000006
上記式３、４、５を式１に代入して数学的操作を施すと次式（式６）を得ることができる。
【００１６】
【数５】
2004046380.tif 000007
図２は４要素に分割された平面板の例である。この例では要素１１の左辺で変位ｕ，ｖが固定され、要素２４の右辺にｘ方向の分布荷重バーｔが分布する。この例に上記式６を摘要すると次式（式７）が得られる。なお、要素２１〜２４は無節点要素である。
【００１７】
【数６】
2004046380.tif 000008
ここに、Ｍ１１、Ｍ１２はマトリックス、ａｉは要素１１の未定係数ベクトルである。
【００１８】
ここで、式７をマトリックス−ベクトル表示すると次式（式８）になる。
【００１９】
【数７】
2004046380.tif 000009
この式１２が計算すべき最終的な連立一次方程式である。環境条件は導入済みであるからコンピュータでマトリックス演算を行って連立一次方程式の係数値｛ａ１、ａ２、ａ３、ａ４　｝を算出し、要素２４の右辺に分布するｘ方向の分布荷重バーｔを得ることができる。
【００２０】
（実施例）
１．ハードウエア構成
図３は無節点有限要素法による解析計算を実行可能な数値解析システムの一実施例の構成を示すブロック図であり、この例は数値解析システムとしてパーソナルコンピュータ（以下、パソコン）を例としている。
図３で、パソコン５は制御部５１、マウス等のポインティングデバイス、入力装置５３、作業用メモリ５４、保存記憶用メモリ５５、表示装置５６、プリンタ５７を有している。
制御部５１はＣＰＵ（又はＭＰＵ）およびプログラム格納用ＲＯＭのほか周辺制御回路を備えた構成をなし、装置全体の動作を制御すると共に、保存記憶用メモリ２３に保存記憶されている無節点解析計算プログラム群（無節点解析手段）により本発明の無節点有限要素法に基く無節点解析処理の実行を制御する。
【００２１】
ポインティングデバイス５２はオペレータによる領域指定時の位置や描画パラメータの指定や各種操作指示等の手段であり、操作がなされると位置情報を制御部５１に送出する。また、入力装置５３はキーボードやマウス等のポインティングデバイス等を有しており、ユーザはデータやパラメータを入力したりすることができる。
【００２２】
また、作業用メモリ５４はＲＡＭ等の一時記憶メモリからなり、パソコンの起動時にプログラム格納用ＲＯＭに記憶されている制御プログラムが読み出されて記憶されると共に、固体力学の解析処理が選択されると具体的なモデル毎に対応の解析用プログラムや、上界近似値計算プログラムおよび下界近似値計算プログラムを含む近似解収束プログラムが保存記憶メモリ５５から適時読み出されて記憶される。
【００２３】
また、保存記憶用メモリ５５は磁気ディスクやＦＤ等のリムーバブルな記録媒体からなっており、本発明の無節点有限要素解析手順を実現するために作成されている無節点解析用プログラムモジュール群及び解析支援プログラムや、入力データ、パラメータおよび解析結果（出力データ）、精度評価用データのほか、描画プログラムを含む対象領域設定プログラムや、領域分割プログラム、中抜き演算によるマトリックス方程式計算プログラム、上界近似値計算プログラムおよび下界近似値計算プログラムを含む近似解収束プログラム等のサブプログラム群やサポートプログラム群を保存記録している。
【００２４】
２．無節点有限要素法による解析手順
図４は無節点有限要素法による解析手順を示すフローチャートである。
まず、解析対象とする材料や構造物等（以下、対象領域と記す）を特定すると共に解析対象とする部位（以下、対象領域と記す）と対象領域中の特異領域（例えば、平板に穿たれた穴の周辺、屈曲面の周辺、溶接個所及びその周辺等）を特定し、そのサイズや重量等、特性、荷重等の負荷量、値からの方向等の条件値（パラメータとして変化させる場合は、そのピッチも含む）を設定する。なお、対象領域は通常は２次元（平面）又は３次元（立体）であるがこれに限定されない（ステップＳ１）。
【００２５】
次に、対象領域（特異領域がある場合はそれも含む）を図５（ａ）の有孔板の例に示すようにモデル化する。モデル化にあたってはモデル全体を１とした正規化を行うことが望ましく、更に、図５（ｂ）の例に示すように解析計算が効率よくできるように解析対象領域（図５（ｂ）の例で符号７１’で示す部分）を決定することが計算上望ましい。モデルの作成及び正規化は後述するようにパソコン５を用いて容易に行うことができる（ステップＳ２）。
【００２６】
次に、モデルの解析対象領域が２次元の場合はｎ×ｎ（ｎは２以上の整数）のメッシュ、３次元の場合はｎ×ｎ×ｎの立方格子に分割し（ステップＳ３）、更に、上記ステップＳ１で特定した特異領域の周辺についてステップＳ３で分割した各メッシュを細分して小メッシュに分割し、各小メッシュ、つまり各要素に番号をつけると共に各辺の長さ（＝分割数でもよい）を対応付ける。この際、図６に示すように対象領域全体から特異領域を絞り込むように何段階かに分けてメッシュの細分を行うようにするとよい。このような領域分割は後述するようにパソコン５を用いて簡単な操作で行うことができる。（ステップＳ４）。
【００２７】
まず、第１段階（つまり、ステップＳ３）でメッシュ（２次元の場合は正方格子、３次元の場合は立方格子）に分割した各小領域（要素）毎にその座標値及び材料等の条件値（大きさ、重量、特性値、負荷、力の方向等（モデルが正規化されている場合は条件値も正規化したもの））を与え、前述した式８の連立一次方程式を演算して各要素毎の変位（未定係数ベクトル）を算出する。つまり、無節点有限要素法では要素毎に結果が算出される（ステップＳ５）。
【００２８】
次に、第１段階で得た小領域のうち段階ごとに細分された各小要素にその座標値及び材料等の物性値を与え、前述した式８の多元連立一次方程式（正方マトリックス方程式）を演算して各小要素毎の変位（未定係数ベクトル）を算出する。なお、無節点有限要素法では各要素を独立的に扱えるので、本願出願人が１９９９年５月２５日に出願の発明（マトリックスの中抜き圧縮処理（特開平２０００−３３９２９５））に基くマトリックスの中抜き圧縮処理を行ない、連立一次方程式の演算を行い、未定係数ベクトル（図５の例では各要素ごとの変位）の近似解（上界近似解及び下界近似解）として求めることができる（ステップＳ６）。
【００２９】
次に近似精度の検証を行い、精度が許容誤差内になった場合は解析計算を狩猟する。また、許容誤差に至らない場合はステップＳ４に戻って各要素を更に分割するか、ステップＳ１に戻って条件値の設定を変更する（ステップＳ７）。
【００３０】
なお、ステップＳ７の近似精度の検証は本発明の無節点有限要素法における近似解の収束手段、つまり、下記に述べる上界近似値及び下界近似値による挟み撃ち法によって行う。
【００３１】
（実施例）
以下、本発明の近似解の収束方法（つまり、上界近似値及び下界近似値による挟み撃ち法）について説明する。なお、以下の説明で、総自由度＝小領域数×相関関係式の項数、である。
【００３２】
式１で表されるエネルギー原理は、Ｌｓｇｒａｎｇｅ乗数（Ｌｓｇｒａｎｇｅは数学者の名前）を導入したために停留原理ではあるが最小原理ではない。従って真の解への収束性は必ずしも保証されない。これに対して、本発明に係わる新たなエネルギー原意は最小原理である。従って、この新原理を用いた場合収束性が保証される。
【００３３】
ここで、例えば、図４のステップＳ３での領域分割のように粗い要素分割の時に図９（ｂ）に示すように変位と応力が高めの値になり、要素分割数を増加、正確には総自由度（＝要素分割数と次式９を級数展開した式の項数との積）を増加するにつれて単調減少に変位と応力の値をとりながら収束する近似値を下界（Ｌｏｗｅｒ　ｂｏｕｎｄ）近似値をいう。
【００３４】
有限要素法では一般には真の解が分からないので数値解析による近似解を用いて設計を行うが、近似解が下界の場合は安全設計になる。一方、下界近似値とは反対に、総自由度を増加するにつれて単調増加に変位と応力の値をとりながら収束する近似値、つまり上界近似値を用いれば危険側の設計になる。無節点有限要素法に下記式９で表される新エネルギー原理を摘要すると解の下界を得ることができる。
【００３５】
【数８】
2004046380.tif 000010
上記式９の右辺の第１項に部分積分とＧａｕｓｓ（ガウス（数学者名））の発散定理を導入し、Πｉ　の最小条件から次式１０を得る。
【００３６】
【数９】
2004046380.tif 000011
上記式１０を用いてメッシュ分割を細かくするか或いは一要素あたりの自由度を増加すると最小ポテンシャル原理とコンプリメンタリエネルギー原理を同時に限りなく充足度を高めることができ、解の下界に収束する。また、左辺の第３項の平衡条件をアプリオリに満たす許容関数を用いると第３項が消えるので、下記式１１のＴｒｅｆｆｔｚ型辺分原理を得る。
【００３７】
【数１０】
2004046380.tif 000012
図６は片持平面板解析モデルの一実施例を示す図であり、片持平面板６は８個の正方メッシュ６１で分割されている。また、左辺は固定されており、右辺に荷重ｘが働く。また、ＵＡ　は荷重辺の中点であるＡ点の変位（計算値）、σＢ　はＢ点の応力（計算値）である。
【００３８】
図７は図６で示した片持平面板の荷重辺中点の変位の上界近似値及び下界近似値の分布を示す図であり、図７（ａ）は分布図を、図７（ｂ）は近似値の例を示す。図７から明らかなように真の解は上界近似値７１及び下界近似値７２の間にあって、総自由度の増加に追従して収束する上界近似値７１及び下界近似値７２の間で収束する。
【００３９】
ここで、図７（ａ）の上界近似値７１は下記式１２の応力関数（Ｇｏｕｒｓａｔ）法に基いてパソコン５等の数値解析装置（上界近似値取得手段）によって算出し、取得したものであり、下界近似値７２は新エネルギー原理を表す式９から導出される下記式１３の変位関数を基にパソコン５等の数値解析装置（下界近似値取得手段）によって算出し、取得したものである。また、総自由度数は一要素当たりの自由度数×メッシュ数であり、例えば、総自由度＝１２８＝１６×８（縦の分割数×横の分割数）である。
【００４０】
また、図８は図６で示した片持平面板の荷重辺中点の応力の上界近似値及び下界近似値の分布を示す図であり、図８（ａ）は分布図を、図８（ｂ）は近似値の例を示す。図８から明らかなように真の解は上界近似値８１及び下界近似値８２の間にあって、総自由度の増加に追従して収束する上界近似値８１及び下界近似値８２の間で収束する。
【００４１】
ここで、図８（ａ）の上界近似値８１は下記式１２の応力関数（Ｇｏｕｒｓａｔ）法に基いてパソコン５等の数値解析装置（上界近似値取得手段）によって変位を算出し、ひずみ応力関係式を基に変換して取得したものであり、下界近似値７２は新エネルギー原理を表す式９から導出される下記式１３の変位関数を基にパソコン５等の数値解析装置（下界近似値取得手段）によって変位を算出し、ひずみ応力関係式を基に変換して取得したものである。
【００４２】
【数１１】
2004046380.tif 000013
但し、１６自由度数／要素とする。
【００４３】
【数１２】
2004046380.tif 000014
但し、１６自由度数／要素とする。
【００４４】
上記説明では無節点有限要素法の場合を例として説明したが、本発明の適用範囲は無節点有限要素法による数値解析に限定されない。つまり、従来の節点有限要素法を含む有限要素法に摘要でき、更に、弾性モデルを含む構造工学の諸モデルの数値解析についても摘要可能である。
【００４５】
一般に節点法を用いて数値解析を行うと解の上界のみが得られる。これに対し、変位関数法を用いると下界が得られるので構造工学上必要な問題について上界近似値と下界近似値を計算して真の値を挟み打ちすることができる。また、設計時には上界近似値と下界近似値を得てその差が所定の閾値より小となる下界近似値又は上界近似値（或いはそのような下界近似値又は上界近似値の中間値）を前記対象領域全体の変位及び応力値として取得することができる。これにより、危険設計側及び安全設計側からのアプローチが可能となる。
【００４６】
以上、本発明の一実施例について説明したが本発明は上記実施例に限定されるものではなく、種々の変形実施が可能であることはいうまでもない。
【００４７】
【発明の効果】
上記説明したように、本発明によれば構造工学上必要な問題について上界近似値と下界近似値を計算して真の値を挟み打ちすることができるので、従来のように上界近似値による危険側設計だけでなく、下界近似解を用いた安全設計ができるので、構造物の剛性を過大に評価することがない。また、真の解は上界近似値と下界近似解の間で収束するので収束性が保証され、実験によって確かめたり、解析解との対比を行う必要がなくなる。
【図面の簡単な説明】
【図１】無節点有限要素法における無節点要素（小領域）の一実施例を示す図である。
【図２】４要素に分割された平面板の例である。
【図３】無節点有限要素法による解析計算を実行可能な数値解析システムの一実施例の構成を示すブロック図である。
【図４】無節点有限要素法による解析手順を示すフローチャートである。
【図５】有孔板の引張荷重モデルの一例である。
【図６】片持平面板解析モデルの一実施例を示す図である。
【図７】図６で示した片持平面板の荷重辺中点の上界近似値及び下界近似値の分布を示す図である。
【図８】図６で示した片持平面板の荷重辺中点の上界近似値及び下界近似値の分布を示す図である。
【図９】上界近似値及び下界近似値の説明図である。
【符号の説明】
５　パソコン（上界近似値取得手段、下界近似値取得手段）
６１　上界近似値
６２　下界近似値
【図１】
2004046380.tif 000015
【図２】
2004046380.tif 000016
【図３】
2004046380.tif 000017
【図４】
2004046380.tif 000018
【図５】
2004046380.tif 000019
【図６】
2004046380.tif 000020
【図７】
2004046380.tif 000021
【図８】
2004046380.tif 000022
【図９】
2004046380.tif 000023