TOP > 国内特許検索 > 解析装置、解析方法及び解析プログラム > 明細書

明細書 :解析装置、解析方法及び解析プログラム

発行国 日本国特許庁(JP)
公報種別 特許公報(B2)
特許番号 特許第5886309号 (P5886309)
登録日 平成28年2月19日(2016.2.19)
発行日 平成28年3月16日(2016.3.16)
発明の名称または考案の名称 解析装置、解析方法及び解析プログラム
国際特許分類 G06F  17/50        (2006.01)
G06F  17/13        (2006.01)
FI G06F 17/50 612H
G06F 17/13
請求項の数または発明の数 16
全頁数 323
出願番号 特願2013-541831 (P2013-541831)
出願日 平成24年10月31日(2012.10.31)
国際出願番号 PCT/JP2012/078262
国際公開番号 WO2013/065764
国際公開日 平成25年5月10日(2013.5.10)
優先権出願番号 2011239812
2012085968
優先日 平成23年10月31日(2011.10.31)
平成24年4月4日(2012.4.4)
優先権主張国 日本国(JP)
日本国(JP)
審査請求日 平成26年4月30日(2014.4.30)
特許権者または実用新案権者 【識別番号】504176911
【氏名又は名称】国立大学法人大阪大学
発明者または考案者 【氏名】林 茂弘
個別代理人の代理人 【識別番号】100104444、【弁理士】、【氏名又は名称】上羽 秀敏
【識別番号】100112715、【弁理士】、【氏名又は名称】松山 隆夫
【識別番号】100125704、【弁理士】、【氏名又は名称】坂根 剛
【識別番号】100120662、【弁理士】、【氏名又は名称】川上 桂子
【識別番号】100132506、【弁理士】、【氏名又は名称】山内 哲文
審査官 【審査官】松浦 功
参考文献・文献 松村清重外2名,2次元未定浸水長問題の自己随伴変分原理について,2004年度秋季造船三学会連合大会 関西造船協会講演概要集[online],公益社団法人日本船舶海洋工学会 ,2004年11月25日,第23号,pp. 95-98,[平成27年5月19日検索],CiNii,URL,http://ci.nii.ac.jp/naid/110003885168/
坂井藤一,有限要素法と差分法の等価性およびある離散化手法,土木学会論文報告集[online],社団法人土木学会,1973年12月,第220号,pp. 39-52,[平成27年5月19日検索],J-STAGE,URL,https://www.jstage.jst.go.jp/article/jscej1969/1973/220/1973_220_39/_article/-char/ja/
CHEN, Chun Hsiung et al.,The Variational Principle for Non-Self-AdjointElectromagnetic Problems,IEEE Transactions on Microwave Theory and Techniques[online],IEEE,1980年 8月,Vol. MTT-28, No. 8,Pages 878-886,[平成27年5月19日検索],IEEE Exploer,URL,http://ieeexplore.ieee.org/xpls/abs_all.jsp?arnumber=1130186&tag=1
調査した分野 G06F 17/50
G06F 17/11 -17/13
IEEE Xplore
CiNii
特許請求の範囲 【請求項1】
処理対象となる系の構造及び系の構成要素の性質を示すデータを読み込んで、読み込んだデータに基づいて、前記系を表す方程式であって、求めるべき物理量を表す変数を含むn個の初期方程式を決定する初期方程式決定部と、
境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定部と、
前記n個の初期方程式を、2n個の変数の方程式に変形し、前記変形した2n個の変数の方程式において、前記境界条件によって既知となる変数を含む既知部と、未知の変数を含む未知部とを決定し、当該未知部について方程式の解を計算する演算部を備え
前記初期方程式決定部は、前記系の構成要素の節点における物理量を示すn次元の変数ベクトルを2つと、n行の行列を含むn個の方程式を前記初期方程式として決定し、
前記境界条件決定部は、前記変数ベクトルにおいて、値が既知となる変数の自由度の数と、値が未知となる変数の自由度の数が異なる境界条件も決定可能とし、
前記演算部は、前記2つの変数ベクトルに基づいて2n次元のベクトルを生成し、当該2n次元のベクトルの変数に基づいて前記n行の行列を2n列の行列に変形した上で、当該2n次元のベクトルの変数のうち、前記境界条件によって既知となる変数を含む既知部分と、未知の変数を含む未知部分であって、必ずしも自由度が同じでない既知部分と未知部分を決定し、前記未知部分の変数を、前記既知部分の変数で表す形式になるよう前記2n列の行列及び前記2n次元のベクトルを変形して、変形した行列を用いて、前記未知部分の変数を計算する、情報処理装置。
【請求項2】
処理対象となる系の構造及び系の構成要素の性質を示すデータを読み込んで、読み込んだデータに基づいて、前記系を表す方程式であって、求めるべき物理量を表す変数を含むn個の初期方程式を決定する初期方程式決定部と、
境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定部と、
前記n個の初期方程式を、2n個の方程式を含む方程式に変形し、前記変形した2n個の方程式を含む方程式において、前記境界条件によって既知となる変数を含む既知部と、未知の変数を含む未知部とを決定し、当該未知部について方程式の解を計算する演算部を備え、
前記初期方程式決定部は、前記物理量を表す変数を含むn個の微分方程式を決定し、
前記演算部は、前記初期方程式決定部が決定した前記n個の微分方程式の微分作用素と、当該微分作用素から決まる随伴微分作用素を用いて、前記n個の主連立微分方程式とn個の双対連立微分方程式を示すデータを生成し、前記n個の主連立微分方程式と前記n個の双対連立微分方程式の解を計算することにより、少なくとも1つの前記物理量を出力する、情報処理装置。
【請求項3】
解析対象の系を表す初期方程式の原初微分作用素及び変数の境界条件を設定する設定部と、
前記境界条件から定まる随伴境界条件を得る随伴境界条件演算部と、
下記式において、解析対象の解uを、非同次境界条件を満たす項uBjと同次境界条件を満たす項uHjとの和(uBj+uHj)で表し、uHjに前記境界条件及び前記随伴境界条件から定まる主境界条件を満たす関数群を代入し、δuに前記境界条件及び前記随伴境界条件から定まる双対境界条件を満たす関数群を代入して得られる連立方程式を解くことで、前記解析対象の解uを計算する演算部を備える、情報処理装置。
下記式において、Lijは、前記原初微分作用素を表す。

【数1】
JP0005886309B2_000756t.gif

【請求項4】
請求項1~3いずれか1項に記載の情報処理装置で計算された前記解が複数存在する場合、下記式の汎関数Πにおいて、変分がゼロとなるように、同次解に対するモード係数を決定し、決定したモード係数を用いて解を求める、情報処理装置。
下記式において、Sは系の内部領域、Lijは前記系が満たすべき微分方程式の微分作用素、f、ujは、物理量を表す変数を表す。

【数753】
JP0005886309B2_000757t.gif

【請求項5】
前記決定したモード係数近傍の値の入力をユーザから受け付け、入力された値のモード係数を用いて解を計算し、前記解または前記解から得られる情報を出力する、請求項4に記載の情報処理装置。
【請求項6】
処理対象となる系の構造及び系の構成要素の性質を示すデータを読み込んで、読み込んだデータに基づいて、前記系を表す微分方程式であって、求めるべき物理量を表す変数を含む微分方程式を初期方程式として決定する初期方程式決定部と、
境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定部と、
前記初期方程式から定まる主連立微分式の主変数と同じ数の双対変数と、双対連立微分方程式を定義した場合に、前記主連立微分方程式の微分作用素を前記変数へ作用させたものと前記双対変数との積分和すなわち内積を部分積分して得られる境界項を計算し、前記境界条件の下で、当該境界項をゼロにする双対変数の条件を、随伴境界条件として決定する随伴境界条件決定部と、
前記随伴境界条件と前記境界条件とを比較した結果を出力する判定部とを備え、
前記主連立微分方程式の微分作用素を前記変数へ作用させたものと前記双対変数との内積は、前記変数と前記双対連立微分方程式の微分作用素を前記双対変数へ作用させたものとの内積に等しい関係にあり、
前記判定部は、前記随伴境界条件と前記境界条件とが一致しないと判断した場合、解析対象が非自己随伴である旨のメッセージを出力するか、又は、請求項1~3のいずれか1項に記載の情報処理装置の演算部に、前記初期方程式に基づく演算を実行させる、情報処理装置。
【請求項7】
前記境界条件と前記随伴境界条件とが一致するか否かは、下記式を用いて、節点力Fおよび節点変位Uそれぞれの既知部Fb、Ubの組と、双対節点力Fおよび節点変位Uそれぞれの既知部Fb、Ubの組が一致するか否かによって判定する、請求項6に記載の情報処理装置。
下記式において、Rは、前記境界項を、Gは剛性率を表す。

【数754】
JP0005886309B2_000758t.gif

【請求項8】
処理対象となる系の複数の要素における物理量を表す2つのn次元の変数と、n行の行列で表されるn個の初期方程式を決定する設定部と、
境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定部と、
前記2つのn次元の変数のうち、前記境界条件によって既知となる変数を含む既知部と、未知の変数を含む未知部を決定し、1つの変数の未知部の自由度が、他の変数の既知部の自由度と等しくない場合に、前記境界条件が非自己随伴であると判断する判定部を備え
前記判定部は、前記境界条件が非自己随伴であると判断した場合、解析対象が非自己随伴である旨のメッセージを出力するか、又は、請求項1~3のいずれか1項に記載の情報処理装置の演算部に、前記初期方程式に基づく演算を実行させる、情報処理装置。
【請求項9】
解析対象の原初微分作用素及び変数の境界条件を設定する設定部と、
前記境界条件から定まる随伴境界条件を得る随伴境界条件演算部と、
原初微分作用素から定まる主微分作用素及び双対微分作用素を得て、主連立微分方程式及び双対連立微分方程式、並びに、前記境界条件及び前記随伴境界条件を用いて、主固有関数及び双対固有関数を求めることにより、連立微分方程式の解を計算する演算部と、を備える請求項に記載の情報処理装置。
【請求項10】
前記境界条件と前記随伴境界条件が一致するか否かを判定する自己随伴判定部を備え、
前記演算部は、
前記境界条件と前記随伴境界条件が一致すると判定された場合、原初微分作用素から自己随伴問題の自己随伴固有関数を求めることにより、前記自己随伴問題の解を計算する自己随伴演算部と、
前記境界条件と前記随伴境界条件が一致しないと判定された場合、
前記主連立微分方程式及び前記双対連立微分方程式、並びに、前記境界条件及び前記随伴境界条件を用いて、主固有関数及び双対固有関数を求めることにより、連立微分方程式の解を計算する非自己随伴演算部とを含む、請求項9に記載の情報処理装置。
【請求項11】
下記式(23)で表される主問題の連立微分方程式の解uは、非同次境界条件を満たす項uBjと同次境界条件を満たす項uHjとの和(uBj+uHj)で表し、
下記式(65)で表される双対問題の連立微分方程式の解uは、非同次随伴境界条件を満たす項uBjと同次随伴境界条件を満たす項uHjとの和(uBj+uHj)で表し、
Hjを、下記式(69)に示すように主固有関数φの和で表し、uHjを、下記式(70)に示すように双対固有関数φの和で表した場合に、
前記演算部は、下記式(105)で表される前記主連立微分方程式、及び下記式(106)で表される前記双対連立微分方程式を満たす主固有関数φ及び双対固有関数φを求める、請求項9又は10に記載の情報処理装置。
下記式において、c、cは、k次の係数、λは固有値となる定数、wは、重みを示す定数である。

【数4】
JP0005886309B2_000759t.gif

【請求項12】
前記境界条件決定部又は前記設定部は、ユーザから前記物理量を表す変数で値が未知の部分を示す情報の入力を受け付け、当該情報を用いて前記境界条件を決定する、請求項1~3、6~11のいずれか1項に記載の情報処理装置。
【請求項13】
コンピュータが、処理対象となる系の構造及び系の構成要素の性質を示すデータを読み込んで、読み込んだデータに基づいて、前記系を表す方程式であって、求めるべき物理量を表す変数を含むn個の初期方程式を決定する初期方程式決定工程と、
コンピュータが、境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定工程と、
コンピュータが、前記n個の方程式を、2n個の変数の方程式に変形し、前記変形した2n個の変数の方程式において、前記境界条件によって既知となる変数を含む既知部と、未知の変数を含む未知部とを決定し、当該未知部について方程式の解を計算する演算工程とを有し、
前記初期方程式決定工程では、前記系の構成要素の節点における物理量を示すn次元の変数ベクトルを2つと、n行の行列を含むn個の方程式を前記初期方程式として決定し、
前記境界条件決定工程では、前記変数ベクトルにおいて、値が既知となる変数の自由度の数と、値が未知となる変数の自由度の数が異なる境界条件も決定可能とし、
前記演算工程では、前記2つの変数ベクトルに基づいて2n次元のベクトルを生成し、当該2n次元のベクトルの変数に基づいて前記n行の行列を2n列の行列に変形した上で、当該2n次元のベクトルの変数のうち、前記境界条件によって既知となる変数を含む既知部分と、未知の変数を含む未知部分であって、必ずしも自由度が同じでない既知部分と未知部分を決定し、前記未知部分の変数を、前記既知部分の変数で表す形式になるよう前記2n列の行列及び前記2n次元のベクトルを変形して、変形した行列を用いて、前記未知部分の変数を計算する、情報処理方法。
【請求項14】
処理対象となる系の構造及び系の構成要素の性質を示すデータを読み込んで、読み込んだデータに基づいて、前記系を表す方程式であって、求めるべき物理量を表す変数を含むn個の初期方程式を決定する初期方程式決定処理と、
境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定処理と、
前記n個の方程式を、2n個の変数の方程式を含む方程式に変形し、前記変形した2n個の変数の方程式を含む方程式において、前記境界条件によって既知となる変数を含む既知部と、未知の変数を含む未知部とを決定し、当該未知部について方程式の解を計算する演算処理とをコ ンピュータに実行させ
前記初期方程式決定処理は、前記系の構成要素の節点における物理量を示すn次元の変数ベクトルを2つと、n行の行列を含むn個の方程式を前記初期方程式として決定する処理を含み、
前記境界条件決定処理では、前記変数ベクトルにおいて、値が既知となる変数の自由度の数と、値が未知となる変数の自由度の数が異なる境界条件も決定可能とし、
前記演算処理は、前記2つの変数ベクトルに基づいて2n次元のベクトルを生成し、当該2n次元のベクトルの変数に基づいて前記n行の行列を2n列の行列に変形した上で、当該2n次元のベクトルの変数のうち、前記境界条件によって既知となる変数を含む既知部分と、未知の変数を含む未知部分であって、必ずしも自由度が同じでない既知部分と未知部分を決定し、前記未知部分の変数を、前記既知部分の変数で表す形式になるよう前記2n列の行列及び前記2n次元のベクトルを変形して、変形した行列を用いて、前記未知部分の変数を計算する処理を含む、情報処理プログラム。
【請求項15】
コンピュータが、処理対象となる系の構造及び系の構成要素の性質を示すデータを読み込んで、読み込んだデータに基づいて、前記系を表す微分方程式であって、求めるべき物理量を表す変数を含むn個の微分方程式を初期方程式として決定する初期方程式決定工程と、
コンピュータが、境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定工程と、
コンピュータが、前記初期方程式から定まる主連立微分式の主変数と同じ数の双対変数と、双対連立微分方程式を定義した場合に、前記主連立微分方程式の微分作用素を前記変数へ作用させたものと前記双対変数との積分和すなわち内積を部分積分して得られる境界項を計算し、前記境界条件の下で、当該境界項をゼロにする双対変数の条件を、随伴境界条件として決定する随伴境界条件決定工程と、
コンピュータが、前記随伴境界条件と前記境界条件とを比較した結果を出力する出力工程を有し、
前記主連立微分方程式の微分作用素を前記変数へ作用させたものと前記双対変数との内積は、前記変数と前記双対連立微分方程式の微分作用素を前記双対変数へ作用させたものとの内積に等しい関係にあり、
前記出力工程では、コンピュータが、前記随伴境界条件と前記境界条件とが一致しないと判断した場合、解析対象が非自己随伴である旨のメッセージを出力するか、又は、請求項1~3のいずれか1項に記載の情報処理装置の演算部に、前記初期方程式に基づく演算を実行させる、情報処理方法。
【請求項16】
処理対象となる系の構造及び系の構成要素の性質を示すデータを読み込んで、読み込んだデータに基づいて、前記系を表す微分方程式であって、求めるべき物理量を表す変数を含むn個の微分方程式を初期方程式として決定する初期方程式決定処理と、
境界条件を示すデータとして、前記物理量を表す値を読み込んで、境界条件を決定する境界条件決定処理と、
前記初期方程式から定まる主連立微分式の主変数と同じ数の双対変数と、双対連立微分方程式を定義した場合に、前記主連立微分方程式の微分作用素を前記変数へ作用させたものと前記双対変数との積分和すなわち内積を部分積分して得られる境界項を計算し、前記境界条件の下で、当該境界項をゼロにする双対変数の条件を、随伴境界条件として決定する随伴境界条件決定処理と、
前記随伴境界条件と前記境界条件とを比較した結果を出力する出力処理をコンピュータに実行させ、
前記主連立微分方程式の微分作用素を前記変数へ作用させたものと前記双対変数との内積は、前記変数と前記双対連立微分方程式の微分作用素を前記双対変数へ作用させたものとの内積に等しい関係にあり、
前記出力処理は、前記随伴境界条件と前記境界条件とが一致しないと判断した場合、解析対象が非自己随伴である旨のメッセージを出力するか、又は、請求項1~3のいずれか1項に記載の情報処理装置の演算部に、前記初期方程式に基づく演算を実行させる、情報処理プログラム。
発明の詳細な説明 【技術分野】
【0001】
本発明は、連続体等の解析対象の運動を微分方程式で表し、境界条件に応じた解を計算することにより、連続体の運動又は状態を解析する技術に関する。
【背景技術】
【0002】
連続体の運動は連立偏微分方程式で記述されるのが通例で,各種境界条件に応じた解を得る事が目的となる。工学的には自己随伴問題を取り扱う事が多く,それに適した様々な解法が研究されており,なかでもヒルベルトの展開定理に基づく固有関数法が有用で,応用例も豊富にある。これに比べて,非自己随伴問題についての研究例は極めて少ないが,シュミットの展開定理(例えば、非特許文献1、2参照)に基づく固有関数法が存在する。さらに,ミフリン(例えば、非特許文献3、4参照)は,この固有関数を用いれば最小2乗法が正解を与える事を示している。しかし,これらの定理は積分方程式論の範疇にあって境界条件が積分核に埋め込まれる為,汎用法としては利用しづらい面がある。また,複数の未知関数を解く連立方程式の形ではない事も一因となり,必然的に,応用例も限られている。
【先行技術文献】
【0003】

【非特許文献1】R.クーラン,D.ヒルベルト;齋藤利弥,丸山滋弥 共訳:数理物理学の方法1,東京図書,初版1959.
【非特許文献2】Erhard Schmidt: Zur Theorie der linearen und nichtlinearen Integralgleichungen. I., Math. Ann. Bd. 63, 1907, pp.433-476.
【非特許文献3】B.A.フィンレイソン;鷲津久一郎,山本善之,川井忠彦 共訳:重みつき残差法と変分原理,培風館,1974.
【非特許文献4】S.G.Mikhlin; translated by T.Boddington: Variational Methods in Mathematical Physics, Pergamon Press,1964.
【発明の概要】
【発明が解決しようとする課題】
【0004】
弾性学,振動学,材料力学,構造力学などの分野では,理論体系の根幹に仮想仕事の原理が据えられる為,内容が自己随伴問題に偏っているのが現状であり,非自己随伴問題の存在すら認識されていない感がある。しかし,実際には非自己随伴問題は数多く存在するので,仮想仕事の原理の持つ曖昧性を認識して,より高い視点から力学の体系を再構築し,その可能性を幅広く奥深く伸展させ,更なる発展を誘う事が必要である。
【0005】
そこで、本発明は、連立微分方程式で表された非自己随伴問題の解を算出することできる、解析装置、解析方法、及び解析プログラムを提供することを目的とする。
【課題を解決するための手段】
【0006】
本願開示の解析装置は、解析対象に対する原初微分作用素Lijを表すデータを用いて,前記解析対象の微分方程式の解を求める解析装置であって,解析対象へ働く力をf、双対変位uの変分のことを双対変分δuとするとき,下記式によって,解uを計算する。
【0007】
【数1】
JP0005886309B2_000002t.gif

【図面の簡単な説明】
【0008】
【図1】固有関数と微分作用素や境界条件との関係を示す図
【図2】第1の実施形態の解析装置の構成の一例を示す機能ブロック図
【図3】第2の実施形態における解析装置の構成の一例を示す機能ブロック図
【図4】図3に示す解析装置の動作例を示すフローチャート
【図5】第3の実施形態における解析装置の構成の一例を示す機能ブロック図
【図6】図5に示す解析装置10bの動作例を示すフローチャート
【図7】第4の実施形態における解析装置の動作例を示すフローチャート
【図8】16節点四角形の要素の例を示す図
【図9】一様重力下の円環で「変位ゼロかつ表面力ゼロ」を境界条件とする解析結果を示す図
【図10】一様重力下の正方形で、上辺の境界条件を「変位ゼロかつ表面力ゼロ」とした場合の解析結果を示す図
【図11】周辺を固着した正方形の平板の面内変形の解析結果を示す図
【図12】正方形の面外変形平板要素の周辺を固着した場合の自重による変形を解析した結果を示す図
【図13A】左端自由固着の梁のモデルを図
【図13B】図13Aに示す梁の主固有関数φ、φ、φ、及び双対固有関数φ、φ、φの計算結果を示す図
【図13C】固有関数法により求めた解および解析解を示す図
【図14】無限領域に正方形の仮想境界を定義し,その中心に湧き出しを配置した場合の解析結果を示す図
【図15】無限領域に正方形の仮想境界を定義し,その中心に渦度を配置した場合の解析結果を示す図
【図16】ばね減衰系の例を示す図
【図17】固有関数の例を示す図
【図18】(525)式の解析解と,(494)式の固有関数法で得た結果の比較を示す図
【図19】(528)式の解析解と,(494)式の固有関数法で得た結果の比較を示す図
【図20】(531)式の解析解と,(494)式の固有関数法で得た結果の比較を示す図
【図21】質点ばね系の例を示す図
【図22】(610)式でのそれぞれの固有関数を示す図
【図23】(611)式でのそれぞれの固有関数を示す図
【図24】(625)式の解析解と,(578)式の固有関数法で得た結果の比較を示す図
【図25】(632)式の解析解と,(584)式の固有関数法で得た結果の比較を示す図
【図26】弦の静たわみの例を示す図
【図27】固有関数の例を示す図
【図28】(706)式の解析解と、(702)式の固有関数法で得た結果の比較を示す図
【図29】(709)式の解析解と、(702)式の固有関数法で得た結果の比較を示す図
【図30】(713)式の解析解と、(704)式の固有関数法で得た結果の比較を示す図
【図31】(715)式の解析解と、(704)式の固有関数法で得た結果の比較を示す図
【図32】梁の静たわみの例を示す図
【図33】固有関数の例を示す図
【図34】固有関数の例を示す図
【図35】(784)式の解析解と,(780)式の固有関数法で得た結果の比較を示す図
【図36】(787)式の解析解と、(780)式の固有関数法で得た結果の比較を示す図
【図37】(790)式の解析解と、(782)式の固有関数法で得た結果の比較を示す図
【図38】(793)式の解析解と、(782)式の固有関数法で得た結果の比較を示す図
【図39】複数の解のうちの一例を示す解以外の解の一つを示す図
【図40A】有限要素を用いて解析した場合に計算される複数の解の例を示す図
【図40B】有限要素を用いて解析した場合に計算される複数の解の例を示す図
【図41】4節点形状関数のための26節点を示す図
【図42】円環の変形と応力の解析解を示す図
【図43】モードSAの変形を示す図
【図44】モードASの変形を示す図
【図45】モードSSの変形を示す図
【図46】モードAAの変形を示す図
【図47】円環の第1モードを示す図
【図48】円環の第2モードを示す図
【図49】円環の第3モードを示す図
【図50】円環の第4モードを示す図
【図51】2モードを用いた場合の変形と応力分布を示す図
【図52】3モードを用いた場合の変形と応力分布を示す図
【図53】10モードを用いた場合の変形と応力分布を示す図
【図54】30モードを用いた場合の変形と応力分布を示す図
【図55】湧き出しの無次元速度の様子を示す図
【図56】渦度の無次元速度の様子を示す図
【図57】関数(1)の第1モード(m=1,N=1)を示す図
【図58】関数(2)の第1モード(m=1,N=1)を示す図
【図59】関数(3)の第1モード(m=1,N=1)を示す図
【図60】関数(4)の第1モード(m=1,N=1)を示す図
【図61】関数(5)の第1モード(m=1,N=1)を示す図
【図62】関数(6)の第1モード(m=1,N=1)を示す図
【図63】関数(7)の第1モード(m=1,N=1)を示す図
【図64】関数(8)の第1モード(m=1,N=1)を示す図
【図65】固有関数法を用いた湧き出しの速度分布を示す図
【図66】固有関数法と解析解の湧き出しの速度分布の差を示す図
【図67】情報処理装置の構成例を示す図
【図68】情報処理装置の構成例を示す図
【図69】境界条件設定処理の例を示すフローチャート
【図70】境界条件設定処理の例を示すフローチャート
【図71】境界条件設定処理の例を示すフローチャート
【図72】解析処理の例を示すフローチャート
【発明を実施するための形態】
【0009】
1.緒言
ここで,連立偏微分方程式で表された非自己随伴問題を解く為の,新しい固有関数法を示す。次に,本手法と最小2乗法,変分直接法,弾性学のエネルギ原理との関係について論じる。その後,仮想仕事の原理及び有限要素定式化の深部について論じる。さらに、2次元弾性体の有限要素の定式化について述べる。

2. 連立偏微分方程式
2.1 用語と造語
固有関数と微分作用素や境界条件との関係について,必要な用語や造語を図1に示す。

【0010】
解きたい問題を主問題[Primal Problem]と呼ぶ。これに与えられた非同次境界条件[Inhomogeneous Boundary Condition]から同次境界条件[Homogeneous Boundary Condition]が定まる。同じ同次境界条件を持つ問題は多数あり,境界条件の和を算法とする加群を成す。これを表の群[Heads Group]と呼ぶ。即ち,主問題の非同次境界条件が,要素のひとつとして,表の群に含まれる。主問題の解関数を主解[Primal Solution]と呼び,この変分を主変分[Primal Variation]と呼ぶ。

【0011】
内積[Inner Product]の部分積分により,原初微分作用素[Original Differential Operator]からは随伴微分作用素[Adjoint Differential Operator]が定まり,境界項をゼロとする条件からは同次随伴境界条件[Homogeneous Adjoint Boundary Condition]が定まる。同じ同次随伴境界条件を持つ問題は多数あり,同様に加群を成す。これを裏の群[Tails Group]と呼ぶ。その要素のひとつを非同次随伴境界条件[Inhomogeneous Adjoint Boundary Condition]と呼び,これを与える問題を双対問題[Dual Problem]と呼ぶ。双対問題の解関数を双対解[Dual Solution]と呼び,この変分を双対変分[Dual Variation]と呼ぶ。内積を通じて,両群は表裏一体に結ばれる。

【0012】
2.2 自己随伴作用素の例 (弾性体の静的釣合)
2次元弾性体の静的釣合方程式は,直交座標系で,
【数5】
JP0005886309B2_000003t.gif
となる。ここで,u,uはx,y方向の変位で,求解により主解となる。b,bはx,y方向の単位体積あたりの物体力,Gは剛性率である。ポアソン比[Poisson’s ratio]をνとして,定数μは,
【数6】
JP0005886309B2_000004t.gif
である。
平面応力,平面ひずみ状態共に,変位ひずみ関係は
【数7】
JP0005886309B2_000005t.gif
であり,応力ひずみ関係は,
【数8】
JP0005886309B2_000006t.gif

【0013】
となる。よって,応力成分は,
【数9】
JP0005886309B2_000007t.gif
となり,境界面の外向き単位法線ベクトルの成分をn,nとして,表面力p,pは,コーシーの公式[Cauchy’s formula]から,
【数10】
JP0005886309B2_000008t.gif

【0014】
である。(1)~(14)式のように定数μを用いた表式により,平面応力状態と平面ひずみ状態を統一的に取り扱えるようになる。(1),(2)式の微分作用素を原初微分作用素Lijと呼び,
【数11】
JP0005886309B2_000009t.gif
と定義すると,次章の部分積分により随伴微分作用素Lijは,
【数12】
JP0005886309B2_000010t.gif
となる場合を自己随伴微分作用素[Self-Adjoint Differential Operator]と呼ぶ。随伴微分作用素Ljiの全ての成分が,その添字j,iを反転させてi,jとした,原初微分作用素Lijに一致すれば,Lijは自己随伴微分作用素である。

【0015】
2.3 非自己随伴作用素の例 (ポテンシャル流れ)
2次元ポテンシャル流れの問題を,速度ポテンシャルを用いずに,速度だけで解こうとすれば,直交座標系で,
【数13】
JP0005886309B2_000011t.gif
の連立偏微分方程式となる。ここで,u,uはx,y方向の速度で,求解により主解となる。(18),(19)式の微分作用素を原初微分作用素Lijと呼び,
【数14】
JP0005886309B2_000012t.gif
と定義すると,次章の部分積分により随伴微分作用素Lijは,
【数15】
JP0005886309B2_000013t.gif
となる場合を非自己随伴微分作用素[Non-Self-Adjoint Differential Operator]と呼ぶ。随伴微分作用素Ljiが,その添字j,iを反転させてi、jとした,原初微分作用素Lijに一致しない成分がひとつでもあれば,Lijは非自己随伴微分作用素である。

【0016】
3.主問題と双対問題
3.1 境界条件の同次化
以降は,前章の弾性体の例に題材を絞って記す。主解となる変位u,uを主変位と呼びu,uと書く。方程式の右辺をf,fと書けば,連立偏微分方程式(1),(2)は,
【数16】
JP0005886309B2_000014t.gif
となる。非同次境界条件を満たす項に添字を付けてuBjとおき,同次境界条件を満たす項に添字Hを付けてuHjとおく。主変位ujをその和で,
【数17】
JP0005886309B2_000015t.gif
と表す。(10)~(14)式および上式に従って,応力や表面力もuBjあるいはuHjから成る項に分かれて,
【数18】
JP0005886309B2_000016t.gif
となる。uHj,pHjは同次境界条件を満たす未知関数であるが,pHjはuHjで表現されるから,最終的な未知関数はuHjのみである。同様に,uBj,pBjは非同次境界条件を満たす既知関数であるが,pBjはuBjで表現されるから,最終的な既知関数はuBjのみである。uBjが(23)式を満たすとは限らず,(24)式の形で(23)式を満たす。(24)式を(23)式に代入して,同次境界条件で表された連立偏微分方程式
【数19】
JP0005886309B2_000017t.gif
である。以上の手順は,微分方程式を解く為の常套手段で,既知関数uBjに応じた未知関数uHjを求める問題にすり替わる。

【0017】
3.2 双対変位
双対解となる変位を双対変位uと呼ぶ。前節同様に,非同次随伴境界条件を満たす項uBjと同次随伴境界条件を満たす項uHjの和で,双対変位uを,
【数20】
JP0005886309B2_000018t.gif
と表す。(10)~(14)式および上式に従って,応力や表面力もuBjあるいはuHjから成る項に分かれて,
【数21】
JP0005886309B2_000019t.gif
となる。uHj,pHjは同次随伴境界条件を満たす未知関数であるが,pHj
Hjで表現されるから,最終的な未知関数はuHjのみである。同様に,uBj,pBjは非同次随伴境界条件を満たす既知関数であるが,pBjはuBjで表現されるから,最終的な既知関数はuBjのみである。

【0018】
3.3 部分積分
境界表面をc,内部領域をsとする。方程式(23)に双対変位uを乗じた積分和,即ち,内積は,
【数22】
JP0005886309B2_000020t.gif
となる。左辺を部分積分して,
【数23】
JP0005886309B2_000021t.gif
となる。ここで,Rは境界項で,前章の弾性体の場合には,
【数24】
JP0005886309B2_000022t.gif
である。同次化された連立偏微分方程式(40)と関数uHjとの内積は,
【数25】
JP0005886309B2_000023t.gif
となる。左辺を部分積分して,
【数26】
JP0005886309B2_000024t.gif
となる。ここで,Rは境界項で,
【数27】
JP0005886309B2_000025t.gif
である。非同次境界条件と非同次随伴境界条件が共にゼロの場合には,(60)式と(63)式は一致する。同次境界条件がuHi,pHiに課された下で,境界項Rをゼロとする条件,即ち,同次随伴境界条件がuHi,pHiに課される。その結果,(62)式は,
【数28】
JP0005886309B2_000026t.gif
となる。
同次境界条件と同次随伴境界条件が一致する場合を自己随伴境界条件[Self-Adjoint Boundary Condition]と呼び,一致しない場合を非自己随伴境界条件[Non-Self-Adjoint Boundary Condition]と呼ぶ。
微分作用素と境界条件が共に自己随伴な問題を自己随伴問題[Self-Adjoint Problem]と呼び,どちらか一方でも非自己随伴な問題を非自己随伴問題[Non-Self-Adjoint Problem]と呼ぶ。

【0019】
3.4 主問題と双対問題の関係
(59),(64)式右辺の作用素に着目して,双対問題の連立偏微分方程式
【数29】
JP0005886309B2_000027t.gif
を定義する。上式右辺のfは,双対問題の外力である。上式に(42)式を代入して,双対問題の,同次化された連立偏微分方程式
【数30】
JP0005886309B2_000028t.gif

【0020】
である。(64)式に(40),(66)式を代入して,関係式
【数31】
JP0005886309B2_000029t.gif
を得る。これは,(主,双対)問題の解関数uHi,uHiと,外力項fHi,fHiとの間に成り立つ関係を示している 。

【0021】
4.連立固有値問題
4.1関数群
同次境界条件を満足する関数群をφとおき,この群のk番目の関数をφjkと書いて,解関数uHjをその和で,
【数32】
JP0005886309B2_000030t.gif
と表す。同様に,同次随伴境界条件を満足する関数群をφとおき,k番目の関数をφjkと書いて,解関数uHjをその和で,
【数33】
JP0005886309B2_000031t.gif
と表す。ここで,c,cはk次の係数である。ルベーグ[Lebesgue]の意味で「ほとんど至る所等しい」解関数を得る為には,関数群φ,φを固有関数とすべきである。また,(68)式に着目すれば,外力項fHjを関数群φから構成し,外力項fHjを関数群φから構成する事により,内積を通じて,係数c,c間に関係が生じる事も判る。

【0022】
4.2 関数群による応力と表面力
応力を表す(27),(30),(33)式に(69)式を代入すれば,
【数34】
JP0005886309B2_000032t.gif
となる。よって,関数群φにより生じる応力に添字Eを付けてσEx,σEy,τExyと表すと,
【数35】
JP0005886309B2_000033t.gif
となる。関数群φにより生じる表面力pEx,pEyは,
【数36】
JP0005886309B2_000034t.gif
となる。これは,関数群φによる表面力pEiに,境界条件を満足させる為の準備である。

【0023】
同様に,応力を表す(45),(48),(51)式に(70)式を代入すれば,
【数37】
JP0005886309B2_000035t.gif
となる。よって,関数群φにより生じる応力をσEx,σEy,τExy,と表すと,
【数38】
JP0005886309B2_000036t.gif
となる。関数群φにより生じる表面力pEx,pEy
【数39】
JP0005886309B2_000037t.gif
となる。これは,関数群φによる表面力pEiに,境界条件を満足させる為の準備である。

【0024】
4.3 (主,双対)連立微分方程式と(主,双対)連立固有値問題
4.1節の知見により,同次化された微分方程式(40)に注目して,主連立微分方程式[Primal Simultaneous Differential Equations]
【数40】
JP0005886309B2_000038t.gif
を定義する。さらに,微分方程式(66)に注目して,双対連立微分方程式[Dual Simultaneous Differential Equations]
【数41】
JP0005886309B2_000039t.gif
を定義する。ここで,w,wは定数で重み,λ,λは定数でいずれ固有値となる。

【0025】
主連立微分方程式(87)とφとの内積は,
【数42】
JP0005886309B2_000040t.gif
となる。左辺を部分積分して,
【数43】
JP0005886309B2_000041t.gif
となる。ここで,REは境界項で,
【数44】
JP0005886309B2_000042t.gif
である。4.1節により,同次境界条件がφ,pEiに課され,同次随伴境界条件がφ
,pEiに課される為,境界項Rはゼロとなる。よって,(88)式により,(90)式は,
【数45】
JP0005886309B2_000043t.gif
となる。(89)式と(92)式の左辺は同じであるから,その右辺に注目して,
λ=λ (93)
となるようにλ,λを定めれば,
【数46】
JP0005886309B2_000044t.gif
を得る。即ち,(93)式の下で関数φ,φの重み付き内積は等しい。

【0026】
一方,主連立微分方程式(87)と双対連立微分方程式(88)を相互に代入すれば,2つの連立固有値問題
【数47】
JP0005886309B2_000045t.gif
が得られる。これらは,共に,固有値がλとなる連立固有値問題である。(98)式を主連立固有値問題[Primal Simultaneous Eigenvalue Problem],(99)式を双対連立固有値問題[Dual Simultaneous Eigenvalue Problem]と呼ぶ。関数φを主固有関数[Primal Eigenfunction],φを双対固有関数[Dual Eigenfunction]と呼ぶ。左辺の作用素をそれぞれ,主微分作用素[Primal Differential Operator],双対微分作用素[Dual Differential Operator]と呼ぶ。(主,双対)微分作用素は共に自己随伴な微分作用素である。

【0027】
さて,関数φには同次境界条件が,関数φには同次随伴境界条件が課されている。主連立微分方程式(87)の右辺には関数φがある為,関数φは同次境界条件を満たすと共に,その微分係数の和Σijφが同次随伴境界条件を満たす必要がある。これらを合わせて主境界条件[Primal Boundary Condition]と呼ぶ。双対連立微分方程式(88)の右辺には関数φがある為,関数φは同時随伴境界条件を満たすと共に,その微分係数の和Σjiφが同次境界条件を満たす必要がある。これらを合わせて双対境界条件[Dual Boundary Condition]と呼ぶ。(主,双対)境界条件は共に自己随伴な境界条件である。

【0028】
これらの境界条件を課された(98)式のφおよび(99)式のφはそれぞれ直交性を有し,ヒルベルト空間[Hilbert space]の基底関数を構成する。例えば,m次の固有値と主固有関数をλ、φimとおき,n次の固有値と主固有関数をλ,φinとおくと,内積の演算から,
【数48】
JP0005886309B2_000046t.gif
となるから,主固有関数φに直交性がある事を知る。同様にして,m次の固有値と双対固有関数をλ,φimとおき,n次の固有値と双対固有関数をλ,φinとおくと,内積の演算から,
【数49】
JP0005886309B2_000047t.gif
となるから,双対固有関数φにも直交性がある事を知る。

【0029】
よって,(94),(97)式を考慮すれば,クロネッカーのデルタ[Kronecker’s delta]をδmnとして,
【数50】
JP0005886309B2_000048t.gif
となるよう,関数φ,φを正規化しておける。

【0030】
4.4 固有関数法への準備
主連立微分方程式(87)は,(93),(97)式より,
【数51】
JP0005886309B2_000049t.gif
となる。双対連立微分方程式(88)は,変形なく,
【数52】
JP0005886309B2_000050t.gif

である。次章の固有関数法に用いる関数φ,φは,(主,双対)連立微分方程式(105),(106)を満たす組合せとする。ただし,φはそのままで,λとφの符号を反転させた組合せも要件を満たすが,この組合せは排除する。

【0031】
原初微分作用素Lijと同次境界条件が共に自己随伴な場合には,(105)式と(106)式は一致して,自己随伴連立固有値問題[Self-Adjoint Simultaneous Eigenvalue Problem]
【数53】
JP0005886309B2_000051t.gif
となる。上式はチモシェンコ梁[Timoshenko beam]やミンドリン平板[Mindlin plate]などの固有振動を表す形式で,固有関数法から直接,(123)式を導く鍵となる。

【0032】
5.固有関数法と最小2乗法
5.1固有関数法
5.1.1 主問題
同次化された微分方程式(40)に(69)式を代入して,

【0033】
【数54】
JP0005886309B2_000052t.gif
となる。上式は,主連立微分方程式(105)により,
【数55】
JP0005886309B2_000053t.gif
となる。ここで,λはk次の固有値である。関数φと内積をとって直交性を利用すれば,
【数56】
JP0005886309B2_000054t.gif
となる。これにより,係数cが定まる(104)式の正規化によって,
【数57】
JP0005886309B2_000055t.gif
となり,外力項fHiと双対固有関数φの内積から係数cが定まる事を知る。これを(69)式に戻して変位uHjを得,uHjを(24)式に戻せば変位uを得る。一方,(109)式に戻せば外力項fHiが再構成され検証に役立つ。即ち,(111)式の解法は,変位uHjを主固有関数φで,外力項fHiを双対固有関数φで表現する,ヒルベルト空間での固有関数法である。

【0034】
5.1.2 双対問題
同次化された微分方程式(66)に(70)式を代入して,
【数58】
JP0005886309B2_000056t.gif

となる。上式は,双対連立微分方程式(106)により,
【数59】
JP0005886309B2_000057t.gif

となる。関数φと内積をとって直交性を利用すれば,
【数60】
JP0005886309B2_000058t.gif
となる。これにより,係数cが定まる。(104)式の正規化によって,
【数61】
JP0005886309B2_000059t.gif

【0035】
となり,外力項fHiと主固有関数φの内積から係数cが定まる事を知る。これを(70)式に戻して変位uHjを得,uHjを(42)式に戻せば変位uを得る。一方,(113)式に戻せば外力項fHiが再構成され検証に役立つ。即ち,(115)式の解法は,変位uHjを双対固有関数φで,外力項fHiを主固有関数φで表現する,ヒルベルト空間での固有関数法である。

【0036】
5.1.3 双対原理
ここで改めて図1を熟視すれば,あたかもコインの如く表裏一体の構造は,原初微分作用素がその材質に,同次境界条件が表裏紋様の礎に当り,これに内在する固有値と(主,双対)固有関数が解法の主役を演ずる事が判る。

【0037】
また,双対問題から出発すると逆向きに同様な関係が得られる為,その求解に,(主,双対)固有関数を反転利用できる事も判る。即ち,図1は連立偏微分方程式の双対原理[Duality Principle]を示している。「双対」の名を冠した根拠は,上記の特徴「表裏一体・反転利用」にある。

【0038】
5.2 変分直接法と最小2乗法
主変分[Primal Variation]δuは(24)式,双対変分[Dual Variation]δuは(42)式より,それぞれ,
【数62】
JP0005886309B2_000060t.gif
となる。非同次の境界条件を満たす項は既知関数の為,上式では,
【数63】
JP0005886309B2_000061t.gif
を用いた。結局,(69),(116)式及び,(70),(117)式より,
【数64】
JP0005886309B2_000062t.gif
を得る。これは,主変分δuを主固有関数φで,双対変分δuを双対固有関数φで表現した式に他ならない。よって,(110)式を得る過程で用いた,関数φと内積をとる手順を複数組み合わせる事は,(121)式の双対変分δuと内積をとる事と等価である。即ち,
【数65】
JP0005886309B2_000063t.gif
が成立する。逆に言えば,様々なδuに対して上式を作成し,連立方程式を代数的に解く事で,係数cが定められる。

【0039】
また,自己随伴問題の場合には(107)式を通じて(主,双対)固有関数φ,φは一致するから,(120),(121)式も同等になり,(122)式は,
【数66】
JP0005886309B2_000064t.gif
となる。
一方,(105),(120)式により,主変分δuの微分係数の和について,

【0040】
【数67】
JP0005886309B2_000065t.gif
を得る。成分ごとに同じ重みwを用いた(主,双対)固有関数では,
【数68】
JP0005886309B2_000066t.gif
と定義できるから,(121),(124),(125)式より,双対変分δuを,
【数69】
JP0005886309B2_000067t.gif
と認識できる。(122)式に上式を代入して,
【数70】
JP0005886309B2_000068t.gif
を得る。fを既知外力として扱う場合,上式は,汎関数Πを,
【数71】
JP0005886309B2_000069t.gif
とおいて,その変分がゼロになる事と等価である。即ち,最小2乗法[Least Squares Method]が変分原理[Variational Principle]となる事を示している。なお,(122),(123),(127)式は共に,変分直接法[Direct Variational Method]を示しており,(127)式は,非自己随伴問題で(122)式に等価,自己随伴問題で(123)式に等価になる。

【0041】
6.弾性学のエネルギ原理について
物体力が指定される問題を扱う。(128)式より,基本的には全ての問題に対して,最小2乗法が変分原理となる。

【0042】
(91)式の境界項Rに着目すれば,境界表面で
JP0005886309B2_000070t.gifとなる条件は自己随伴な境界条件である事が判る。即ち,表面力および変位のx,y方向の成分について,各成分ごとにどちらか一方だけを与える問題は自己随伴問題となる。従って,この場合には,(123)式を用いても良い。単位厚さ当りの全ひずみエネルギUは,
【数72】
JP0005886309B2_000071t.gif
で,その変分δUは,
【数73】
JP0005886309B2_000072t.gif
である。(123)式を変形して,上式を用いれば,
【数74】
JP0005886309B2_000073t.gif
を得る。これは,仮想仕事の原理[Principle of Virtual Work]と同じ形式であるが,(123)式そのものであるから,変分直接法として認識すべきである。境界表面から変位指定領域を取り除いた領域をcとすれば,c上で表面力は指定値となるから,上式は,汎関数Πを,
【数75】
JP0005886309B2_000074t.gif
とおいて,その変分がゼロになる事と等価である。これは,最小ポテンシャルエネルギの原理[Principle of Minimum Potential Energy]である。

【0043】
エネルギ保存則を表す,クラペイロンの定理[Clapeyron’s theorem]
【数76】
JP0005886309B2_000075t.gif
において変分をとれば,
【数77】
JP0005886309B2_000076t.gif
となる。上式と(132)式より,
【数78】
JP0005886309B2_000077t.gif
となる。これは,補仮想仕事の原理[Principle of Complementary Virtual Work]と同じ形式であるが,(123)式にエネルギ保存則を適用した結果得られたものであるから,変分直接法の別形式として認識すべきである。境界表面から表面力指定領域を取り除いた領域をCuとすれば,Cu上で変位は指定値となるから,上式は,汎関数Πを,
【数79】
JP0005886309B2_000078t.gif
とおいて,その変分がゼロになる事と等価である。これは,最小仕事の原理[Principle of Minimum Work]である。


【0044】
7.まとめ
(1)連立偏微分方程式(23)は,(111)式の固有関数法によって解ける。(128)式より,最小2乗法は変分原理である。非自己随伴問題では(122),(127)式が変分直接法となり,自己随伴問題では(123),(127)式が変分直接法になる。

【0045】
(2) (123)式の変分δuは,自己随伴固有関数で構成されるから,幾何学的境界条件と力学的境界条件の両方を満足する変位に対する変分である。一方,偶然にも(123)式と同
一形式の仮想仕事の原理では,δuは「幾何学的境界条件を満たす,すべての関数」で
仮想変位と称される。力学的境界条件を満たすとは限らない関数に対して変分記号δを使う事は不適切であるから,仮想変位をδと書くと,仮想仕事の原理とは,
【数80】
JP0005886309B2_000079t.gif
である。これは,重み付き残差法に他ならない。一般に,上式においてδとδuを同一視した変分演算が行われ,力学的境界条件と微分方程式が得られるが,実は,(123
)式による恩恵である。

【0046】
(3) (123)式の恩恵を仮想仕事の原理による効果と捉えて,ひとたび,これを力学体系の根幹に据えてしまうと,自己随伴問題しか扱えなくなる。非自己随伴問題であっても,重み付き残差法としての機能で解けてしまう事がある為,問題点に気づかない。仮想仕事の原理の持つ曖昧な部分を再考し,本章(1),特に(122)式を論拠とする事で,弾性学など各種力学は,より良い体系へと再構築され新たな発展へ向かう事が期待できる。

【0047】
(4) 一般に,線形弾性有限要素は(132)式により定式化されるが,h法要素では変分の演算は行われていない。(132)式から変分記号δを除去して1/2を乗じれば,(134)式のエネルギ保存則と一致するが,h法要素はこれに従っている。(123)式にて変分の演算を行
えばp法と同等の要素が得られ,(122),(127)式からは新たな有限要素を構成できる。

【0048】
8.新しい有限要素
8.1 形状関数
例えば,16節点四角形要素の節点番号を図8のように表すと,その形状関数N(i=1~16)は,

【0049】
【数81】
JP0005886309B2_000080t.gif
の不等号の範囲で境界内部を表し,等号にて境界を表す。

【0050】
8.2 要素内座標
形状関数Ni(i=1~n)を用いて要素内部の座標xjを,
【数82】
JP0005886309B2_000081t.gif

【0051】
8.3 要素内変位
主解となる変位を,境界変位を表現する項と,境界で変位がゼロとなる項とに分けて,
【数83】
JP0005886309B2_000082t.gif
と表す。添字Aが境界変位を表現する項を,添字oが境界で変位がゼロとなる項を表す。まとめて,変位関数u
【数84】
JP0005886309B2_000083t.gif

【0052】
となる。
内挿関数χAi(i=1~n)を用いて要素内部の変位uAj
【数85】
JP0005886309B2_000084t.gif
とかく。Uijは節点変位である。ここでは,アイソパラメトリック要素と同様に,内挿関数(補間関数)χAiを形状関数Nと一致させる。変位の詳細は,
【数86】
JP0005886309B2_000085t.gif
と表す。uAjは要素境界のみならず要素内部の変位をも表現するが,これが微分方程式を満足する保証は無い。要素内部におけるuAjの変位表現能力には限界があるため,要素境界での変位表現能力のみに期待する。その意味でuAjのことを境界関数と呼ぶ。uAjでは表現しきれない要素内部の変位を,uojが補正することを期待する。その意味でuojのことを補正関数と呼ぶ。即ち、(162)式は、変位関数uを境界関数uAjと補正関数uojの和で表す式である。

【0053】
補正関数uojを,境界での変位がゼロとなる試行関数ψok(k=1~1)の和で表現して,
【数87】
JP0005886309B2_000086t.gif
である。
【数88】
JP0005886309B2_000087t.gif
となる。

【0054】
一方,(24)式では,非同次境界条件を満たす項uBjと,同次境界条件を満たす項uHjとの和で,主解となる変位uを表している。行列で表現して,
【数89】
JP0005886309B2_000088t.gif
となる。境界条件を与えた段階で,その詳細は定まる。

【0055】
8.4 試行関数

【数90】
JP0005886309B2_000089t.gif

【0056】
8.5 双対変位
双対解となる双対変位uを,境界変位を表現する項と,境界で変位がゼロとなる項とに分けて,
【数91】
JP0005886309B2_000090t.gif
とかく。Uijは双対節点変位である。ここでは,アイソパラメトリック要素と同様に,内挿関数(補間関数)χAiを形状関数Nと一致させる。変位の詳細は,
【数92】
JP0005886309B2_000091t.gif

【0057】
補正関数uojを試行関数ψok(k=1~l)の和で表現して,
【数93】
JP0005886309B2_000092t.gif
一方,(42)式では,非同次随伴境界条件を満たす項uBjと,同次随伴境界条件を満たす項uHjの和で,双対変位uを表している。行列で表現して,
【数94】
JP0005886309B2_000093t.gif
となる。境界条件を与えた段階で,その詳細は定まる。

【0058】
8.6 等価節点力
ひずみ成分は,(4)~(6)式から,

【0059】
【数95】
JP0005886309B2_000094t.gif
となる。

【0060】
変位の微分係数は,
【数96】
JP0005886309B2_000095t.gif
となる。右辺の行列を,
【数97】
JP0005886309B2_000096t.gif
となる。

【0061】
応力は
【数98】
JP0005886309B2_000097t.gif

【0062】
【数99】
JP0005886309B2_000098t.gif

【0063】
8.7 双対等価節点力
ひずみ成分は,(4)~(6)式から,
【数100】
JP0005886309B2_000099t.gif
と表す。応力成分は,(7)~(9)式から
【数101】
JP0005886309B2_000100t.gif
となる。これは,(10)~(12)式の応力成分を表している。

【0064】
変位の微分係数は,
【数102】
JP0005886309B2_000101t.gif

【0065】
応力は、
【数103】
JP0005886309B2_000102t.gif
である。

【0066】
要素境界面の変位は、
【数104】
JP0005886309B2_000103t.gif

【0067】
8.8 境界上での積分
各辺での積分を実行する具体的な手順は,以下の通り。
【数105】
JP0005886309B2_000104t.gif

【0068】
8.9 ひずみエネルギ
要素内部のひずみエネルギengは,
【数106】
JP0005886309B2_000105t.gif

【0069】
8.10 双対ひずみエネルギ
要素内部のひずみエネルギengは,
【数107】
JP0005886309B2_000106t.gif

【0070】
8.11 随伴境界条件
境界表面での変位は(208),双対変位は(229)式であり,表面力は(207),双対表面力は(228)式である。これらを用いて,(60)式の境界項Rは,
【数108】
JP0005886309B2_000107t.gif
となる。等価節点力(216),(235)式を用いて,
【数109】
JP0005886309B2_000108t.gif

【0071】
となる。節点力Fを既知部Fと未知部Fに分け,節点変位Uを既知部Uと未知部Uに分ける。これに応じて,双対節点変位U*が未知部Uv*と既知部Ub*に分かれ,双対節点力Fが未知部Fと既知部Fに分かれて,
【数110】
JP0005886309B2_000109t.gif
となる。そして,
【数111】
JP0005886309B2_000110t.gif
である。既知部Fb,をゼロに差し替えれば,同次境界条件を与えたことになる。このとき,境界項Rをゼロとする条件より,同次随伴境界条件が得られ,
【数112】
JP0005886309B2_000111t.gif
となる。既知部Fb,bと既知部Fb*, Ub*の組み合わせが一致する場合には自己随伴境
界条件であり,一致しない場合には非自己随伴境界条件である。

【0072】
8.12 主試行関数と主固有関数
等価節点力の定義式(216)を変形して,方程式
【数113】
JP0005886309B2_000112t.gif
を得る.ここで,[I]は単位行列である。上式は,節点変位U、節点力F、係数Cの間に成立せねばならない関係を表している。上式に、(173)式の要素内変位uを組み合わせて、

【数114】
JP0005886309B2_000113t.gif
となる.節点変位Uの既知部Uと、節点力Fの既知部Fを合わせて,節点既知部s
【数115】
JP0005886309B2_000114t.gif
と定義する.その個数nbは,
【数116】
JP0005886309B2_000115t.gif
である。節点変位Uの未知部Uvと,節点力Fの未知部Fvを合わせて,節点未知部sv
【数117】
JP0005886309B2_000116t.gif
と定義する.その個数nvは,
【数118】
JP0005886309B2_000117t.gif
である.また,(253),(259),(261)式より,
【数119】
JP0005886309B2_000118t.gif
である.(257)式のU、Fの順番を、(258)、(260)式に従って並べ直すと、行列部の列も入れ替わって、
【数120】
JP0005886309B2_000119t.gif
となる.ここで,行列[K]は、既知部U、Fの並びに従って行列[K]、[-I]から該当列を抽出して得られるもので,同様に、未知部U,Fの並びに従って行列[K]、[-I]から該当列を抽出すれば,行列[K]が得られる.また,行列[Y]、[Y]は、次の様にして得られる.単位行列[I]とゼロ行列[0]を用いて,(257)式の行列[χ]、[0]を
【数121】
JP0005886309B2_000120t.gif
と認識すれば,行列[χA]から列を抽出する操作は[I]の列を抜き出すことで,同様
に,行列[0]から列を抽出する操作は[0]の列を抜き出すことで表現できる.既知部Ubと未知部Uvの並びに従って、[I]から該当列を抽出して[IUb]、[IUv]を構成する。同様に、既知部Fbと未知部Fvの並びに従って,「0」から該当列を抽出して[0Fb]、[0Fv]を構成する。

【0073】
(264)式に対して、
【数122】
JP0005886309B2_000121t.gif
を定義し,(265)式に対して,
【数123】
JP0005886309B2_000122t.gif
を定義する.さらに,係数行列
【数124】
JP0005886309B2_000123t.gif
を定義すれば、行列[Yb]、[Yv]は、
【数125】
JP0005886309B2_000124t.gif
となる.以上の操作によって,(263)式の行列は定まり,要素内変位uは,
【数126】
JP0005886309B2_000125t.gif
となる.上式で変位uは,節点既知部sb、節点未知部sv、未知係数coで表される状態
である。そして、未知部sv、coと既知部sbの間には、
【数127】
JP0005886309B2_000126t.gif
が成立することを示している。従って、未知部sv、coはもう少し整理できて、未知部の総数を減らせることが判る。

【0074】
そこで、未知部sv、coを並べ直して未知部se、shとし,それに応じて行列[Kv
、[K0]の列を入れ替えて[Ke]、[Kh]とする。同様に,行列[Yv]、[ψ0]の列を入れ替えて[Ye]、[Yh]とする.その結果,(263)式は,
【数128】
JP0005886309B2_000127t.gif
となる.この手続きは,[Ke]が逆行列を持つようにして、未知部seを消去することが狙いである.よって,[Ke]は正方行列とする.未知部shの個数nhは,(263),(276)式より、
【数129】
JP0005886309B2_000128t.gif
である。(262)式により、
【数130】
JP0005886309B2_000129t.gif
でもある。未知部se、shを次の手順で構成する。
(1)未知部svのうち未知部seに入るものnve個を、未知部seの前半に集める。
(2)未知部coのうち未知部seに入るものnoe個を、未知部seの後半に集める。
(3)未知部svのうち未知部shに入るものnvh個を、未知部shの前半に集める。
(4)未知部coのうち未知部shに入るものnoh個を、未知部shの後半に集める。
よって、未知部の個数は、
【数131】
JP0005886309B2_000130t.gif
となる。手順(1)、(3)からは、
【数132】
JP0005886309B2_000131t.gif
となり、手順(2)、(4)からは、
【数133】
JP0005886309B2_000132t.gif
となる。 (277),(279),(280),(281)式から、
【数134】
JP0005886309B2_000133t.gif
となる。これは、未知部seの個数がnとなることを示している。

【0075】
自己随伴境界条件では、
【数135】
JP0005886309B2_000134t.gif
であるから、上式を(277)式あるいは(278)式に代入して、
【数136】
JP0005886309B2_000135t.gif
となる。さらに、(279),(281),(284)式より、
【数137】
JP0005886309B2_000136t.gif
となる。

【0076】
行列[Yu]は(273)式で表されるから,行列[Yv]から列を抽出する操作は,行列[
v]の列を抜き出すことで表現できる。[cv]から該当列を抽出して[cve]、[cvh]を構成して、
【数138】
JP0005886309B2_000137t.gif
を定義する。単位行列[I]を用いて、(263)式の行列[ψ0]を
【数139】
JP0005886309B2_000138t.gif
と認識すれば,行列[ψ0]から列を抽出する操作は,行列[I]の列を抜き出すことで
表現できる.行列[I]から該当列を抽出して[Ioe]、[Ioh]を構成して、
【数140】
JP0005886309B2_000139t.gif
を定義する、よって、行列[Ye]、[Yh]は、
【数141】
JP0005886309B2_000140t.gif
となる。さらに、関数行列
【数142】
JP0005886309B2_000141t.gif
および、係数行列
【数143】
JP0005886309B2_000142t.gif
を定義すれば、行列[Ye]、[Yh]は、
【数144】
JP0005886309B2_000143t.gif
となる。

【0077】
(276)式を分離して、
【数145】
JP0005886309B2_000144t.gif
となる。(298)式の未知部seを解いて、
【数146】
JP0005886309B2_000145t.gif
となり、上式を(299)式に代入して、
【数147】
JP0005886309B2_000146t.gif
となる.第1項は既知部sbを備えた関数であり,非同次境界条件を満足する変位[uB]となる.第2項は未知部shを備えた関数であり,同次境界条件を満足する変位[uH]となる。即ち、
【数148】
JP0005886309B2_000147t.gif
と定義すれば、
【数149】
JP0005886309B2_000148t.gif
であり、
【数150】
JP0005886309B2_000149t.gif
となる。上式は(174)式の詳細を示すもので,変位[u]が、非同次境界条件を満足する変
位[uB]と、同次境界条件を満足する変位[uH]とで表現された事が判る。

【0078】
(272),(296),(302)式より,行列[χB]は
【数151】
JP0005886309B2_000150t.gif
となる。上式第1項は、(293)式より、
【数152】
JP0005886309B2_000151t.gif
となることに注意すれば、
【数153】
JP0005886309B2_000152t.gif
とおいて、
【数154】
JP0005886309B2_000153t.gif
となる。

【0079】
(296),(297),(303)式より,行列[ψh]は
【数155】
JP0005886309B2_000154t.gif
となる。
【数156】
JP0005886309B2_000155t.gif
とおいて、
【数157】
JP0005886309B2_000156t.gif
となる。
(310)、(313)式により、行列[χB]、[ψh]が関数部[Γ]と係数部[Hb]、[Hh]とに分離された事が判る。また、変位{u}が、
【数158】
JP0005886309B2_000157t.gif
と表現された事から,任意の節点既知部sbによって非同次境界条件を満たす関数が構成
され,任意の未知部shによって同次境界条件を満たす関数が構成されることを知る。未
知部shに対して変位[uH]が同次境界条件を満たすのは、行列[ψh]による効果である
から,係数ベクトルを{eh}として,主固有関数{φ}を、
【数159】
JP0005886309B2_000158t.gif
とおくことができる。同時にその変分{δφ}を、
【数160】
JP0005886309B2_000159t.gif
とおくことができる。重要な役割を演ずる行列[ψh]を、主試行関数と呼ぶ。


【0080】
8.13 双対試行関数と双対固有関数
等価節点力の定義式(235)を変形して,方程式
【数161】
JP0005886309B2_000160t.gif
を得る.上式は,双対節点変位U、双対節点力F、係数Cの間に成立せねばならない関係を表している。上式に、(190)式の要素内変位uを組み合わせて、
【数162】
JP0005886309B2_000161t.gif
となる。双対節点変位Uの既知部Uと、双対節点力Fの既知部Fを合わせて,節点既知部s
【数163】
JP0005886309B2_000162t.gif
と定義する.その個数nbは,
【数164】
JP0005886309B2_000163t.gif
である。双対節点変位Uの未知部Uと,双対節点力Fの未知部Fを合わせて,節点未知部s
【数165】
JP0005886309B2_000164t.gif
と定義する.その個数nは,
【数166】
JP0005886309B2_000165t.gif
である.また,(253),(320),(322)式より,
【数167】
JP0005886309B2_000166t.gif
である.(318)式のU、Fの順番を、(319)、(321)式に従って並べ直すと,行列部の列も入れ替わって、
【数168】
JP0005886309B2_000167t.gif
となる.ここで,行列[K]は、既知部U、Fの並びに従って行列[K]、[-I]から該当列を抽出して得られるもので,同様に、未知部U,Fの並びに従って行列[K]、[-I]から該当列を抽出すれば,行列[K]が得られる.また,行列[Y]、[Y]は、次の様にして得られる。単位行列[I]とゼロ行列[0]を用いて,(318)式の行列[χ]、[0]を
【数169】
JP0005886309B2_000168t.gif
と認識すれば,行列[χ]から列を抽出する操作は[I]の列を抜き出すことで,同様に,行列[0]から列を抽出する操作は[0]の列を抜き出すことで表現できる.既知部Uと未知部Uの並びに従って、[I]から該当列を抽出して[IUb*]、[IUv*]を構成する。同様に、既知部Fと未知部Fの並びに従って,[0]から該当列を抽出して[0Fb*]、[0Fv*]を構成する。

【0081】
(325)式に対して、
【数170】
JP0005886309B2_000169t.gif
を定義し,(326)式に対して,
【数171】
JP0005886309B2_000170t.gif
を定義する.さらに,係数行列
【数172】
JP0005886309B2_000171t.gif
を定義すれば、行列[Yb*]、[Yv*]は、
【数173】
JP0005886309B2_000172t.gif
となる.以上の操作によって,(324)式の行列は定まり,要素内双対変位u*は,
【数174】
JP0005886309B2_000173t.gif
となる.上式で双対変位u*は,節点既知部sb*、節点未知部sv*、未知係数co*で表さ
れる状態である。そして、未知部sv*、co*と既知部sb*の間には、
【数175】
JP0005886309B2_000174t.gif
が成立することを示している。従って、未知部sv*、co*はもう少し整理できて、未知部の総数を減らせることが判る。

【0082】
そこで、未知部sv*、co*を並べ直して未知部se*、sh*とし,それに応じて行列[Kv*]、[K0]の列を入れ替えて[Ke*]、[Kh*]とする。同様に,行列[Yv*]、[ψ0]の列を入れ替えて[Ye*]、[Yh*]とする.その結果,(324)式は,
【数176】
JP0005886309B2_000175t.gif
となる.この手続きは,[Ke*]が逆行列を持つようにして、未知部se*を消去することが狙いである.よって,[Ke*]は正方行列とする.未知部sh*の個数nh*は,(324),(337)式より、
【数177】
JP0005886309B2_000176t.gif
である。(323)式により、
【数178】
JP0005886309B2_000177t.gif
でもある。未知部se*、sh*を次の手順で構成する。
(1)未知部sv*のうち未知部se*に入るものnve*個を、未知部se*の前半に集める。
(2)未知部co*のうち未知部se*に入るものnoe*個を、未知部se*の後半に集める。
(3)未知部sv*のうち未知部sh*に入るものnvh*個を、未知部sh*の前半に集める。
(4)未知部co*のうち未知部sh*に入るものnoh*個を、未知部sh*の後半に集める。
よって、未知部の個数は、
【数179】
JP0005886309B2_000178t.gif
となる。手順(1)、(3)からは、
【数180】
JP0005886309B2_000179t.gif
となり、手順(2)、(4)からは、
【数181】
JP0005886309B2_000180t.gif
となる。 (338),(340),(341),(342)式から、
【数182】
JP0005886309B2_000181t.gif
となる。これは、未知部se*の個数がnとなることを示している。

【0083】
自己随伴境界条件では、
【数183】
JP0005886309B2_000182t.gif
であるから、上式を(338)式あるいは(339)式に代入して、
【数184】
JP0005886309B2_000183t.gif
となる。さらに、(340),(342),(345)式より、
【数185】
JP0005886309B2_000184t.gif
となる。

【0084】
行列[Yu*]は(334)式で表されるから,行列[Yv*]から列を抽出する操作は,行列
[cv*]の列を抜き出すことで表現できる。[cv*]から該当列を抽出して[cve*]、
[cvh*]を構成して、
【数186】
JP0005886309B2_000185t.gif
を定義する。単位行列[I]を用いて、(324)式の行列[ψ0]を
【数187】
JP0005886309B2_000186t.gif
と認識すれば,行列[ψ0]から列を抽出する操作は,行列[I]の列を抜き出すことで
表現できる.行列[I]から該当列を抽出して[Ioe*]、[Ioh*]を構成して、
【数188】
JP0005886309B2_000187t.gif
を定義する、よって、行列[Ye*]、[Yh*]は、
【数189】
JP0005886309B2_000188t.gif
となる。さらに、係数行列
【数190】
JP0005886309B2_000189t.gif
を定義すれば、行列[Ye*]、[Yh*]は、
【数191】
JP0005886309B2_000190t.gif
となる。

【0085】
(337)式を分離して、
【数192】
JP0005886309B2_000191t.gif
となる。(358)式の未知部se*を解いて、
【数193】
JP0005886309B2_000192t.gif
となり、上式を(359)式に代入して、
【数194】
JP0005886309B2_000193t.gif
となる.第1項は既知部sb*を備えた関数であり,非同次随伴境界条件を満足する変位[
B*]となる.第2項は未知部sh*を備えた関数であり,同次随伴境界条件を満足する変
位[uH*]となる。即ち、
【数195】
JP0005886309B2_000194t.gif
と定義すれば、
【数196】
JP0005886309B2_000195t.gif
であり、
【数197】
JP0005886309B2_000196t.gif
となる。上式は(191)式の詳細を示すもので,変位[u*]が、非同次随伴境界条件を満足する変位[uB*]と、同次随伴境界条件を満足する変位[uH*]とで表現された事が判る。

【0086】
(333),(356),(362)式より,行列[χB*]は
【数198】
JP0005886309B2_000197t.gif
となる。上式第1項は、(293)式より、
【数199】
JP0005886309B2_000198t.gif
となることに注意すれば、
【数200】
JP0005886309B2_000199t.gif
とおいて、
【数201】
JP0005886309B2_000200t.gif
となる。

【0087】
(356),(357),(363)式より,行列[ψh*]は
【数202】
JP0005886309B2_000201t.gif
となる。
【数203】
JP0005886309B2_000202t.gif
とおいて、
【数204】
JP0005886309B2_000203t.gif
となる。
(370)、(373)式により、行列[χB*]、[ψh*]が関数部[Γ]と係数部[Hb*]、[Hh*]とに分離された事が判る。また、変位{u*}が、
【数205】
JP0005886309B2_000204t.gif
と表現された事から,任意の節点既知部sb*によって非同次随伴境界条件を満たす関数が構成され,任意の未知部sh*によって同次随伴境界条件を満たす関数が構成されることを知る。未知部sh*に対して変位[uH*]が同次随伴境界条件を満たすのは、行列[ψh*]による効果であるから,係数ベクトルを{eh*}として,双対固有関数{φ*}を、
【数206】
JP0005886309B2_000205t.gif
とおくことができる。同時にその変分{δφ*}を、
【数207】
JP0005886309B2_000206t.gif
とおくことができる。重要な役割を演ずる行列[ψh*]を、双対試行関数と呼ぶ。

【0088】
8.14 個数の関係
(259),(320),(252)式より,
【数208】
JP0005886309B2_000207t.gif
となる。(278)、(339)、(377)式より、
【数209】
JP0005886309B2_000208t.gif
となる。

【0089】
8.15 主固有関数と双対固有関数の決定
主連立微分方程式(105)、双対連立微分方程式(106)を行列で表現して、それぞれ、
【数210】
JP0005886309B2_000209t.gif
となる。そこで、微分作用素を
【数211】
JP0005886309B2_000210t.gif
とおき、重みを
【数212】
JP0005886309B2_000211t.gif
とおいて、(主、双対)連立微分方程式は、
【数213】
JP0005886309B2_000212t.gif
となる。(122)式の変分直接法に従って、変分{δφ*}、{δφ}を乗じて積分すれば、
【数214】
JP0005886309B2_000213t.gif

となる。(315),(375)式の固有関数{φ}、{φ*}及び(316)、(376)式の変分{δφ}、{δφ*}を代入して、
【数215】
JP0005886309B2_000214t.gif
となる。左辺の積分結果を、
【数216】
JP0005886309B2_000215t.gif
とおき、右辺の積分結果を、
【数217】
JP0005886309B2_000216t.gif
とおく。[Bh]、[Bh*]は、共に対称行列となる。また、主試行関数[ψh]は同次境界条件を満たし、双対試行関数[ψh*]は同次随伴境界条件を満たすから、部分積分により、
【数218】
JP0005886309B2_000217t.gif
となる。即ち、
【数219】
JP0005886309B2_000218t.gif
である。(390)~(393)式を用いて(388),(389)式は
【数220】
JP0005886309B2_000219t.gif
となる。これが、任意の変分に対して成立するためには、
【数221】
JP0005886309B2_000220t.gif
でなくてはならない。(398)、(399)式を相互に代入して連立をはずせば、
【数222】
JP0005886309B2_000221t.gif
となる.これを解くことにより,固有値λと固有ベクトル{eh}、{eh*}を得ること
ができる。逆に,(398),(399)式をひとつにまとめると、
【数223】
JP0005886309B2_000222t.gif
【数224】
JP0005886309B2_000223t.gif
である.なお,(403)式の等号については,(395)式による。

【0090】
これを解くことにより,固有値λと固有ベクトル{eh}、{eh*}を得ることができ
る。[As]、[Bs]は共に対称行列であるから、固有値λは実数となり,固有ベクトル{es}は[As]、[Bs]を挟んで直交する。

【0091】
(400),(401)式で計算すれば,固有値問題のサイズは小さくて済む.その反面,固有ベ
クトル{eh}、{eh*}の組み合わせについて,注意深く取り扱う必要が生じる.一方
,(402)式で計算すれば,固有値問題のサイズは大きくなるが,固有ベクトル{eh}、{eh*}の組み合わせも同時に得られるので、組み合わせの確認は不要となる。どちらの方法を用いるかは、メモリサイズや計算速度に応じて決めれば良いが、定式化の都合からは、(402)式の方が扱い易いので、以降は、これを考え方の基準にする。

【0092】
(402)式からは,2l個の固有値λと固有ベクトル{es}の組み合わせが得られる。正の固有値と負の固有値の個数は一致して、これをmp個とする.これは,(398),(399)式で{eh}はそのままで、λと{eh*}の符号を反転させた組み合わせも連立方程式を満た
す事によるもので、4.4 節で述べた事でもある。[As]がゼロ空間を持つ場合には,固
有値がゼロとなる。[Ah]のゼロ空間をm0 個とする。

【0093】
固有値λに属する固有ベクトル{es}の成分{eh}、{eh*}を、(315),(375)式に
戻すことにより,固有関数{φ}、{φ*}が得られる。

【0094】
正の固有値に属するmp個の固有関数{φ}、{φ*}及び固有ベクトル{eh}、{eh*}を列方向に並べて、それぞれ、[φp]、[φp*]及び[ep]、[ep*]と表せば、
【数225】
JP0005886309B2_000224t.gif
となる。mp個の固有値を対角上に並べて、
【数226】
JP0005886309B2_000225t.gif
を定義すれば、(主、双対)連立微分方程式は、
【数227】
JP0005886309B2_000226t.gif
となる。固有関数の直交性により、[ep]、[ep*]が、
【数228】
JP0005886309B2_000227t.gif
を満たすように正規化しておけば、
【数229】
JP0005886309B2_000228t.gif
となる。

【0095】
ゼロ固有値に属するm0個の固有関数{φ}及び固有ベクトル{eh}を列方向に並べて、それぞれ、[φ0][e0]と表せば、
【数230】
JP0005886309B2_000229t.gif
となる。固有関数の直交性により、
【数231】
JP0005886309B2_000230t.gif
となる。

【0096】
8.16 行列の詳細
(390)式の[Ah]を、(293),(313),(373)式により変形すれば、
【数232】
JP0005886309B2_000231t.gif
となる.(391)式の[Ah]を,(293),(313),(373)式により変
形すれば,
【数233】
JP0005886309B2_000232t.gif
となる。(392)式の[Bh]を、(293)、(313)式により変形すれば、
【数234】
JP0005886309B2_000233t.gif
となる。[Bh]が対称行列となる事を如実に示している。(393)式の[Bh*]を、(293)
、(373)式により変形すれば、
【数235】
JP0005886309B2_000234t.gif
となる。[Bh*]が対称行列となる事を如実に示している。

【0097】
8.17 固有関数法
主問題の連立偏微分方程式(23)を行列で表現して、
【数236】
JP0005886309B2_000235t.gif
である。変位{u}は(314)式で表されるから、上式は、
【数237】
JP0005886309B2_000236t.gif
となる。一方,変位{uH}を固有関数[φp]、[φ0]で表すと、その係数ベクトルを
{ap}、{a0}と定義して、(406),(415)式を利用すれば、
【数238】
JP0005886309B2_000237t.gif
となる。上式と(305)式を比較して、未知部shは、
【数239】
JP0005886309B2_000238t.gif
となる。即ち、係数ベクトル{ap}、{a0}が定まれば、未知部shが決まることを示
している。
(422)、(423)式より、連立偏微分方程式は、
【数240】
JP0005886309B2_000239t.gif
となる。両辺に[φp*Tを乗じて内積をとれば、(416)式より、左辺第2項はゼロとなり、
【数241】
JP0005886309B2_000240t.gif
となる。これは、(340)、(407)、(413)式により、
【数242】
JP0005886309B2_000241t.gif
となる。整理して、
【数243】
JP0005886309B2_000242t.gif
となる。(293)、(310)、(373)式を用いて、右辺の積分結果を、
【数244】
JP0005886309B2_000243t.gif
と計算しておけば
【数245】
JP0005886309B2_000244t.gif
である。既知部sbに応じて,係数ベクトル{ap}が定まることがわかる。

【0098】
上式を(424)式に戻せば,ゼロ固有値が存在しない場合には未知部shが定まるので、問題が完全に解けたことになる。ゼロ固有値が存在する場合には、未知部shは係数ベクトル
{a0}の関数となるので,節点情報の他に,辺上の境界条件を加えることで{a0}を定めることができる。

【0099】
このようにして,未知部shが定まれば,(300)式により未知部seが定まる。よって、(299)式あるいは(301)式より、変位{u}が完全に定まる。変位が定まれば、要素内部の
応力分布、ひずみエネルギ分布なども得られる。

【0100】
8.18 ゼロ固有値の取り扱い
固有値がゼロになるのは、(384)、(385)式の(主、双対)連立微分方程式にて、
【数246】
JP0005886309B2_000245t.gif
か、あるいは、
【数247】
JP0005886309B2_000246t.gif
のどちらかが成立する場合である。

【0101】
(415)式や(423)式では、(432)式が成立するものとしている。これは、同次の微分方程
式を意味しているので,その解に付随の係数{a0}は任意である。前節では、その任意
性を辺上の境界条件を付加することで決定する方法を示したが,計算を行う設計者が自ら決る方法もある。係数{a0}に応じてどのような変形状態,あるいは応力状態となるか
を図で表示するシステムを具備して,この図を見ながら実現可能な境界条件となるように,設計者が自ら係数{a0}を決めるのである。


【0102】
9.従来有限要素の利用
9.1 全体系の方程式
7 章の(4)項にて,従来の線形弾性有限要素,特にh法要素では変分の演算は行われて
おらず,(134)式のエネルギ保存則と一致することを述べた。h法要素を複数個組み合わ
せて,全体系の運動方程式は、
【数248】
JP0005886309B2_000247t.gif
と表され,各節点の変位と外力のうち、どちらか一方のみを与えること、即ち自己随伴な境界条件を前提として解法が組み立てられている。例えば,全自由度n個のうち、n1
が変位指定領域、n2個が外力指定領域で,
【数249】
JP0005886309B2_000248t.gif
とするとき、剛性行列[K]、節点変位{U}、節点外力{F}の成分を入れ替えて、
【数250】
JP0005886309B2_000249t.gif
と表す。既知量は{U1}、{F2}であり、未知量は{U2}、{F1}である。上式を分解して、
【数251】
JP0005886309B2_000250t.gif
および、
【数252】
JP0005886309B2_000251t.gif
と表す。(438)式を変形して、
【数253】
JP0005886309B2_000252t.gif
となり、未知量{U2}が定まる。これを(437)式に代入して、
【数254】
JP0005886309B2_000253t.gif
となり、未知量{F1}が定まる。即ち、全ての未知量が得られた。このような手順は、
有限要素法の常套手段として確立されている。

【0103】
しかし、変分の演算を含まず、エネルギ保存則から等価な方程式として(434)式が得ら
れていることに着目すれば、何も自己随伴な境界条件に縛られる必要はなく、非自己随伴な境界条件に対して上式を解いても、良いはずである。

【0104】
9.2 非自己随伴な境界条件での計算
全体系の運動方程式(434)を,単位行列[I]を用いて,
【数255】
JP0005886309B2_000254t.gif
と認識すれば、
【数256】
JP0005886309B2_000255t.gif
となる。節点変位Uの既知部Uと,節点外力Fの既知部Fを合わせて,節点既知部s
【数257】
JP0005886309B2_000256t.gif
と定義する.節点変位Uの未知部Uと,節点外力Fの未知部Fを合わせて、節点未知部s
【数258】
JP0005886309B2_000257t.gif
と定義する。このとき、
【数259】
JP0005886309B2_000258t.gif
である.(442)式のU、Fの順番を、(443),(444)式に従って並べ直すと,行列部の列も入れ替わって、
【数260】
JP0005886309B2_000259t.gif
となる.ここで,行列[Kb]は,既知部Ub、Fbの並びに従って行列[K]、[-I]
から該当列を抽出して得られるもので、同様に、未知部Uv、Fvの並びに従って行列[K]、[-I]から該当列を抽出すれば,行列[Kv]が得られる。上式からは、
【数261】
JP0005886309B2_000260t.gif
が得られる。右辺は既知部のみから構成されるから、左辺の節点未知部svを得ることが
できる。即ち、行列[Kv]が、解の複数存在する劣決定系になるのか、解が存在しない
可能性のある優決定系になるのかに応じて、線形代数学の知見に基づき、未知部svを得
ることができる。解が複数存在する場合には、8.18節のゼロ固有値の場合に相当するので、設計者自身によって解を決定することができる。

【0105】
9.3 重力荷重のある場合
全体系の運動方程式(434)の節点外力Fを、重力荷重など物体力を表す項Fgと表面力を表す項Fsとに分けて、
【数262】
JP0005886309B2_000261t.gif
とおくと、(434)式は、
【数263】
JP0005886309B2_000262t.gif
となる。表面力Fsを移項して、
【数264】
JP0005886309B2_000263t.gif
となる。(446)式と同じ手順により上式は、
【数265】
JP0005886309B2_000264t.gif
となる。左辺に未知部を残すよう変形すると、
【数266】
JP0005886309B2_000265t.gif
となる。前節と同様に、線形代数学の知見に基づき、未知部svを得ることができる。

【0106】
9.4 解が複数存在する場合
解が複数存在する場合には,(452)式を満足する1つの特解{sp}と、(452)式の右辺
をゼロとおいた同次解{s0}が得られる。同次解の個数n0は、
【数267】
JP0005886309B2_000266t.gif
となる。同次解{s0}を列方向に並べた行列を[ψ0]とする。
【数268】
JP0005886309B2_000267t.gif
同次解に対するモード係数を{a0}とおいて、解{sv}は
【数269】
JP0005886309B2_000268t.gif
となる。

【0107】
さて、これで全ての節点変位Uと節点外力Fが定まったので,非自己随伴な問題も解けたことになるが、いずれもモード係数{a0}を含んだ形式で、任意の{a0}に対して解となる。

【0108】
しかし,モード係数{a0}を適切に選ばないと、従来型の有限要素では要素が折れ曲
がった状態となって,不正な解析結果と見なされることが良くある。これは、全体系として,エネルギ保存則と等価な方程式(434)を満たすものの、要素内部で微分方程式を満足
しないことに起因する。

【0109】
そこで,全体系の変形がモード係数{a0}によって支配されることに注目して、要素
の微分方程式に等価な変分原理、即ち、(122)式~(128)式の最小自乗法を用いることで、適切な{a0}を決める.適切な{a0}の周辺の値で要素が折れ曲がった状態にならない場合には、まずまず正確な解析結果が得られていると考えることができる。

【0110】
9.5 最小自乗法の利用
従来型の有限要素では,要素内部の変位は(166)式の{uA}で表現される。{uA}は
節点変位Uで表され,Uは(455)式のようにモード係数{a0}によって表される。

【0111】
要素内部は微分方程式で支配されるから,(421)式の{u}に{uA}を代入することで,モード係数{a0}を決めることができれば望ましい。(128)式の汎関数を全体系について作成し,モード係数{a0}で微分することで(127)式を構成できる。その結果、モード係数{a0}を決定できる。

【0112】
モード係数{a0}をスライダーに連動させて、設計者がスライダーによってモード係
数を調節するのに合わせて、変形図が表示されるように構成すれば、どのような支持方法があり得るのか検討できて、設計者のイメージを湧き立たせることが可能となる。実現可能な変形状態を絞り込んで境界条件を再定義することで、より良い構造を得るための検討が可能となる。


【0113】
以下、図面を参照して、本発明の実施の一形態を詳細に説明する。

【0114】
[第1の実施形態]
図2は、本実施形態の解析装置の構成の一例を示す機能ブロック図である。図2に示す解析装置10(解析システム10)は、演算部1を備えるコンピュータにより構成することができる。演算部1は、解析対象物の設計データ、境界条件データ及び解析対象に対する微分作用素を含む微分方程式データを入力し、解析対象の微分方程式の解を計算して、解析結果データとして出力です。ここで、演算部1は、双対変位ui*の変分のことを双対変分δui*とするとき,下記式によって解ujを計算する。下記式は、上述した式(12
2)と同じである。

【0115】
【数1】
JP0005886309B2_000269t.gif

【0116】
例えば、演算部1は、設計データとして、解析対象物の形状、材質を示すデータを入力し、設計データに基づいて、計算に用いる原初微分作用素Lijのデータを読み込む。原初微分作用素Lijを用いて、上記式で表される微分方程式の解を計算する。これにより非自己随伴問題の解を計算することが可能になる。

【0117】
上記式(数1)において、一例として、fiは解析対象へ働く力、uiは解析対象の変位
とすることができるが、fi及びuiは、これに限られない。解析対象及び微分作用素に応じて、fi及びuiは力及び変位以外のものを示す場合がありうる。例えば、構造解析、振動解析、電磁場解析、材料解析、流体解析、温度解析、音場解析、電磁気学解析、又は回路シミュレータ等の様々な解析において、解析対象に応じた適切な値がfi及びuiに設定することができる。

【0118】
(変形例1)
演算部1は、下記式によって解ujを計算してもよい。下記式は上述した式(123)
と同様である。この場合、自己随伴問題の解を計算することができる。

【0119】
【数2】
JP0005886309B2_000270t.gif

【0120】
(変形例2)
演算部1は、変分直接法により下記式において解析対象の解ujを計算することもでき
る。下記式は、上記(127)と同様である。

【0121】
【数3】
JP0005886309B2_000271t.gif

【0122】
前記変分直接法による解ujの計算は、下記式の汎関数Πにおいて、変分がゼロとなる
解を計算することができる。

【0123】
【数4】
JP0005886309B2_000272t.gif

【0124】
[第2の実施形態]
図3は、第2の実施形態における解析装置の構成の一例を示す機能ブロック図である。図3に示す解析装置10aは、設定部11、随伴境界条件演算部12、及び非自己随伴演算部13を備える。

【0125】
設定部11は、解析対象の原初微分作用素及び変位の境界条件を設定する。例えば、設定部11は、ユーザから、設計データとして解析対象の構成または材料を表すデータの入力を受け付けて、解析対象の構成または材料に基づいて原初微分作用素を決定することができる。また、設定部11は、境界条件データの入力をユーザから受け付けることができる。

【0126】
随伴境界条件演算部12は、設定部11で受け付けた境界条件から随伴境界条件を計算する。例えば、随伴境界条件演算部12は、設定部11が読み込んだ微分方程式データに含まれる原初微分作用素を含む微分方程式に、双対固有関数φi*または双対変位ui*を乗じた積分和すなわち内積を、部分積分して得られる境界項から随伴境界条件を計算することできる。

【0127】
例えば、解析対象が2次元弾性体の場合、上記式(91)を用いて随伴境界条件を計算することができる。この場合、上記式(91)において、uHi、pHiに課された同次境界条件をφi、pEiに課した場合に、境界項がゼロになるようなφi*、pEi*の条件を随伴境界条件として求めることができる。一例として、解析対象の境界表面に、力がゼロ(pEi=0)の境界条件が課された場合、境界項REがゼロになるには、φi*は任意でよく、pEi*及びφiのいずれかがゼロにならなければいけない。ここで、解析対象の境界表面に変
位もゼロ(φi=0)の境界条件が課されていた場合、pEi*は任意となる。この場合、随伴境界条件は、境界表面の変位、力とも任意となる。このような場合、境界条件と随伴境界条件が一致しないため、非自己随伴問題となるただし、上記式(106)の作用により,Lji*φj*がφiに課された境界条件を満たさねばならない。これは前記双対境界条件に含まれる。非自己随伴問題の場合には、前記主境界条件と前記双対境界条件とを利用することにより、適切に問題が解ける。

【0128】
非自己随伴演算部13は、原初微分作用素から主連立微分方程式及び双対連立微分方程式を計算し、主連立微分方程式及び双対連立微分方程式、並びに、前記境界条件及び前記随伴境界条件を用いて、主固有関数及び双対固有関数を求めることにより、連立微分方程式の解を計算する。

【0129】
例えば、原初微分作用素をLijとすると、主連立微分方程式は上記式(105)、双対連列微分方程式は上記式(106)と定義することができる。また、上記式(90)により随伴微分作用素Lij*を計算することができる。非自己随伴演算部13は、主連立微分
方程式(上記式(105))、双対連立微分方程式(上記式(106))、境界条件、随伴境界条件を用いて、主連立固有値問題(上記式(98))と双対連立固有値問題(上記式(99))を定義して、主固有関数φ、及び双対固有関数φ*を計算することができる


【0130】
なお、非自己随伴演算部13の解法は、主連立固有値問題(上記式(98))と双対連立固有値問題(上記式(99))を用いたものに限定されない。例えば、上記式(122)及び式(123)を用いて主固有関数φ、及び双対固有関数φ*を計算することもでき
る。その場合、主固有関数φ、及び双対固有関数φ*は、解析的にも近似的にも計算する
ことができる。例えば、後述するように、上記式(122)及び式(123)によって近似固有関数を計算することができる。

【0131】
非自己随伴演算部13は、主固有関数φ、及び双対固有関数φ*を用いて、主問題の連
立方程式の解を計算する。例えば、上記式(111)で示される係数ckを算出し、これ
を上記式(69)に代入して、解関数uHjを計算することができる。また、非自己随伴演算部13は、上記式(115)により双対問題の係数ck*を計算でき、これを上記(70)に代入することで双対問題の解関数uHj*を計算することができる。

【0132】
(動作例)
図4は、図3に示す解析装置10aの動作例を示すフローチャートである。図4に示す例では、まず、設定部11は、ユーザからの設計データの入力を受け付ける(S1)。設計データは特定の対象や形式に限定されない。例えば、解析対象が固体の場合は、解析対象の形状や材料を表すデータの入力を受け付けることができる。解析対象が流体の場合は
、流体の物性を表すデータの入力を受け付けることができる。また、有限要素法による解析の場合は、要素ごとの構成や物性を示すデータの入力を受け付けてもよい。また、設定部11は、ユーザからの入力のみならず、例えば、アクセス可能な記録媒体から設計データを読み込む、あるいはネットワークを介するダウンロード等により、設計データを取得することができる。

【0133】
設定部11は、S1で取得した設計データに基づいて、原初微分作用素を決定し、主問題の微分方程式を決定する(S2)。例えば、設定部11は、予め記録されたデータから、設計データで示される解析対象の形状に応じた微分作用素を抽出し、抽出した微分作用素を含む微分方程式を主問題の微分方程式として設定することができる。このように、設計データに応じて主問題の微分方程式を決定することができる。ここでは、一例として、原初微分作用素がLijである場合について説明する。

【0134】
設定部11は、さらに、境界条件データを取得する(S3)。境界条件データは、解析対象の境界条件に関するデータであり、特定の対象や形式に限定されない。例えば、解析対象物の境界表面における変位及び力の少なくとも一方の条件を示すデータを境界条件として入力することができる。本願開示の発明では、従来は不可能であった、境界表面における変位及び力の双方を設定することができる。

【0135】
随伴境界条件演算部12は、S2で設定された微分方程式と、S3で取得された境界条件データを用いて、随伴境界条件データを生成する。随伴境界条件は、例えば、上述したように、上記式(91)を用いて計算することもできる。また、例えば、有限要素法を用いた例では、上記式(60)で表される境界項を、境界変位、双対変位、表面力、双対表面力をそれぞれ示す上記式(208)、式(229)、式(207)、式(228)を用いて表すことにより計算することができる。

【0136】
非自己随伴演算部13は、原初微分作用素Lijから主連立微分方程式(例えば、上記式(105))及び双対連立微分方程式(上記式(106))を示すデータを生成する(S5)。そして、非自己随伴演算部13は、S5で計算した主連立微分方程式及び双対連立微分方程式、並びに、S2で取得した境界条件及び、S4で計算した随伴境界条件を用いて、主固有関数φ、及び双対固有関数φ*を計算する(S6)。

【0137】
非自己随伴演算部13は、S6で計算した主固有関数φ、及び双対固有関数φ*を用い
て、固有関数法により、主問題の連立方程式の解関数uHjを、例えば、上記式(111)及び式(69)により、計算する(S7)。また、非自己随伴演算部13は、主固有関数φ、及び双対固有関数φ*を用いて、双対問題の連立方程式の解関数uHj*を、例えば、上記式(115)及び式(70)により、計算することができる。このようにして、計算された解は、解析結果をコンピュータのディプレイ(図示せず)などに出力される(S8)。解析結果の具体例については、後述する。

【0138】
上記処理によれば、境界条件に基づいて随伴境界条件を求め、さらに、原初微分作用素から主微分作用素及び双対微分作用素を計算して、これらの主固有関数及び双対固有関数を求めることが可能になる。これらの固有関数を用いて、主問題の連立方程式の解を算出することにより、非自己随伴問題の解の算出が可能になる。なお、自己随伴問題であっても、上記処理によって解を算出することは可能である。

【0139】
以上、解析装置10aの動作例について説明したが、解析装置10aの処理は上記例に限定されない。なお、本実施形態は、上記第1の実施形態の具体例と一つと言える。

【0140】
[第3の実施形態]
図5は、第3の実施形態における解析装置の構成の一例を示す機能ブロック図である。図5に示す解析装置10bは、設定部11、随伴境界条件演算部12、自己随伴判定部15、非自己随伴演算部13、及び、自己随伴演算部14を備える。設定部11、随伴境界条件演算部12、及び非自己随伴演算部13は、上記第1の実施形態の同様に構成することができる。

【0141】
解析装置10bは、境界条件と、随伴境界条件演算部が計算した随伴境界条件とが一致するか否かを判定する自己随伴判定部15を備える。自己随伴判定部15が境界条件と前記随伴境界条件が一致すると判定した場合、自己随伴演算部14は、原初微分作用素から自己随伴問題の自己随伴固有関数を求めることにより、前記自己随伴問題の解を計算する。境界条件と前記随伴境界条件が一致しないと判定された場合、非自己随伴演算部13が、主連立微分方程式及び双対連立微分方程式、並びに、前記境界条件及び前記随伴境界条件を用いて、主固有関数及び双対固有関数を求めることにより、連立微分方程式の解を計算する。

【0142】
自己随伴判定部15は、境界条件と随伴境界条件が一致するか否かを判定することにより、境界条件及び随伴境界条件が自己随伴か非自己随伴かを判定することができる。すなわち、原初微分作用素が自己随伴であって、境界条件と随伴境界条件が一致する場合は自己随伴であり、一致しない場合は非自己随伴であると判断することができる。また、自己随伴判定部15による、境界条件と随伴境界条件が一致するか否かを判定については、原初微分作用素Lijと、随伴微分作用素Ljiが一致するか否か(Lij=Ljiか否か)の判定も含まれる。例えば、随伴微分作用素Lijは、原初微分作用素を含む微分方程式に双
対固有関数φi*または双対変位ui*を乗じた積分和すなわち内積を、部分積分することにより計算することができる。例えば、上記式(62)又は式(90)により、随伴微分作用素Lij*を計算することができる。

【0143】
自己随伴演算部14は、原初微分作用素から自己随伴問題の連立微分方程式を生成し、この連立微分方程式の解を計算する。例えば、上記式(107)又は式(123)を自己随伴問題の連立微分方程式として、自己随伴固有関数φiを計算することができる。自己
随伴演算部14は、例えば、自己随伴固有関数φikを上記式(111)のφik*と置き換
えた式により、ckを計算し、このckを上記式(69)に代入することにより、自己随伴問題の解関数uHjを求めることができる。

【0144】
(動作例)
図6は、図5に示す解析装置10bの動作例を示すフローチャートである。図6において、解析装置10bは、S1~S4の処理を、上記第2の実施形態におけるS1~4の処理と同様に実行することができる。S21において、自己随伴判定部15は、S1で取得された境界条件を示すデータと、S4で生成された随伴境界条件データとを比較することにより、境界条件が随伴境界条件と一致するか否かを判定する。S21で自己随伴判定部15が、一致しないと(すなわち、非自己随伴である)と判定すると、非自己随伴演算部13は、S5~S7の処理を実行する。S5~S7の処理は、上記第2の実施形態におけるS4~S7と同様にしてもよい。

【0145】
S21で自己随伴判定部15が、境界条件と非境界条件が一致する(すなわち、自己随伴である)と判定すると、自己随伴演算部14は、原初微分作用素から自己随伴微分方程式データを生成する(S22)。例えば、自己随伴演算部14は、上記式(107)を自己随伴微分方程式とすることができる。この場合、自己随伴演算部14は、S3で取得した原初微分作用素Lijを自己随伴微分作用素とし、S2で設定した境界条件を自己随伴境界条件として、自己随伴固有関数φiを計算する。そして、自己随伴演算部14は、自己
随伴固有関数φiを用いて、自己随伴連立微分方程式の解関数ujを計算する(S24)
。S8の結果出力は、上記第2の実施形態におけるS8と同様にしてもよい。

【0146】
以上の処理により、解析対象の問題が自己随伴か非自己随伴を判断し、判断結果に応じて適切に処理を切り替えることができる。なお、本実施形態は、上記第1又は第2の実施形態の具体例と一つと言える。

【0147】
[第4の実施形態]
本実施の形態は、解析装置を、有限要素法(FEM)を用いた解析に適用した場合の例である。本実施形態では、原初微分作用素及び境界条件は、有限要素に対して設定される。例えば、上記の第1~第3の実施形態で説明した解析装置の構成及び解析処理を、有限要素に対しても、同様に、適用することができる。

【0148】
ここで、上記第1~第3の実施形態で示した解析処理を有限要素に用いた形態である本実施形態が、従来の有限要素法とどのように異なるかについて説明する。まず、従来のh法では、上述したように、変分演算を行っていないのに対して、本実施形態では、変分直接法により、解を計算している(例えば、上記第1の実施形態)。p法は、高次の関数により、要素内変位を近似する技術であり、要素内変位を適切に近似する関数が手探りで見つけ出されて用いられる。これに対して、本実施形態では、主固有関数及び双対固有関数が計算され、これらをもって、要素内の変位が適切に説明しようとするものである。

【0149】
図7は、本実施形態における解析装置の動作例を示すフローチャートである。本実施形態における解析装置の構成は、図2または図5と同様とすることができる。なお、図7に示すフローチャートは一例であり、本願発明は、これに限られない。図7に示す例では、まず、設定部11が、解析対象の形状を示すデータの入力を受け付ける(S31)。例えば、有限要素の形状を示すデータが入力される。また、設定部11は重み定数を作成する(S32)。重み定数wは、固有値λが無次元量として定まるように使うのが好ましい。例えば、解析対象物が円環の場合、微分作用素の階数が2階なので、代表長さをrとすれ
ば、w=1/r2とするのが好ましい。ここで、代表長さには、物体のどこか特徴的な部分の長さが採用される。円環の場合には、半径を代表長さにとるのが好ましい。この場合、相似図形の固有値λが一致するので、相似図形の物体については、計算を節約することができる。

【0150】
また、設定部11は、節点変位と節点力を入力する(S33)。例えば、節点の変位u及び双対変位u*は、要素境界の変位を表す境界関数uAと補正関数u0で表し、補正関数
0は、試行関数ψ0の和で表すことができる(例えば、上記式(173)及び式(190))。設定部11は試行関数を設定する(S34)。試行関数は、例えば、上記式(179)のように設定することができる。例えば、図8に示すような16節点四角形の要素の場合、境界表面での変位は上記式(208)、双対変位は上記式(229)、表面力は上記式(207)、双対表面力は上記式(228)で表すことができる。

【0151】
随伴境界条件演算部12は、境界条件から随伴境界条件を計算する。例えば、原初微分作用素を含む微分方程式に、双対変位ui*を乗じた積分和すなわち内積を、部分積分して得られる境界項Rに、要素の変位及び表面力を代入して、随伴境界条件を計算することできる。

【0152】
例えば、境界項Rを上記式(250)で表すことができる。この場合、式(250)において、節点力Fの行列要素を既知部Fbと未知部Fvに分けて並べかえ、節点変位Uを既知部分Ubと未知部Uvに分けて並べかえる。この並べかえた、節点力Fb、Fv及び節点変位Ub、Uvのベクトル要素に対応するように、双対節点力F*も既知部Fb*と未知部Fv*
に、双対変位U*も既知部Ub*と未知部Uv*に分けてならべかえると式(251)のよう
になる。Rをゼロとするには、未知部Fv、Uvに対応するUb*、Fb*がゼロでなければならないので、式(254)の随伴境界条件が得られる。

【0153】
非自己随伴演算部13は、例えば、主試行関数[ψh]及び双対試行関数[ψh*]を生成し
(S37、S38)、近似主固有関数を{φ}、近似双対固有関数を{φ*}、それぞれ
の係数ベクトルを{eh}、{eh*}とおいて、
{φ}=[ψh]{eh
{φ*}=[ψh*]{eh*
と表し、これらを主連立微分方程式(上記式(105))及び双対連立微分方程式(上記式(106))へ代入して、上記式(122)及び式(123)によって係数ベクトル{eh}および{eh*}を得た後、近似固有関数{φ}および{φ*}を得ることができる(S39)。

【0154】
一例として、主試行関数[ψh]は、上記式(313)のように表すことができる。この
場合、関数行列Γを上記式(293)で定義し、係数行列[ce][ch]を上記式(294)(295)で定義している。双対試行関数[ψh*]は、例えば、上記式(373)のように表すことができる。この場合、係数行列[ce*][ch*]を上記式(354)(355)で定義している。

【0155】
非自己随伴演算部13は、S39で計算した主固有関数及び双対固有関数を用いて、固有関数法により、主問題の連立方程式の解関数を、計算する(S40)。このようにして、計算された解は、解析結果をコンピュータのディスプレイ(図示せず)などに出力される(S41)。解析結果の具体例については、後述する。

【0156】
S40においては、例えば、(数1)に従って、S37で計算した近似固有関数{φ}の変分{δφ}と、S38で計算した近似双対固有関数{φ*}の変分{δφ*}とを主連立微分方程式及び双対連立微分方程式に乗じて積分して得られる式(例えば、上記式(400)(401)又は、上記式(402)を解くことにより、固有値と固有ベクトルを得ることができる。これらを用いて、主問題の変位{u}が求められる。

【0157】
以上のように、有限要素法を用いることで、形状が複雑で、固有関数を直接計算するのが困難な場合であっても、近似の固有関数を求めることができる。

【0158】
[第5の実施例]
図67に示す情報処理装置10cは、初期方程式決定部21、境界条件決定部22、及び演算部23を備える。初期方程式決定部21は、処理対象となる系の構造及び構成要素の性質を示すデータを読み込んで、系を表す方程式を決定する。この方程式は、例えば、求めるべき物理量を表す変数を含むn個の方程式とすることができる。方程式は、例えば、行列のデータ又は微分方程式のデータで表すことができ、データの表現形式は特に限定されない。変数は、例えば、ベクトル又は関数等のデータで表すことができるが、表現形式に特に限定されない。

【0159】
初期方程式決定部21は、例えば、ユーザからの処理対象の系の指定を受け付けて、予め記録しておいた連立微分方程式のうちから指定された系に応じたものを選択することができる。あるいは、連立微分方程式自体をユーザから受け付ける構成であってもよい。

【0160】
境界条件決定部22は、境界条件を示すデータとして、物理量を現す値を読み込んで境界条件を決定する。境界条件決定部22は、例えば、1つの点の1方向、すなわち1自由度に対して外力と変位を指定するような境界条件も受け付ける。後述する演算部23の機能により、1自由度に対して複数の物理量を境界条件として定義しても、対応可能となっている。また、境界条件として、既知部だけでなく、未知部の指定を受け付けることもできる。例えば、ある1つの節点に対して、変位も外力も未知とする設定を受け付けることができる。

【0161】
演算部23は、前記n個の初期方程式を、2n個の変数の方程式又は2n個の方程式を含む方程式に変形し、変形した2n個の変数の方程式又は2n個の方程式を含む方程式において、境界条件によって既知となる変数を含む既知部と、未知の変数を含む未知部とを決定し、未知部について方程式の解を計算する。例えば、初期方程式として連立微分方程式が決定された場合、式(963)で示すように、初期方程式の2倍の数の方程式が含まれる形式に変形することができる。初期方程式として行列と変数ベクトルで表される方程式が決定された場合、例えば、式(442)で示すように、変数ベクトルの2倍の自由度を持った変数ベクトルを含む行列に変形することができる。

【0162】
演算部23は、計算した解又は解から得られる情報を出力することができる。例えば、処理対象の系において物体が力を受けて変位する場合、解から得られる情報として、変形状態、応力分布、それらの誤差分布等を出力することができる。

【0163】
(随伴微分作用素を用いた演算例)
初期方程式決定部21は、例えば、前記系の構成要素である物体に作用する力と、前記物体の物理量を表す変数との関係を表す連立微分方程式を決定することができる。物体の物理量は、例えば、変位又は速度等とすることができるが、特定のものに限定されない。

【0164】
演算部23は、初期方程式決定部21が決定した連立微分方程式の微分作用素と、当該微分作用素から決まる随伴微分作用素を用いて、前記2n個の方程式を示すデータを生成し、2n個の方程式の解を計算する。演算部23は、対象となる系を表す主問題の連立微分方程式の微分作用素及び変数に対して、双対変数を含む双対連立方程式(n個の方程式を含むもの)を新たに導入する。これにより、2n個の方程式の解を計算することになる。このように、微分作用素から決まる随伴微分作用素を用いて、前記2n個の方程式を生成し、2n個の方程式の解を計算することで、1自由度に対して複数の物理量が境界条件として設定された場合にも解を計算することができる。

【0165】
なお、例えば、式(64)に示すように、連立微分方程式の微分作用素を変数の作用させたものと双対変数との内積は、変数と双対連立微分方程式の微分作用素を双対変数に作用させたものとの内積に等しい関係にある。本願明細書では、様々な連立微分方程式について、双対微微分方程式を求める例を開示している。なお、予め、双対微分方程式に関する情報を、微分方程式と対応付けて記録しておくことができる。演算部23は、初期方程式決定部21で決定した微分方程式に対応する双対微分方程式に関する情報を読み出して演算に用いることができる。

【0166】
(有限要素を用いた演算例)
また、演算部23は、随伴作用素を計算せずに、初期方程式を2n個の変数を含む方程式に変形し、方程式の解を計算することができる。例えば、初期方程式決定部21が、2つのn次元の変数ベクトルとn行の行列を決定することで、系の構成要素の節点における物理量を示す初期方程式を決定することができる。例えば、節点の変位を表す変数のn次元ベクトルと、節点の力を表す変数のn次元ベクトルとの関係を表すn行の行列を生成することができる。これは、系が満たすべき微分方程式を、n次元ベクトルとn行の行列を用いて表したものと言える。例えば、式(441)や式(449)等に示すような形で初期方程式を決定することができる。式(441)の例では、n行n列の行列Kと、n次元ベクトルの変数ベクトル{U},{F}でn個の方程式が表される。
境界条件決定部22は、方程式の変数において、値が既知となる変数の自由度の数と、値が未知となる変数の自由度の数が異なる境界条件も決定可能に構成される。例えば、ユーザから、境界条件として、1自由度に対して複数の物理量が未知または既知であると指定された場合も、ユーザの設定が不適切である旨等のエラーは出さないよう構成することができる。

【0167】
この場合、演算部23は、2つの変数ベクトルに基づいて2n次元のベクトルを生成する。さらに、この2n次元のベクトルの変数に基づいて行列を2n列の行列に変形した上で、当該2n次元のベクトルの変数のうち、境界条件によって既知となる変数を含む既知部分と、未知の変数を含む未知部分とを決定する。この時、未知部分の自由度と、既知部分の自由度は同じでなくてもよい。例えば、演算部23は、未知部分の自由度と、既知部分の自由度を同じにするという制限を加えずに、未知部分と既知部分を決定することができる。演算部23は、未知部分の変数を、既知部分の変数で表す形式になるよう前記2n列の行列及び前記2n次元のベクトルを変形して、変形した行列を用いて、前記未知部分の変数を計算する。例えば、演算部23は、変数の既知部と未知部をそれぞれまとめるように変数ベクトルの要素を並べ替え、これに応じて行列の要素も並べ替えることができる。なお、この具体例は、9.2節の式(441)~式(447)に示されている。
上記演算部23は、上記処理を実行することにより自己随伴な問題と、非自己随伴な問題の両方について解を計算することができる。例えば、1つの自由度、例えば、ある節点のx方向については変位も力もゼロであるという境界条件が設定された場合や、ある節点のy方向については変位も力も未知であるといった境界条件も受け入れることができる。上記のような、変数をn倍にした方程式で、未知部と既知部の自由度の数に制限をかけない演算により、非自己随伴境界条件の下でも解を計算することが可能になる。

【0168】
解が複数存在する場合、例えば、式(128)の汎関数Πにおいて、変分がゼロとなるように、同次解に対するモード係数を決定し、決定したモード係数を用いて解を求めることができる。また、決定したモード係数近傍の値の入力をユーザから受け付け、入力された値のモード係数を用いて解を計算し、前記解または前記解から得られる情報を出力する構成とすることもできる。なお、これらの構成は、演算部23が設けられた情報処理装置に設ける必要は必ずしもなく、演算部23で計算された複数の解について、上記モード係数を用いて解を求める処理を、他の情報処理装置が実行してもよい。また、ユーザからモード係数近傍の値の入力をユーザから受け付け、入力された値のモード係数を用いて解を計算して出力する処理を他の情報処理装置が実行してもよい。

【0169】
例えば、演算部23が計算した方程式の解として、1つの特解と、複数の同次解が得られた場合、解は、特解と同次解に任意のモード係数をかけたもので表すことができる。このモード係数を用いて表される解を、式(128)の汎関数Πにおいて、変分がゼロとする式に当てはめたときに得られるモード係数により、構造全体系として近似的に微分方程式を満足する解が得られる。このような解が複数存在する場合の処理の例は、例えば、9.4節、9.5節及び11.3.10節に記載されている。
本実施形態では、上記演算部23では、自由度を2倍にした方程式を、未知部と既知部の自由度が必ずしも同じでない状況で方程式の解を計算する。この演算部23では、複数の解が算出される場合がある。その場合にも、上記例のように、複数解のうち代表的な解を計算し、この解を含む所定の範囲で解を出力することで適切な結果を出力することができる。また、ユーザの操作により解を変動させて表示することで、複数解の結果を認識可能な形態で伝えることができるユーザインタフェースが提供される。

【0170】
[第6の実施形態]
図68に示す情報処理装置10dは、初期方程式決定部21、境界条件決定部22、及び判定部33を備える。初期方程式決定部21及び境界条件決定部22は、図67の初期方程式決定部21及び境界条件決定部22と同様に構成することができる。判定部33は、初期方程式決定部21で決定された連立微分方程式の境界条件が、自己随伴条件か否かを判定する。

【0171】
(随伴境界条件を用いて判定する例)
判定部33は、例えば、連立微分方程式の変数と同じ数の双対変数と、双対連立微分方程式とを定義した場合に、連立微分方程式と双対変数との積分和すなわち内積を部分積分して得られる境界項を計算し、境界条件決定部22で決定された境界条件の下で、当該境界項をゼロにする双対変数の条件を、随伴境界条件として決定することができる。判定部33は、随伴境界条件が境界条件と一致する場合は、自己随伴条件であると判断できる。

【0172】
随伴境界条件を算出するため、例えば、予め、微分方程式に対応付けて境界項を計算するための式を情報処理装置10dがアクセス可能な記録部に記録しておくことができる。判定部33は、初期方程式決定部21が決定した微分方程式に対応する境界項の式を記録部から読み出して、境界項を計算することができる。

【0173】
具体例として、四角形要素モデルにおいて、式(250)又は式(1814)に示したように、境界項Rを、節点外力F及び変位Uの転置行列と,双対節点外力F及び双対変位Uの行列との積で表すことができる。このような場合に、式(251)又は式(1814)のように、節点外力F、変位U、双対節点外力F及び双対変位Uそれぞれを未知部と既知部に分けて並べ替えることで、随伴境界条件が得られる。このように、微分方程式の変数の転置行列と双対変数の行列との積で表される境界項において、それぞれの変数の未知部と既知部を分けて並べ替えることで、随伴境界条件を得ることができる。式(251)又は式(1814)は、例えば、初期微分方程式の作用素及び境界条件を、有限要素に設定して解析を実行する場合に用いることができる。

【0174】
判定部は、例えば、自己随伴条件であるか否かを判定結果に基づく情報をユーザに対して出力することができる。例えば、境界条件と随伴境界条件が一致せず、非自己随伴条件であると判定した場合、解析対象が境界条件化において非自己随伴である旨のメッセージを出力することができる。例えば、判定部33は、解析には非自己随伴問題を扱える機能が必要であること、又は、非自己随伴問題を解析できるソフトウエアの入手先または入手方法など、非自己随伴問題の解析のガイドとなる情報を出力することができる。

【0175】
或いは、判定部33は、随伴境界条件と境界条件が一致するか否かの情報に基づいて、その後の演算を制御することもできる。例えば、随伴境界条件と境界条件が一致する場合は、従来の有限要素法を用いた解析処理を実行させ、随伴境界条件と境界条件が一致しない場合は、上記のように、上記したような、非自己随伴を扱うことが可能な処理を情報処理装置に実行させることができる。非自己随伴を扱うことができる処理は、例えば、初期方程式を2n個の方程式を含む方程式に変形し、変形した2n個の方程式を含むにおいて、既知の変数を含む既知部と、未知の変数を含む未知部とを決定し、当該未知部について方程式の解を計算するといった処理である。

【0176】
(随伴境界条件を計算せずに自己随伴条件を判定する例)
なお、上記例のように、随伴境界条件を計算しなくても、自己随伴条件の判定ができる場合がある。例えば、11.7.1に記載されているように、上記の式(251)において、自己随伴境界条件となる場合、節点変位Uの既知部Ubの個数n1は、節点外力の未知部Fvに等しく、節点外力Fの既知部Fbの個数n2は、節点変位の未知部Uvの個数に等しくなる。したがって、判定部33は、1つの変数の未知部の自由度が、他の変数の既知部の自由度と等しいか否かによって、自己随伴条件か否かを判定することができる。

【0177】
(判断部と演算部との組み合わせ例)
なお、上記図68の判定部33を図68の情報処理装置に追加した構成も可能である。例えば、判定部33の判定結果に基づいて、演算部23が、2n次元の方程式を用いた演算を行うか否かを制御することができる。これにより、演算部23は、2n次元の方程式を用いた非自己随伴の演算と、n次元の方程式を用いた自己随伴の演算とを適切に切り替えて、処理効率を向上させることができる。

【0178】
具体例として、境界条件決定部22は、ユーザから前記変数のうち既知のものと未知のものの指定を受け付けることができる。判定部33は、前記ユーザから受け付けた指定に基づいて前記変数の未知部と既知部を決定して、1つの変数の未知部の自由度が、他の変数の既知部の自由度と等しいか否かを出力する。判断部33において、1つの変数の未知部の自由度が、他の変数の既知部の自由度と等しいと判断された場合、演算部23は、前記微分作用素を表す前記n行の行列を用いて、前記変数の未知部を計算する。判断部33において、1つの変数の未知部の自由度が、他の変数の既知部の自由度と等しくないと判断された場合、前記微分方程式を、変数が2n個で2n個の方程式を含む形式に変形し、変形した2n個の方程式のうち、前記未知部について方程式の解を計算する。

【0179】
[動作例]
(自己随伴条件に対応する設定例)
図69は、従来の有限要素を用いた解析において、境界条件が自己随伴条件を満たすことを前提として、節点変位ベクトル、節点外力ベクトルを設定する処理例を示すフローチャートである。図69に示す例では、情報処理装置は、まず、節点変位の入力を受け付ける(S101)。例えば、変位がわかっている自由度とその数値がユーザから入力される。あるいは、節点変位の自由度と数値を示すデータが読み込まれてもよい。自由度は、例えば、節点の位置と方向により特定される。節点変位が入力されると、節点変位の既知部Ubの自由度と、その数値が確定する(S102)。節点変位の既知部Ubの自由度及び数値が確定すると自動的に節点外力の未知部Fvの自由度を確定することができる(S103)。自己随伴条件を満たす場合、未知部Fvの自由度は、既知部Ubに等しいからである。

【0180】
節点変位の既知部Ubの自由度が決まると、全体自由度から既知部Ubの自由度を引くことで、節点変位の未知部Uvの自由度が求められる(S104)。自己随伴条件を満たす場合は、既知部Fbの自由度は、未知部Uvの自由度に等しいので、既知部Fbの自由度が確定し、既知部Fbの数値はゼロと仮定される(S105)。S107において、ユーザから節点外力の自由度と数値が入力されると、情報処理装置は、節点外力の既知部Fbが自由度に反しないか確認した上で既知部Fbの値を修正する。

【0181】
このように、図69に示すフローでは、節点変位の既知部Ubが決まると未知部Fv,未知部Uv、及び既知部Fbが自動的に計算される。値が指定されていない節点外力の自由度にはゼロが設定される。このように、自由度が自己随伴境界条件を満たすので、設定された節点変位及び節点外力に基づいて、従来型有限要素法で解析可能である(S108)。

【0182】
(非自己随伴条件にも対応した境界条件の設定例1)
図70は、有限要素を用いた解析において、境界条件が非自己随伴条件である場合にも対応可能な節点変位ベクトル、節点外力ベクトルを設定する処理例を示すフローチャートである。図70に示す例では、情報処理装置の境界条件決定部22は、節点変位の自由度及び数値の入力を受け付け(S201)、節点変位の既知部Ubの自由度と数値を決定し記録する(S202)。既知部Ubが決まると節点外力の未知部Fv(S203)、節点変位の未知部Uv(S204)、及び既知部Fb(S205)が自動的に計算される。S205では、値が指定されていない節点外力の自由度にはゼロが設定される。

【0183】
境界条件決定部22は、節点外力の自由度と数値の入力を受け付け(S206)、節点外力Fbの自由度及び数値を入力された自由度及び数値に修正する(S207)。さらに、境界条件決定部22は、節点変位の未知部の入力も受け付ける(S208)。本例では、変位が未知である節点の自由度が入力される。例えば、ユーザが画面上で節点及び方向を指定することで入力できる構成でもよいし、未知部の自由度を示すデータを読み込む構成でもよい。境界条件決定部22は、未知部Uvの自由度を、入力された情報を元に修正する(S209)。

【0184】
同様に、境界条件決定部22は、節点外力の未知部の入力も受け付け(S210)、未知部Fvの自由度を、入力された情報を元に修正する(S211)。S208~S211において、同じ自由度に対して変位と外力の両方を未知部にする指定も受け入れることができる。すなわち、境界条件決定部22は、同じ自由度に対して異なる変数を未知とする条件を設定することができる。

【0185】
判定部33は、自由度が自己随伴境界条件を満たすか否かを判定し(S212)、判定結果を表示する(S213)。判定部33は、例えば、節点変位の既知部Ubの自由度の数と節点外力の未知部Fvの数が等しいか否かを判定することができる。等しい場合は、自己随伴境界条件を満たし、等しくない場合は、自己随伴境界条件を満たさないと判断できる。

【0186】
例えば、ある1つの節点の1つの方向、すなわち1つの自由度に対して、複数の変数(本例では、節点変位が節点外力)を既知とする境界条件や、1つの自由度に対して、複数の変数の未知とする境界条件は、自己随伴境界条件を満たさないと判断される。自己随伴境界条件を満たさない境界条件が設定された場合、従来は、解析を実行することができなかったが、本願明細書に開示の技術(例えば、新型有限要素法を用いた技術)により、非自己随伴境界条件の解析が可能になる(S215)。自己随伴境界条件の解析は、従来型有限要素法で解析可能である(S214))。

【0187】
なお、S214又はS215の解析処理は、別の装置で実行されてもよい。例えば、判定部33は、S212で自己随伴境界条件を満たすと判断された場合は、S214の従来型有限要素法での解析処理を実行し、自己随伴境界条件を満たさないと判断された場合は、S212において、S215の解析処理を可能にするため情報を出力することができる。例えば、S215の処理を実行するためのソフトウエアの入手方法や案内等の情報を出力することができる。

【0188】
(非自己随伴条件にも対応した境界条件の設定例2)
図71は、境界条件が非自己随伴条件である場合にも対応可能な節点変位ベクトル、節点外力ベクトルを設定する処理の他の例を示すフローチャートである。図71に示す例では、節点変位および節点外力の自由度と数値の入力を受け付け、入力がなかった自由度は未知部として設定する。境界条件決定部22は、節点変位の自由度と数値の入力を受け付け(S301)、入力に基づいて、節点変位の既知部Ubの自由度と数値を確定し(S302)、これに伴って未知部Uvの自由度も確定する(S303)。同様にして、節点外力の自由と数値の入力を受け付け(S304)、入力に基づいて、節点外力の既知部Fbの自由度と数値を確定し(S305)、これに伴って未知部Fvの自由度も確定する(S306)。判定処理(S307)、表示処理(S308)、解析処理(S309、310)は、図70の場合と同様に実行することができる。

【0189】
図71に示す例では、特に入力がない節点変位及び節点外力の自由度は、未知部として設定される。ユーザは、未知部を明示的に指定しなくてもよいが、例えば、変位をゼロにしたい節点については、明示的にゼロの値を設定することを入力する必要がある。これに対して図70に示す例では、ユーザが未知部の指定も既知部の指定もしなかった自由度の変数は、既知部としてなんらかの値が設定される。そのため、ユーザがすべての自由度に対して未知か既知かの設定をする必要がなくなる。なお、情報処理装置では、例えば、ユーザが図70又は図71のいずれの処理で入力するかを選択できる構成とすることができる。なお、上記の図70、及び図71の処理例は、11.7.1に記載の2通りの処理方法(1)、(2)に相当する。

【0190】
(有限要素で非自己随伴問題の解析をする例)
図72は、有限要素を用いた、非自己随伴問題の解析の実行例を示すフローチャートである。図72に示す例では、まず、情報処理装置10cは、計算モデルを入力する(S401)。例えば、処理対象となる構造の有限要素モデル、形状、材料、要素、節点等が入力される。

【0191】
初期方程式決定部21は、S401で入力されたデータを用いて、要素剛性行列を作成し、要素剛性行列を重ね合わせて全体剛性行列を作成する(S405)。ここで作成される全体剛性行列は、行列で表現される初期方程式の一例である。S401では、例えば、自由度の数がn個である系の節点における物理量を示すn次元の変数ベクトルを2つと、n行の要素剛性行列が決定される。具体的には、式(441)や式(449)の行列ようにすることができる。

【0192】
境界条件決定部22は、重力荷重のデータを入力し(S403)、重力荷重を作成する(S405)。S403では、例えば、節点と方向で特定される自由度それぞれに対して加わる重力の値又は変位の値を示すデータが入力される。また、境界条件決定部22は、節点変位の未知の部分と、節点外力が未知の部分を指定するデータも入力する(S404)。S403及びS404で入力されたデータによって、節点変位の既知部と未知部、及び節点外力の既知部と未知部が作成される(S405)。なお、S403、S404及びS405の処理は、例えば、図70で示した処理と同様に実行することができる。

【0193】
演算部23は、剛性要素行列を変形して演算に用いる行列を作成する(S406)。例えば、自由度nの2つの変数ベクトルに基づいて、演算部23は、2n次元のベクトルを生成し、当該2n次元のベクトルの変数に応じて前記n行の行列を2n列の行列に変形する。例えば、式(450)のような形に変形できる。演算部23は、当該2n次元のベクトルの変数のうち、境界条件によって既知となる変数を含む既知部分と、未知の変数を含む未知部分とを決定する。未知部分の変数を、前記既知部分の変数で表す形式になるよう前記2n列の行列及び前記2n次元のベクトルを変形する。例えば、式(452)のように変形した行列を用いて、未知部分の変数を計算することができる。

【0194】
演算部23は、S406で作成した行列の方程式の解を計算する(S407)。計算の結果、解が1つだけ存在する場合、演算部23は、節点変位の特解部分及び節点外力の特解部分を作成する(S408)。例えば、式(936)のように、特解の未知部と既知部を表すベクトルを作成することができる。複数解が存在する場合、演算部23は、同次モード係数を定義する(S409)。例えば、同次解の個数をnとすると、n次元ベクトル{a}をモード変数として定義することができる。演算部23は、節点変位の同次解部分と、節点外力の同次解部分を作成する(S410)。例えば、式(937)に示すように、同次解の未知部と既知部を表すベクトルが作成される。

【0195】
S411において、演算部23は、節点変位の解の未知部を、例えば、式(938)を用いて決定し、節点外力の解の未知部を、例えば、式(939)を用いて決定する。S411の処理の結果、解が一意に決定できる場合は、その結果を出力して解析処理を終了する(S412)。

【0196】
解が一意に決定できない場合、最小二乗法により、{a}の最適値を計算する。例えば、演算部23は、式(960)、式(961)のような行列を作成する(S413)。さらに、演算部23は、例えば、式(962)を、モード係数{a}について解くことで、モード係数を決定する(S414)。演算部23は、決定したモード係数を中心値とする所定の範囲のモード係数に対応する解を出力する(S415)。例えば、決定したモード係数近傍の値の入力をユーザから受け付け、入力された値のモード係数を用いて解を計算し、この計算した解またはこの解から得られる情報を出力する。解から得られる情報として、例えば、変形状態、応力分布、誤差分布等が出力される。また、ユーザが、スライダーにより{a}を操作して、S414で決定されたモード係数を含む所定範囲の値を指定できるようにユーザインタフェースを設けることができる。

【0197】
なお、S408~S410の、解が複数のある場合の処理を実行するためのプログラムと、S401~S407の処理を実行するためのプログラムは、それぞれ、別に提供されてよい。すなわち、S408~S410の少なくとも1つのステップを実行するプログラム又はプログラムを記録した記録媒体も、本願発明の実施形態に含まれる。


【0198】
[解析結果例]
図9は、一様重力下の円環で,内辺の境界条件を「変位がゼロかつ表面力もゼロ」とするには,外辺を如何に支持すべきかを求めた結果を示す図である。従来は求解不能とされてきた問題であるが,新たに解析解を得て,これと固有関数法による計算結果とを比較し,解析解に近い事を確認した。円環の上部と下部に大きな垂直応力σrが作用し,左右部
には剪断応力τrθが作用する事で,重力荷重を受け持っている様子が判る。内辺は真円
に保たれ,表面力も作用していない事が確認できる。なお、図9において、左図は応力σr、中央図は応力σθ、右図は剪断応力τrθを示す。図9に示す結果は、円環を有限要素とせず、解析的に計算した結果である。以下に、解析に用いた具体的な主連立微分方程式の一例を示す。

【0199】
【数270】
JP0005886309B2_000273t.gif

【0200】
図10は、一様重力下の正方形で、上辺の境界条件を「変位がゼロかつ表面力もゼロ」とした場合の解析結果を示す図である。すなわち、一様重力下で、正方形の上辺の変位と表面力の両方をゼロとするためには、他辺をどのように支持すればよいかを求める問題の結果である。図10においてそれぞれ、左図は垂直応力σ、中央図は垂直応力σ、右図は剪断応力τxy、の分布を示す。図10の左図より、正方形の下部に大きな応力σとσが作用することによって上辺の変位がゼロに保たれていることがわかる。ここでは、正方形を有限要素として定式化しており、要素はひとつだけである。ひとつの要素で解析できる方法としてp法有限要素がある。図10に示す解析は、非自己随伴な境界条件の例であり、現代のp法有限要素法ではこのような境界条件での解析はできない。また、境界条件を辺に与えるか、節点に与えるかは定式化の都合によるが、ここでは、節点に境界条件を与えるよう定式化している。

【0201】
図11は、周辺を固着した正方形の平板の面内変形の解祈結果を示す図である。図11においてそれぞれ、左図は垂直応力σ、中央図は垂直応力σ、右図は剪断応力τxy、の分布を示す。図11の中央図より、正方形の上部と下部に大きな応力σが作用することがわかる。これは、自己随伴な境界条件である。図11に示す解析結果は、p法有限要素と同等の解析を行ったものとなる。上記発明の解析方法は、p法を包含するものと考えられる。すなわち、本願発明を基に有限要素法をコード化しておけば、シームレスにp法も使うことができる。

【0202】
図12は、正方形の面外変形平板要素の周辺を固着した場合の自重による変形を解析した結果を示す図である。正方形を有限要素として定式化しており、要素はひとつだけである。

【0203】
図13Aは、左端自由固着の梁のモデルを図であり、図13Bは、図13Aに示す梁の主固有関数φ1、φ2、φ3、及び双対固有関数φ1*、φ2*、φ3*の計算結果を示す。図1
3Cは、固有関数法により求めた解(1モード、2モード及び3モード)および解析解(Analytical Solution)を示す図である。3モードの解は、解析解とほぼ重なっている。

【0204】
図14及び図15は、無限領域に正方形の仮想境界を定義し,その中心に湧き出し(図14)或いは渦度(図15)を配置した場合の解析結果を示す図である。図14及び図15に示す解析結果は、上記(12)式或いは(13)式の右辺にディラックのδ関数を与え,固有関数法で計算したポテンシャル流れの様子を示している。これらの図から、内部の流れの様子は解析解に近い事が確認できる。

【0205】
[解析例]
10.1 ばね減衰系
10.1.1 微分方程式
図16のばね減衰系にて、減衰係数をc、ばね定数をkとする。変位をu( t)、外力をf(t)とすると、作用素Lを
【数271】
JP0005886309B2_000274t.gif

【0206】
10.1.2 随伴境界条件と随伴微分作用素
【数272】
JP0005886309B2_000275t.gif

【0207】
【数273】
JP0005886309B2_000276t.gif

【0208】
10.1.3 境界条件の同次化
【数274】
JP0005886309B2_000277t.gif

【0209】
10.1.4 連立固有値問題と固有関数法
【数275】
JP0005886309B2_000278t.gif
【数276】
JP0005886309B2_000279t.gif

【0210】
10.1.5 固有関数のセット
【数277】
JP0005886309B2_000280t.gif
【数278】
JP0005886309B2_000281t.gif

【0211】
10.1.6 固有関数
【数279】
JP0005886309B2_000282t.gif
【数280】
JP0005886309B2_000283t.gif

【0212】
【数281】
JP0005886309B2_000284t.gif
となる。それぞれの固有関数を、最大値が1になるように正規化して図17に示す。

【0213】
10.1.7 解析例1
【数282】
JP0005886309B2_000285t.gif
として,(525)式の解析解と,(494)式の固有関数法で得た結果の比較を図18に示す.図中,破線が解析解を,実線が固有関数法を示している.50Mは50モードまで計算したことを示している.モード数が多いほど,固有関数法の精度は良くなり,10モード程度で変位u^は良く一致していることがわかる。

【0214】
10.1.8 解析例2
【数283】
JP0005886309B2_000286t.gif
として,(528)式の解析解と、(494)式の固有関数法で得た結果の比較を図19に示す。図中、破線が解析解を、実線が固有関数法を示している。50Mは50モードまで計算したことを示している。モード数が多いほど、固有関数法の精度は良くなり、10モード程度で変位u^は良く一致していることがわかる。

【0215】
10.1.9 解析例3
【数284】
JP0005886309B2_000287t.gif
として、(531)式の解析解と、(494)式の固有関数法で得た結果の比較を図20に示す。図中、破線が解析解を、実線が固有関数法を示している。50Mは50モードまで計算したことを示している。モード数が多いほど,固有関数法の精度は良くなり、10モード程度で変位u^は良く一致していることがわかる.

【0216】
10.2 質点ばね系
10.2.1 微分方程式
図21の質点ばね系にて、質点の質量をm、ばね定数をkとする。変位をu(t)、外力をf(t)とすると、作用素Lを
【数285】
JP0005886309B2_000288t.gif
【数286】
JP0005886309B2_000289t.gif

【0217】
10.2.2 随伴境界条件と随伴微分作用素
【数287】
JP0005886309B2_000290t.gif

【0218】
【数288】
JP0005886309B2_000291t.gif

【0219】
【数289】
JP0005886309B2_000292t.gif

【0220】
10.2.3 境界条件の同次化

【数290】
JP0005886309B2_000293t.gif

【0221】
【数291】
JP0005886309B2_000294t.gif

【0222】
10.2.4 連立固有値問題と固有関数法1
[非自己随伴境界条件の例1]の場合には、(554)式の境界関数uBを用いた、非自己随伴問題の解法となる。連立固有値問題
【数292】
JP0005886309B2_000295t.gif
【数293】
JP0005886309B2_000296t.gif
となる。

【0223】
10.2.5 連立固有値問題と固有関数法2
【数291】
JP0005886309B2_000297t.gif

【0224】
10.2.6 固有関数のセット
【数295】
JP0005886309B2_000298t.gif

【0225】
【数296】
JP0005886309B2_000299t.gif
【数297】
JP0005886309B2_000300t.gif

【0226】
これらの組み合わせは、連立微分方程式(596)を満足するだけで、境界条件の制約は課
していない。

【0227】
10.2.7 固有関数
【数298】
JP0005886309B2_000301t.gif
【数299】
JP0005886309B2_000302t.gif
【数300】
JP0005886309B2_000303t.gif
となる。(610)式でのそれぞれの固有関数を、最大値が1になるように正規化して図22に示す。

【0228】
同様に、(611)式最初の2つの固有値は、
【数301】
JP0005886309B2_000304t.gif
となる。(612)式でのそれぞれの固有関数を、最大値が1になるように正規化して図23に示す。

【0229】
10.2.8 解析例1
【数302】
JP0005886309B2_000305t.gif
として、 (625)式の解析解と、(578)式の固有関数法で得た結果の比較を図24に示す。
図中、破線が解析解を、実線が固有関数法を示している。56M は56 モードまで計算した
ことを示している。モード数が多いほど、固有関数法の精度は良くなり、10 モード程度
で変位 u^は良く一致していることがわかる。

【0230】
10.2.9 自己随伴固有関数
【数303】
JP0005886309B2_000306t.gif

【0231】
10.2.10 解析例2
【数304】
JP0005886309B2_000307t.gif
として、(632)式の解析解と、(584)式の固有関数法で得た結果の比較を図25に示す。図中、破線が解析解を、実線が固有関数法を示している。50Mは50モードまで計算したことを示している。モード数が多いほど、固有関数法の精度は良くなり、10モード程度で変位u^は良く一致していることがわかる。

【0232】
10.2.11 ハミルトニアンの原理について

【0233】
【数305】
JP0005886309B2_000308t.gif
【数306】
JP0005886309B2_000309t.gif
【数307】
JP0005886309B2_000310t.gif

【0234】
10.3 弦の静たわみ
10.3.1 微分方程式
図26の弦の静たわみにて、弦の長さをl、張力をTとする。変位をu(x)、分布荷重をf(x)とすると、作用素Lを
【数308】
JP0005886309B2_000311t.gif
【数309】
JP0005886309B2_000312t.gif
【数310】
JP0005886309B2_000313t.gif

【0235】
10.3.2 随伴境界条件と随伴微分作用素
【数311】
JP0005886309B2_000314t.gif
【数312】
JP0005886309B2_000315t.gif

【0236】
10.3.3 境界条件の同次化
非同次境界条件を満たす項に添字Bを付けてuBとおき、同次境界条件を満たす項に添
字Hを付けてuHとおく。変位uをその和で、
【数313】
JP0005886309B2_000316t.gif


【数314】
JP0005886309B2_000317t.gif

【0237】
10.3.4 固有関数のセット
【数315】
JP0005886309B2_000318t.gif
【数316】
JP0005886309B2_000319t.gif

【0238】
10.3.5 固有関数
【数317】
JP0005886309B2_000320t.gif
【数318】
JP0005886309B2_000321t.gif
となる。それぞれの固有関数を、最大値が1になるように正規化して図27に示す。

【0239】
【数319】
JP0005886309B2_000322t.gif

【0240】
10.3.6 固有関数法
【数320】
JP0005886309B2_000323t.gif
【数321】
JP0005886309B2_000324t.gif

【0241】
10.3.7 非自己随伴解析例1
【数322】
JP0005886309B2_000325t.gif
として、(706)式の解析解と、(702)式の固有関数法で得た結果の比較を図28に示す。図中、破線が解析解を,黒色実線が固有関数法を示している。30M は30 モードまで計算し
たことを示している。モード数が多いほど、固有関数法の精度は良くなり、3 モード程度で変位 u^は良く一致していることがわかる。

【0242】
10.3.8 非自己随伴解析例2
【数323】
JP0005886309B2_000326t.gif
として、(709)式の解析解と、(702)式の固有関数法で得た結果の比較を図29に示す。図中、破線が解析解を、黒色実線が固有関数法を示している。30M は30 モードまで計算し
たことを示している。モード数が多いほど、固有関数法の精度は良くなり、3 モード程度で変位 u^は良く一致していることがわかる。

【0243】
10.3.9 自己随伴解析例1
【数324】
JP0005886309B2_000327t.gif
として、(713)式の解析解と、(704)式の固有関数法で得た結果の比較を図30に示す。図中、破線が解析解を、黒色実線が固有関数法を示している。30M は30 モードまで計算し
たことを示している。モード数が多いほど、固有関数法の精度は良くなり、3 モード程度で変位 u^は良く一致していることがわかる。

【0244】
10.3.10 自己随伴解析例2
【数325】
JP0005886309B2_000328t.gif
として、(715)式の解析解と、(704)式の固有関数法で得た結果の比較を図31に示す。図中、破線が解析解を、黒色実線が固有関数法を示している。30M は30 モードまで計算し
たことを示している。モード数が多いほど、固有関数法の精度は良くなり、3 モード程度で変位 u^は良く一致していることがわかる。

【0245】
10.4 梁の静たわみ
10.4.1 微分方程式
図32の梁の静たわみにて,梁の長さをl、曲げ剛性をEIとする。剪断力F(x)を、曲げモーメントをM(x)とする。変位をu(x)、分布荷重を f(x)とすると、
作用素Lを
【数326】
JP0005886309B2_000329t.gif
【数327】
JP0005886309B2_000330t.gif
【数328】
JP0005886309B2_000331t.gif

【0246】
10.4.2 随伴境界条件と随伴微分作用素
【数329】
JP0005886309B2_000332t.gif
【数330】
JP0005886309B2_000333t.gif
【数331】
JP0005886309B2_000334t.gif

【0247】
10.4.3 境界条件の同次化
【数332】
JP0005886309B2_000335t.gif
【数333】
JP0005886309B2_000336t.gif
【数334】
JP0005886309B2_000337t.gif
【数335】
JP0005886309B2_000338t.gif

【0248】
10.4.4 固有関数のセット
【数336】
JP0005886309B2_000339t.gif
【数337】
JP0005886309B2_000340t.gif

【数338】
JP0005886309B2_000341t.gif

【0249】
10.4.5 固有関数
【数339】
JP0005886309B2_000342t.gif
となる。ぞれぞれの固有関数を、最大値が1となるように正規化して図33に示す。

【0250】
[自己随伴境界条件の例1]
同次境界条件が(751)式となる場合には、
【数340】
JP0005886309B2_000343t.gif

となる。それぞれの固有関数を、最大値が1になるように正規化して図34に示す。

【0251】
10.4.6 固有関数法
【数341】
JP0005886309B2_000344t.gif
となる。
【数342】
JP0005886309B2_000345t.gif

【0252】
10.4.7 非自己随伴解析例1
【数343】
JP0005886309B2_000346t.gif
として,(784)式の解析解と,(780)式の固有関数法で得た結果の比較を図35に示す。曲線の向きを図32の方向に合わせるため、各関数値には負号を乗じて作図している。縦軸の座標値のみ,正負反転して読まれたし。図中、破線が解析解を,黒色実線が固有関数法を示している。32M は32 モードまで計算したことを示している。モード数が多いほど、
固有関数法の精度は良くなり、3 モード程度で変位 u^ は良く一致していることがわか
る。

【0253】
10.4.8 非自己随伴解析例2
【数344】
JP0005886309B2_000347t.gif
として、(787)式の解析解と、(780)式の固有関数法で得た結果の比較を図36に示す。曲線の向きを図32の方向に合わせるため、各関数値には負号を乗じて作図している。縦軸の座標値のみ、正負反転して読まれたし。図中、破線が解析解を、黒色実線が固有関数法を示している。32M は32 モードまで計算したことを示している。モード数が多いほど、
固有関数法の精度は良くなり、3 モード程度で変位 u^は良く一致していることがわかる。


【0254】
10.4.9 自己随伴解析例1
【数345】
JP0005886309B2_000348t.gif
として、(790)式の解析解と、(782)式の固有関数法で得た結果の比較を図37に示す。曲線の向きを図32の方向に合わせるため,各関数値には負号を乗じて作図している。縦軸の座標値のみ、正負反転して読まれたし。図中、破線が解析解を、実線が固有関数法を示している。32Mは32モードまで計算したことを示している。モード数が多いほど、固有関数法の精度は良くなり、3モード程度で変位u^は良く一致していることがわかる。

【0255】
10.4.10 自己随伴解析例2
【数346】
JP0005886309B2_000349t.gif
として、(793)式の解析解と、(782)式の固有関数法で得た結果の比較を図38に示す。曲線の向きを図32の方向に合わせるため、各関数値には負号を乗じて作図している。縦軸の座標値のみ、正負反転して読まれたし。図中,破線が解析解を、黒色実線が固有関数法を示している。32M は32 モードまで計算したことを示している。モード数が多いほど、
固有関数法の精度は良くなり、3 モード程度で変位 u^は良く一致していることがわかる。

【0256】
[他の解析結果例]
図39は、図10に示した解析において、固有値がゼロになることにより複数の解が発生する場合の、複数の解のうちの図10に示す解以外の解の一つを示す図である。このような複数の解を設計者に提示することにより、設計者は、設計的に実現可能な方法を検討することができる。

【0257】
図40A及び図40Bは、一様重力下の正方形板の上面の変位及び表面力をゼロにする変位を、有限要素を用いて解析した場合に計算される複数の解の例である。図40A及び図40Bに示す例では、縦横6x6 個の有限要素で正方形を構成している。丸は節点を意味している。合計10個のモードが示されている。これらの任意の組み合わせが解となる。例えば、上記9.2章に示した計算により、複数解を求めることが可能になる。また、モード係数をスライダーに連動させて、設計者がスライダーによってモード係数を調節するのに合わせて、変形図が表示されるように構成すれば、設計者は、どのような支持方法があり得るのか検討でき、設計者のイメージを湧き立たせることが可能となる。また、設計者は、実現可能な支持方法を絞り込んだ後、再度計算させることも可能になる。

【0258】
[その他]
上記第1~4の実施形態における演算部1、設定部11、随伴境界条件演算部12、非自己随伴演算部13、自己随伴判定部15、自己随伴演算部14は、コンピュータのプロセッサがメモリから、プログラム及びデータを読み込んでプログラムに従ってデータの処理を実行する。演算部1は、設計データ、境界条件データを、入力デバイス及びユーザインターフェースを介して、ユーザからの入力を受け付けることにより取得してもよいし、メモリから読み込むことで取得してもよい。微分方程式データは、例えば、演算部1がアクセス可能なメモリに予め記録しておくことが好ましい。

【0259】
演算部1の機能は、1台のコンピュータのCPUで実現することもできるし、互いにデータ通信が可能な複数のコンピュータで実現することもできる。また、上記演算部1の処理をコンピュータに実行させるプログラムやそれを記録した非一時的(non-transitrory
)な記録媒体も本願発明の実施形態に含まれる。

【0260】
なお、本願発明は、上記第1~4の実施形態に限定されない。本願発明は、例えば、構造解析、振動解析、電磁場解析、材料解析、流体解析、温度解析、音場解析、電磁気学解析、回路シミュレータなど、微分方程式を用いたあらゆる解析に適用することができる。

【0261】
11.補足
11.1 節5.2の「変分直接法と最小2乗法」における式変形の詳細
11.1.1 主問題の変分直接法
同次化された微分方程式(40)式
【数347】
JP0005886309B2_000350t.gif
【数348】
JP0005886309B2_000351t.gif
【数349】
JP0005886309B2_000352t.gif
となる.これは,(111)式を複数組み合わせたものに他ならない.様々なδuに対して,即ち様々なδcに対して上式を作成し,連立方程式を代数的に解く事で,係数cを定めることができる.
以上のことから,(795)式を変形した式
【数350】
JP0005886309B2_000353t.gif
【数351】
JP0005886309B2_000354t.gif
となる.これが,主問題の変分直接法である.

【0262】
また,自己随伴問題の場合には(107)式を通じて(主,双対)固有関数φ、φ は一致するから,(120),(121)式も同等になり,(122)式は,
【数352】
JP0005886309B2_000355t.gif
となる。

【0263】
11.1.2 主問題の最小2乗法
主連立微分方程式(105)
【数353】
JP0005886309B2_000356t.gif
を得る.成分ごとに同じ重みwを用いた(主,双対)固有関数では,
【数354】
JP0005886309B2_000357t.gif
とおいて,その変分がゼロになる事と等価である.即ち,最小2乗法が変分原理となる事を示している.

【0264】
11.1.3 双対問題の変分直接法
同次化された微分方程式(66)式
【数355】
JP0005886309B2_000358t.gif
となる.これにより,係数cが定まる.(104)式の正規化によって,
【数356】
JP0005886309B2_000359t.gif
となる.これは,(115)式を複数組み合わせたものに他ならない.様々なδuに対して,即ち様々なδcに対して上式を作成し,連立方程式を代数的に解く事で,係数cを定めることができる.

【0265】
以上のことから,(802)式を変形した式
【数357】
JP0005886309B2_000360t.gif
となる.これが,双対問題の変分直接法である.

【0266】
また,自己随伴問題の場合には(107)式を通じて(主,双対)固有関数φ,φは一致するから,(120),(121)式も同等になり,(808)式は,
【数358】
JP0005886309B2_000361t.gif
となる.

【0267】
11.1.4 双対問題の最小2乗法
双対連立微分方程式(106)
【数359】
JP0005886309B2_000362t.gif
と定義できるから,上式と(120)式より,主変分δuは,
【数360】
JP0005886309B2_000363t.gif
とおいて,その変分がゼロになる事と等価である.即ち,最小2乗法が変分原理となる事を示している.


【0268】
11.2 章6の「弾性学のエネルギ原理について」における式変形の詳細
11.2.1 部分積分公式
2次元領域に直交座標x,yをとり,物体の境界表面をC,内部領域をSとする.境界表面の外向き単位法線ベクトルの成分をn,nとする.任意の関数f,fをx,yの関数として部分積分すると,
【数361】
JP0005886309B2_000364t.gif
となる.

【0269】
11.2.2 自己随伴問題の変分直接法
自己随伴問題の変分直接法は,(123)式で表され,
【数362】
JP0005886309B2_000365t.gif
を得る過程を補足説明する.

【0270】
弾性学の微分方程式は,(1),(2)式で表され,
【数363】
JP0005886309B2_000366t.gif
となる.これらの式を用いて,(123)式の左辺は,
【数364】
JP0005886309B2_000367t.gif
である.これら応力成分を用いて,(819)式を表現すれば,
【数365】
JP0005886309B2_000368t.gif
となる.部分積分公式(817)より,
【数366】
JP0005886309B2_000369t.gif

【数367】
JP0005886309B2_000370t.gif
【数368】
JP0005886309B2_000371t.gif
となる.(123)式の変分直接法は,これがゼロとなることを意味しているから,上式をゼロとおいて,
【数369】
JP0005886309B2_000372t.gif
を得る.(132),(829)式は共に変分直接法に他ならないが,従来は,仮想仕事の原理と呼ばれていた形式である.

【0271】
11.2.3 ひずみエネルギの変分
ひずみエネルギUは(130)式で表され,
【数370】
JP0005886309B2_000373t.gif
となる.上式に含まれる応力の変分を計算するために,(7),(8),(9)式の応力ひずみ関係,
【数371】
JP0005886309B2_000374t.gif
となる.整理して,
【数372】
JP0005886309B2_000375t.gif
を得る.

【0272】
11.2.4 エネルギの保存則
エネルギ保存則を表す(134)式のクラペイロンの定理は,
【数373】
JP0005886309B2_000376t.gif
となる.これは,7章(4)で述べたように,(829)式から変分記号 を除去したものに一致する.

【0273】
11.2.5 h法要素の概略
有限要素法の一般的教科書では,(829)式を仮想仕事の原理と認識して要素の剛性行列が構成されるように記されているが,7章(4)で述べたように,実際のところは,エネルギ保存則を表す(838)式が使われている.次章にて,詳しく述べる.


【0274】
11.3 章9の「従来有限要素の利用」における式変形の詳細
11.3.1 形状関数
8.1節の図8では,16節点四角形要素の節点配置を示した.h法にこれを用いることもできるが,より代表的な要素である4節点四角形要素の節点配置を図41に示す.その形状関数N(i=1~4)は,
【数374】
JP0005886309B2_000377t.gif
の不等号の範囲で境界内部を表し,等号にて境界を表す.

【0275】
11.3.2 要素内座標
形状関数Ni(i=1~n)を用いて要素内部の座標xを,
【数375】
JP0005886309B2_000378t.gif
と表す.

【0276】
11.3.3 要素内変位
8.3節の(163)式では,変位uを行列で,
【数376】
JP0005886309B2_000379t.gif
と表した.ここで,uAjは境界関数,uOjは補正関数である.境界関数uAjは要素境界のみならず,要素内の変位をも表現するが,微分方程式を満足する保証はないので,補正関数uOjを付け加えていた.しかし,従来要素即ちh法の要素では,要素内変位を境界関数uAjのみで表すところに特徴がある.即ち,補正関数uOjを付け加えずに,
【数378】
JP0005886309B2_000380t.gif
とかく. Uijは節点変位である.ここでは,アイソパラメトリック要素として,内挿関数χAiを形状関数Nと一致させる.変位の詳細は,
【数379】
JP0005886309B2_000381t.gif
と比較するとわかりやすい.

【0277】
11.3.4 ひずみと応力
ひずみ成分は,(4)~(6)式から
【数380】
JP0005886309B2_000382t.gif

【0278】
変位の微分係数は,
【数381】
JP0005886309B2_000383t.gif
【数382】
JP0005886309B2_000384t.gif
となる.

【0279】
11.3.5 ひずみエネルギ
ひずみエネルギUは(130)式で表され,
【数383】
JP0005886309B2_000385t.gif
である.これは,エネルギ保存則(837)式の左辺に使われており,2倍して(838)式の左辺を表す.

【0280】
要素内部のひずみエネルギengは,
【数384】
JP0005886309B2_000386t.gif
にて得られる.上式では,エネルギのUと節点変位のUが同じ記号となってしまっているのでまぎらわしいが,左辺のUがひずみエネルギを表し,右辺に2つ現れる{U}が節点変位を表している.上式の積分項を
【数385】
JP0005886309B2_000387t.gif
にて計算できることを示している.

【0281】
11.3.6 等価節点力
外力のなす仕事wrkは,
【数386】
JP0005886309B2_000388t.gif
であり,エネルギ保存則(837)式の右辺に使われており,2倍して(838)式の右辺を表す.表面力の項p,pを境界上で積分して節点に振り分けたものをF,物体力の項b,bを要素内で積分して節点に振り分けたものをF,として,単位厚さあたりの等価節点力 を,
【数387】
JP0005886309B2_000389t.gif
となる.

【0282】
11.3.7 h法の要素剛性行列
エネルギ保存則(837)式に従い,(854)式と(859)式を等置して,
【数388】
JP0005886309B2_000390t.gif
と認識することも有用である.単位厚さあたりの要素剛性行列[Kelm]を
【数389】
JP0005886309B2_000391t.gif
となり,節点変位Uと等価節点力Fとのエネルギ的な等価性を表している.

【0283】
解析対象となる物体を膨大な数の要素に分割した後,個々の要素に対して(865)式を作成し,共有節点に注意しながら全体系の方程式を構成すれば,
【数390】
JP0005886309B2_000392t.gif
となる.ここで,Nは全体系の自由度の総数である.総自由度Nのことをnと表したものが(434)式である.これを力の釣り合いと認識して解く方法がh法による線形弾性有限要素法,即ち従来法である.

【0284】
(865)式を得る手順を熟察すれば,変分演算はどこにも用いられておらず,仮想変位も現れないことがわかる.多くの教科書では(829)式を計算して(865)式が得られるかの如く説明されるが,実際のところはそうではなく,エネルギ保存則(837),(838)式が使われているのである.故に,(866),(434)式のことを,力の釣り合いと理解するよりは,エネルギ的な等価性を表していると理解する方が理にかなっている.

【0285】
従って,9.2節に示したように,自己随伴な境界条件に縛られる必要はなく,非自己随伴な境界条件に対して(866),(434)式を解いても,良いのである.

【0286】
11.3.8 有限要素の特徴
新しい有限要素では,単位厚さあたりの等価節点力Fは,(214)式で,
【数391】
JP0005886309B2_000393t.gif
である.それぞれ,要素境界上で表面力を積分することで得られたものである.
これに対して,従来の有限要素法即ちh法の要素では,
単位厚さあたりの等価節点力Fは,(861)式で,
【数392】
JP0005886309B2_000394t.gif
である.これは,要素内部のひずみエネルギと等価節点力による仕事とを等置して得られたものである.

【0287】
新しい有限要素での剛性行列は[K],[K]の2つから成るのに対し,h法要素での剛性行列は[Keng]の1つから成る.その詳細,即ち(212),(213)式および(853)式の計算法からみて,新旧2つの剛性行列は,全く性質が異なるものであり,補正関数による効果即ち[K]を無視すれば両者が一致するといったものではないことがわかる.

【0288】
11.3.9 長方形要素の剛性行列
図41の節点配置で,節点1,2,3,4の座標をそれぞれ(-a,-b),(a,-b),(a,b),(-a,b)として,アスペクト比κが,

【数393】
JP0005886309B2_000395t.gif
である.その結果,図41の4節点要素での要素内座標xは,

【数394】
JP0005886309B2_000396t.gif

【0289】
要素内変位uは(845)式で表され,
【数395】
JP0005886309B2_000397t.gif
である.アイソパラメトリック要素として,内挿関数χAiを形状関数Nと一致させる.Uijは節点変位で,その並びに応じて各種行列の並びが定まり,(192),(193)式より,
【数396】
JP0005886309B2_000398t.gif
である.よって,
【数397】
JP0005886309B2_000399t.gif
となる.

【0290】
以上により,(870)式が計算できて,
【数398】
JP0005886309B2_000400t.gif
となる.ここで,
【数339】
JP0005886309B2_000401t.gif
とおいて,
【数400】
JP0005886309B2_000402t.gif
である.また,
【数401】
JP0005886309B2_000403t.gif
である.単位厚さあたりの要素剛性行列[Kelm]は,(864)式に(876)式を代入して,平面応力状態,平面ひずみ状態,共に,
【数402】
JP0005886309B2_000404t.gif
と定義するのが一般的である.上式を利用する場合には,平面応力状態と平面ひずみ状態とで,行列要素を書き換える必要があることから,少し面倒である.しかし,教科書で良く見る式を示しておいて,それと(893)式との対応を知る事も必要なので,以下に両者の行列要素を記しておく.結果として,両者を区別する必要のない(893)式の方が扱い易いことがわかるだろう.

【0291】
平面応力状態のとき,ヤング率[Young’s modulus]をEとして,(894)式の係数cは,
【数403】
JP0005886309B2_000405t.gif
とおいて,
【数404】
JP0005886309B2_000406t.gif
である.また,
【数405】
JP0005886309B2_000407t.gif
を代入した要素剛性行列[Kelm]は平面応力状態の(894)式に一致する.

【0292】
同様に,平面ひずみ状態のとき,ヤング率をEとして, (894)式の係数cは,
【数406】
JP0005886309B2_000408t.gif
とおいて、
【数407】
JP0005886309B2_000409t.gif
である.また,
【数408】
JP0005886309B2_000410t.gif
を代入した要素剛性行列[Kelm]は平面ひずみ状態の(894)式に一致する.

【0293】
11.3.10 解の複数存在する劣決定系に関する補足
9.4節では,同次解の個数 は,
【数409】
JP0005886309B2_000411t.gif
である.この場合には,(453)式を(933)式に読み替えれば良い.

【0294】
11.3.11 最小自乗法に関する補足
9.5節で示した変分原理の計算方法を補足説明する.まずは,(251)式と同様に,節点変位Uを既知部Uと未知部Uに分けて,
【数410】
JP0005886309B2_000412t.gif
と定義する.これは,8.11節で述べた内容と同一であるが,ここではnを構造全体をメッシュに分割してh法要素で表したきの総自由度数と捉えているところが異なる.
同次解に対するモード係数を{a}とおいて,解{s}は(455)式で,
【数411】
JP0005886309B2_000413t.gif
および,
【数412】
JP0005886309B2_000414t.gif
となる.未知部UとFはいずれも,モード係数{a}によって表され,問題は完全に解かれたことになる.即ち,任意の{a}に対して解となる.

【0295】
しかし,この解は11.3.7節で述べたように,エネルギ的な等価性を表しているだけで,要素内部でも構造全体でも微分方程式を満足する保証は無い.そこで,(122)式式の変分原理および最小自乗法を用いることで,構造全体系としては近似的に微分方程式を満足する解が得られるよう,変分原理により{a}を求めた後,(938)式で未知変位,(939)式で未知外力を決定する.変分原理により定まった{a}の近傍では,近似的に微分方程式を満足する解が表現できていると考えることができるから,設計者は,スライダーにより を操作して,変形状態や応力分布さらには誤差分布を確認しながら,設計作業を進めることができる.

【0296】
節点変位Uは,(934),(938)式より,
【数413】
JP0005886309B2_000415t.gif
と定義すれば,
【数414】
JP0005886309B2_000416t.gif
および,
【数415】
JP0005886309B2_000417t.gif
である.ここで, nはひとつの要素の自由度数である.同じ記号nが使われているとややこしいので,次のように記号を置き換える.構造全体をメッシュに分割してh法要素で表したきの総自由度数をnとすると,(943)式は,
【数416】
JP0005886309B2_000418t.gif
となる.上式は,要素内の変位分布を表しており,モード係数{a}の関数となっていることがわかる.

【0297】
ある要素の連立偏微分方程式(23)は,
【数417】
JP0005886309B2_000419t.gif
となる.しかし,上式の等号が満たされる保証はなく,適切な{a}を与えることで,近似的に等号が成立することを期待することしかできない.そこで,k番目の要素に対する誤差関数Eを,
【数418】
JP0005886309B2_000420t.gif
とおいて,その変分がゼロになるようにモード係数{a}を決めれば,変分原理を適用したことになる.ただし,この段階では,ひとつの要素について適用したことにしかならない.構造全体を表す要素の総数をn,構造全体での汎関数をΠとすると,
【数419】
JP0005886309B2_000421t.gif
となる.そして, Πの変分がゼロになるようにモード係数{a}を決めれば,構造全体に対して変分原理を適用したことになり,変分直接法を表す(127)式,
【数420】
JP0005886309B2_000422t.gif
である.よって,変分δΠは,
【数421】
JP0005886309B2_000423t.gif
となる.任意の変分{δa}に対して(957)式が成立するためには,
【数422】
JP0005886309B2_000424t.gif
が成立すれば良い.上式は,モード係数{a}について解ける.これにより,構造全体系として近似的に微分方程式を満足する解が得られたと考えることができ,求まった{a}を(938)式に戻して未知変位,(939)式に戻して未知外力が定まる.この値を起点として,設計者は,スライダーにより{a}を操作して,変形状態や応力分布さらには誤差分布を確認しながら,設計作業を進めることができる.

【0298】
11.3.12 構造的特徴の共通点に関する補足
連立偏微分方程式(23)は,
【数423】
JP0005886309B2_000425t.gif
を用いた.2次元弾性体の静的釣合方程式の場合にも当てはまるが,より一般の2次元問題に対して,(105),(106)式をまとめて記すと,
【数424】
JP0005886309B2_000426t.gif
となる.上式で,左辺の行列部分は自己随伴作用素となる.即ち,上式には,非自己随伴な問題を自己随伴化させ,完全系の固有関数を得るための工夫が織り込まれている.2次元問題ゆえに方程式(23)も2つ存在するが,未知関数(変数とも言う)はφ,φ2,φ,φ であって,方程式数の2倍の4つが定義されている.また,8.18節のように,固有値 がゼロとなる関数φ,φ2,φ,φ も用いて解析がなされる.一意の解にも複数の解にも対応できる解法となっており,これが構造的特徴である.

【0299】
一方,従来要素を用いた定式化,(442)式,
【数425】
JP0005886309B2_000427t.gif
のように解が求まる.一意の解にも複数の解にも対応できる解法となっており,これが構造的特徴である.
これらの事項は,両解法に一致した共通の特徴である.

【0300】
11.4 円環の静たわみ
本節は,10章に含める方が良い.
11.4.1 微分方程式
2次元弾性体の静的釣合方程式は,直交座標系で,
【数426】
JP0005886309B2_000428t.gif
となる.ここで,u,uはx,y方向の変位で,b,bはx,y方向の単位体積あたりの物体力,Gは剛性率である.ポアソン比をνとして,定数μは,
【数427】
JP0005886309B2_000429t.gif
と定義し,定数と認識する.

【0301】
11.4.2 座標変換
直交座標系と極座標系との関係は,
【数428】
JP0005886309B2_000430t.gif
である。

【0302】
11.4.3 変位変換
変位の直交座標成分u,uと極座標成分u,uθとの関係は,
【数429】
JP0005886309B2_000431t.gif
となる.

【0303】
11.4.4 表面力変換
表面力の変換は変位と同様である.表面力の直交座標成分p,pと極座標成分p,pθとの関係は,
【数430】
JP0005886309B2_000432t.gif
となる.

【0304】
11.4.5 応力変換
応力の直交座標成分σ,σ,τxyと極座標成分σ,σθ,τrθとの関係は,
【数431】
JP0005886309B2_000433t.gif
【数432】
JP0005886309B2_000434t.gif
となる.

【0305】
11.4.6 微分作用素の変換
(967)式を微分して,
【数433】
JP0005886309B2_000435t.gif
となる.もう一度微分して,
【数434】
JP0005886309B2_000436t.gif
となる.(989)式と(990)式を加えて,
【数435】
JP0005886309B2_000437t.gif
となる.これは,ラプラス作用素[Laplacian]▽の変数変換を表している.

【0306】
11.4.7 重調和方程式
微分方程式(1),(2)の連立をはずして,重調和方程式
【数436】
JP0005886309B2_000438t.gif
となる.作用素を展開すると,
【数437】
JP0005886309B2_000439t.gif
である.
極座標系の微分方程式(996),(997)の一般解をψ(r,θ)とおいて,
【数438】
JP0005886309B2_000440t.gif
とおく.

【0307】
11.4.8 変数分離解
(1000)式を(999)式に代入して,
【数439】
JP0005886309B2_000441t.gif
となる.上記の一般解は,
【数440】
JP0005886309B2_000442t.gif
が,変数分離解となる.
なお,上式ではm=0の場合に不都合が生じるので,
【数441】
JP0005886309B2_000443t.gif
が,変数分離解となる.
さらに,(1006)式ではm=1の場合にも不都合が生じるので,
【数442】
JP0005886309B2_000444t.gif
を得るので,これらを(1001)式に代入して,整理すると,
【数443】
JP0005886309B2_000445t.gif
が,変数分離解となる.
以上のことから,(1006),(1011),(1016)式を微分方程式(1001)の解の形式とすれば良い.

【0308】
11.4.9 境界条件
2次元円環の外径をR ,内径をγR とする.内径比率γの範囲は,
【数444】
JP0005886309B2_000446t.gif
である.一様重力下の円環で,内辺の境界条件を,表面力がゼロ,変位もゼロとするには,外辺をどのように支持すれば良いか求める.従来は求解不能とされた問題である.
内辺の変位をゼロとする境界条件は,
【数445】
JP0005886309B2_000447t.gif
となる.

【0309】
11.4.10 解析解
微分方程式(1),(2)と境界条件式(1018),(1019)を満たす解では,(1016)式の形式は不要となるので,微分方程式(996),(997)の一般解u,u を,それぞれ(1006),(1011)式で表す.
計算の結果,関数を,
【数446】
JP0005886309B2_000448t.gif
にある.ここで, b,bはx,y方向の単位体積あたりの物体力,Gは剛性率である. f,fは長さの逆数の次元を持つから,代表長さを外径Rとして,無次元の荷重係数c,cを,
【数447】
JP0005886309B2_000449t.gif
と定義すると,問題を扱いやすい.

【0310】
平面応力状態で,材料定数を,
【数448】
JP0005886309B2_000450t.gif
として,一様重力荷重fを作用させた場合の,変形と応力分布の様子を図42に示す.ここで,荷重係数c,cを,
【数449】
JP0005886309B2_000451t.gif
とおいて,変位を外径Rで無次元化している.

【0311】
11.4.11 部分積分公式
11.2.1節で示したように,任意の関数f,fをx,yの関数として部分積分すると,
【数450】
JP0005886309B2_000452t.gif
が得られる.また,
【数451】
JP0005886309B2_000453t.gif
が得られる.

【0312】
11.4.12 随伴微分作用素と随伴境界条件
微分方程式(1),(2)の右辺を(964)式で表して,
【数452】
JP0005886309B2_000454t.gif
となる.ここで,原初微分作用素Lijは,(15)式で表され,
【数453】
JP0005886309B2_000455t.gif
となる.

【0313】
双対変位uに関する準備として,(4),(5),(6)式に相当する変位ひずみ関係は,平面応力,平面ひずみ状態共に,
【数454】
JP0005886309B2_000456t.gif
となる.(10),(11),(12)式に相当する応力成分は,
【数455】
JP0005886309B2_000457t.gif
となる.

【0314】
では,(58)式の左辺を部分積分しよう.
・|i=1,j=1| のとき,(1029),(1030)式を利用して,
【数456】
JP0005886309B2_000458t.gif
となる.上式右辺の境界項をR11,微分作用素をL11と書き,
【数457】
JP0005886309B2_000459t.gif
とする.L11は(16)式に示した随伴微分作用素である.
・|i=1,j=2| のとき,(1031)式を利用して,
【数458】
JP0005886309B2_000460t.gif
とする.L12は(16)式に示した随伴微分作用素である.
・ |i=2,j=1|のとき,(1032)式を利用して,
【数459】
JP0005886309B2_000461t.gif
とする.L21は(16)式に示した随伴微分作用素である.
・ |i=2,j=2|のとき,(1029)式,(1030)式を利用して,
【数460】
JP0005886309B2_000462t.gif
となる.上式右辺の境界項をR22 ,微分作用素をL22と書き,
【数461】
JP0005886309B2_000463t.gif
となる場合を自己随伴な微分作用素と呼ぶ.

【0315】
境界項を加え合わせて, Rと表し,
【数462】
JP0005886309B2_000464t.gif
となる.(10),(11),(12)式,及び,(1043),(1044),(1045)式を利用して,
【数463】
JP0005886309B2_000465t.gif
となる.

【0316】
なお,(1059)式の境界項Rを極座標成分で表すために,(971),(972)式を利用して,
【数464】
JP0005886309B2_000466t.gif
となる.上式が,極座標成分で書いた境界項である.
境界項Rをゼロとする条件から,随伴境界条件が定まる.
以下にその例を示すが,円環の外辺の半径をr=R ,内辺の半径をr=γRとしており,境界項Rと同じ記号が現れるのでまぎらわしいが,両者の意味合いが大きく異なるので判別は可能である.

【0317】
[非自己随伴境界条件の例1]
円環にて,内辺(r=γR )で表面力がゼロかつ変位もゼロとする境界条件は,(1018),(1019)式であり,
【数465】
JP0005886309B2_000467t.gif
となる.これは,外辺(r=R)で双対表面力p,pがゼロかつ双対変位u,uもゼロとする条件であり,上式が随伴境界条件となる.(1069)式と(1070)式では条件が異なるので,非自己随伴境界条件である.

【0318】
なお,条件式(1018),(1019)は,(1069)式に一致する.即ち,11.4.9節で与えた条件は,非自己随伴境界条件である.

【0319】
[自己随伴境界条件の例1]
円環にて,内辺(r=γR)で変位がゼロ,外辺(r=R )でも変位がゼロとする境界条件は,
【数466】
JP0005886309B2_000468t.gif
を得る.(1071)式と(1072)式では条件が一致するので,自己随伴境界条件である.

【0320】
11.4.13 境界条件の同次化と境界項
非同次境界条件を満たす項に添字Bを付けてuBjとおき,同次境界条件を満たす項に添字Hを付けてuHjとおく.主変位uをその和で,
【数467】
JP0005886309B2_000469t.gif
を得る.これと,同次随伴境界条件を満たす関数uHiとの内積は,
【数468】
JP0005886309B2_000470t.gif
となる.

【0321】
では,(61)式の左辺を部分積分しよう.
・|i=1,j=1|のとき,(1029),(1030)式を利用して,
【数469】
JP0005886309B2_000471t.gif
とする.L11 は(16)式に示した随伴微分作用素である.
・ |i=1,j=2|のとき,(1031)式を利用して,
【数470】
JP0005886309B2_000472t.gif
とする. L12 は(16)式に示した随伴微分作用素である.
・ |i=2,j=1|のとき,(1032)式を利用して,
【数471】
JP0005886309B2_000473t.gif
とする. L21は(16)式に示した随伴微分作用素である.
・ |i=2,j=2|のとき,(1029),(1030)式を利用して,
【数472】
JP0005886309B2_000474t.gif
とする. L22は(16)式に示した随伴微分作用素である.
この例のように,
【数473】
JP0005886309B2_000475t.gif
となる.上式を変形して,境界積分の関係式(1034)を利用すると,
【数474】
JP0005886309B2_000476t.gif
となる.(27),(30),(33)式,及び,(45),(48),(51)式を利用して,
【数475】
JP0005886309B2_000477t.gif
である.極座標成分で表すために,(1068)式を得たのと同様の変形をして,
【数476】
JP0005886309B2_000478t.gif
となる.

【0322】
11.4.14 連立固有値問題の連立解除
代表長さを外径Rとして,重み定数を
【数477】
JP0005886309B2_000479t.gif
となる.
主連立固有値問題は,(98)式
【数478】
JP0005886309B2_000480t.gif
となる.上式の連立をはずして整理すると,φ,φは共に微分方程式,
【数479】
JP0005886309B2_000481t.gif
【数480】
JP0005886309B2_000482t.gif
を満たす関数φとなる.

【0323】
11.4.15 固有関数のセット
(1094)式と(1098)式は同じ微分方程式であるから,結局,φ,φとφ,φは,同じ関数セットから構成されることがわかる.よって,(1094)式を解く事に照準を合わせて,
【数481】
JP0005886309B2_000483t.gif
となる.
さて,θの関数には2πの周期性を期待してmを整数とおき,変数分離解を求めることにして,微分方程式
【数482】
JP0005886309B2_000484t.gif
である.ここで,Im’は第1,2種変形ベッセル関数[modified Bessel function of the first(second)kind]である.
同様にして,微分方程式
【数483】
JP0005886309B2_000485t.gif
である.

【0324】
11.4.16 主連立固有値問題を満たす解の組み合わせ
主連立固有値問題の(1092),(1093)式を変形して,
【数484】
JP0005886309B2_000486t.gif
となる.解φ,φを,(1105)式あるいは(1107)式を満たす関数系から選んだ場合には,φ,φは共に,
【数485】
JP0005886309B2_000487t.gif
を満たさねばならない.これは,括弧内の項が定数になることを要求している.しかし, φ,φの解形式が(1106),(1108)であることから,ゼロ以外の定数になることはないので, φ,φの組み合わせが満たすべき方程式として,
【数486】
JP0005886309B2_000488t.gif
となる.解φ,φを,(1109)式あるいは(1111)式を満たす関数系から選んだ場合には, φ,φは共に,
【数487】
JP0005886309B2_000489t.gif
を満たさねばならない.これは,括弧内の項が定数になることを要求している.しかし, φ,φの解形式が(1110), (1112)であることから,ゼロ以外の定数になることはないので, φ,φの組み合わせが満たすべき方程式として,
【数488】
JP0005886309B2_000490t.gif
が得られる.

【0325】
11.4.17 双対連立固有値問題を満たす解の組み合わせ
双対連立固有値問題の(1096),(1097)式を変形して,
【数489】
JP0005886309B2_000491t.gif
となる.(1099)式のωを用いて,
【数490】
JP0005886309B2_000492t.gif
を満たさねばならない.これは、括弧内の項が定数になることを要求している.しかし、φx*y*の解形式が(1106),(1108)であることから,ゼロ以外の定数になることはないので,φxy*の組み合わせが満たすべき方程式として,
【数491】
JP0005886309B2_000493t.gif
が得られる.
一方,双対連立固有値問題の(1096),(1097)式は,
【数492】
JP0005886309B2_000494t.gif
となる.解 φ,φを,(1109)式あるいは(1111)式を満たす関数系から選んだ場合には, φ,φは共に,
【数493】
JP0005886309B2_000495t.gif
を満たさねばならない.これは,括弧内の項が定数になることを要求している.しかし, φ,φの解形式が(1110), (1112)であることから,ゼロ以外の定数になることはないので, φ,φの組み合わせが満たすべき方程式として,
【数494】
JP0005886309B2_000496t.gif
が得られる.

【0326】
11.4.18 主連立微分方程式を満たす解の組み合わせ
主連立微分方程式の(1087),(1088)式を変形して,
【数495】
JP0005886309B2_000497t.gif
となる.
解 φ,φを,(1105)式
【数496】
JP0005886309B2_000498t.gif
を得る.これを主連立微分方程式(1087),(1088)に代入して,
【数497】
JP0005886309B2_000499t.gif
となる.これは,(1106)式の関数系から成り,(1120)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している.
【数498】
JP0005886309B2_000500t.gif
【数499】
JP0005886309B2_000501t.gif
となる.これは,(1108)式の関数系から成り,(1120)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している.
一方,主連立微分方程式の(1087),(1088)式は,
【数500】
JP0005886309B2_000502t.gif
となる.
【数501】
JP0005886309B2_000503t.gif
【数502】
JP0005886309B2_000504t.gif
となる.これは,(1110)式の関数系から成り,(1128)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している.
【数503】
JP0005886309B2_000505t.gif
を満たさねばならないから,(1163),(1164)式は,
【数504】
JP0005886309B2_000506t.gif
となる.これは,(1112)式の関数系から成り,(1128)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している.

【0327】
11.4.19 双対連立微分方程式を満たす解の組み合わせ
双対連立微分方程式の(1089),(1090)式を変形して,
【数505】
JP0005886309B2_000507t.gif
【数506】
JP0005886309B2_000508t.gif


を満たさねばならないから,(1179),(1180)は,
【数507】
JP0005886309B2_000509t.gif
となる.これは,(1106)式の関数系から成り,(1136)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している.
解φ,φを,(1107)式
【数508】
JP0005886309B2_000510t.gif
を得る.これを双対連立微分方程式(1089),(1090)に代入して,
【数509】
JP0005886309B2_000511t.gif
となる.これは,(1108)式の関数系から成り,(1136)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している.
一方,双対連立微分方程式の(1089),(1090)式は,
【数510】
JP0005886309B2_000512t.gif
なる変形もできる.(1100)式のω2を用いて,
【数511】
JP0005886309B2_000513t.gif
【数512】
JP0005886309B2_000514t.gif
を得る.これを双対連立微分方程式(1089),(1090)に代入して,
【数513】
JP0005886309B2_000515t.gif
となる.これは,(1110)式の関数系から成り,(1144)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している.
解φ,φを,(1111)式
【数514】
JP0005886309B2_000516t.gif
を満たさねばならないから,(1195),(1196)式は,
【数515】
JP0005886309B2_000517t.gif
となる.これは,(1112)式の関数系から成り,(1144)式を満足するφ,φの組み合わせは,上式の連立微分方程式を満足することを示している。

【0328】
11.4.20 解関数の様子
(1106),(1108),(1110),(1112)式の解を用いるにあたって,整数mを偶数と奇数の場合に分けるのが便利で、nを整数として、
【数516】
JP0005886309B2_000518t.gif
が成立することから,上下対称(x軸に対称),左右対称(y軸に対称)な関数となる.即ち、x方向に偶関数(even)で、y方向に偶関数(even)となる.頭文字をとって,添字eeを付ける。

【0329】
(2) 関数φoeを、
【数517】
JP0005886309B2_000519t.gif
が成立することから,上下対称(x軸に対称),左右反対称(y軸に反対称)な関数となる.即ち、x方向に奇関数(odd)で,y方向に偶関数(even)となる. 頭文字をとって,添字oeを付ける。

【0330】
(3) 関数φooを,
【数518】
JP0005886309B2_000520t.gif
が成立することから,上下反対称(x軸に反対称),左右反対称(y軸に反対称)な関数となる.即ち,x方向に奇関数(odd)で,y方向に奇関数(odd)となる.頭文字をとって,添字ooを付ける。

【0331】
(4) 関数φeoを,
【数519】
JP0005886309B2_000521t.gif
が成立することから,上下反対称(x軸に反対称),左右対称(y軸に対称)な関数となる.即ち, x方向に偶関数(even)で、y方向に奇関数(odd)となる.頭文字をとって,添字eoを付ける。

【0332】
さて今度は,x方向の変位をφ,y方向の変位をφとする.上記関数の組み合わせによって,以下の代表的な4つの変形モードが生じる。

【0333】
[SA]
x方向の変位が(1210)式の形の関数φxeeで,y方向の変位が(1214)式の形の関数φyooで表されるとき,変形は図43の様にx軸に対称(symmetry),y軸に非対称(asymmetry)になる.頭文字をとって,このような変形をモードSAと呼ぶ.この変形を表す関数の組み合わせφSAを、
【数520】
JP0005886309B2_000522t.gif
のように,添字SAを付けて定義する。

【0334】
[AS]
x方向の変位が(1214)式の形の関数φxooで、y方向の変位が(1210)式の形の関数φyeeで表されるとき,変形は図の44の様にx軸に非対称(asymmetry),y軸に対称(symmetry)になる.頭文字をとって,このような変形をモードASと呼ぶ.この変形を表す関数の組み合わせφASを,
【数521】
JP0005886309B2_000523t.gif
のように,添字ASを付けて定義する。

【0335】
[SS]
x方向の変位が(1212)式の形の関数φxoeで,y方向の変位が(1216)式の形の関数φyeoで表されるとき,変形は図45の様にx軸に対称(symmetry),y軸に対称(symmetry)になる.頭文字をとって,このような変形をモードSSと呼ぶ.この変形を表す関数の組み合わせφssを,
【数522】
JP0005886309B2_000524t.gif
のように,添字SSを付けて定義する。

【0336】
[AA]
x方向の変位が(1216)式の形の関数φxeoで,y方向の変位が(1212)式の形の関数φyoeで表されるとき,変形は図46の様にx軸に非対称(asymmetry),y軸に非対称(asymmetry)になる.頭文字をとって,このような変形をモードAAと呼ぶ.この変形を表す関数の組み合わせφAAを,
【数523】
JP0005886309B2_000525t.gif
のように,添字AAを付けて定義する.

【0337】
11.4.21 第1種ベッセル関数を用いた解関数
固有値ω1、ωを総称してωと表し,(1106),(1110)式の第1種ベッセル関数Jを用いた解を,(1209)式に従い,
【数524】
JP0005886309B2_000526t.gif
と表す.これらの関数の主変数はr、θであり,n、ωは助変数である.これら変数を全て記すと式が長くなるので,以下では,注目する変数だけを記す.例えば,等式の両辺でnは変化しても,r、θ、ωは変化しない場合には,引数を(n)と表現する.このとき,次のように導関数
【数525】
JP0005886309B2_000527t.gif
が得られる.さらに,(1226)~(1229)式からは,関係式
【数526】
JP0005886309B2_000528t.gif
【数527】
JP0005886309B2_000529t.gif
が得られる。

【0338】
ここで,(1234)と(1239)式,(1235)と(1238)式,(1237)と(1240)式,(1236)と(1241)式の組合せに注目して,
【数528】
JP0005886309B2_000530t.gif
の4組を定義すると,それぞれ,(1120)式
【数529】
JP0005886309B2_000531t.gif
が成立する.よって,(1242)~(1245)式でωを(1099)式のωに差し替えた組合せは,主連立固有値問題(1092),(1093)式を満足することがわかる.即ち、φ、φは主固有関数である.さらに,(1151),(1152)式から、φ、φの符号を反転させたものは双対固有関数φ*、φ*となる.以上により,主連立微分方程式(1087),(1088)を満たす関数セットφ、φ、φ*、φ* が4組得られたことになる。

【0339】
また,(74)~(76)式を用いて,(1242)~(1245)式のそれぞれに応力が計算できて,
【数530】
JP0005886309B2_000532t.gif
となる.よって、ωを(1099)式のω1に差し替えた組合せは,主固有関数に基づく応力σx、σy、τxyを表す。なお,(82)~(84)式及び,(1151),(1152)式から、σx、σy、τxyの符号を反転させれば,双対固有関数に基づく応力σ*σ*τxy*を表す。以上により,応力を表す関数セットσx、σy、τxy、σ*σ*τxy*が4組得られたことになる。

【0340】
改めて定義すると、
【数531】
JP0005886309B2_000533t.gif
の4組である。まだ境界条件は反映されていないが,任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0341】
同様に,(1238)と(1235)式,(1239)と(1234)式,(1241)と(1236)式,(1240)と(1237)式の組合せに注目して,
【数532】
JP0005886309B2_000534t.gif
が成立する。よって,(1254)~(1257)式でωを(1100)式のω2に差し替えた組合せは、主連立固有値問題(1092),(1093)式を満足することがわかる。即ち、φ、φは主固有関数である。さらに、(1167),(1168)式から、φ、φの符号を反転させたものは双対固有関数φ*、φ*となる。以上により、主連立微分方程式(1087),(1088)を満たす関数セットφ、φy、φ*、φ*が4組得られたことになる。

【0342】
また、(74)~(76)式を用いて、(1254)~(1257)式のそれぞれに応力が計算できて、
【数533】
JP0005886309B2_000535t.gif
となる。よって、ωを(1100)式のωに差し替えた組合せは、主固有関数に基づく応力σx、σy、τxyを表す。なお、(82)~(84)式及び,(1167),(1168)式から、σx、σy、τxyの符号を反転させれば,双対固有関数に基づく応力σ*σ*τxy*を表す。以上により、応力を表す関数セットσx、σy、τxy、σ*σ*τxy*が4組得られたことになる。

【0343】
改めて定義すると、
【数534】
JP0005886309B2_000536t.gif
の4組である。まだ境界条件は反映されていないが、任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0344】
11.4.22 第2種ベッセル関数を用いた解関数
固有値ω、ωを総称してωと表し、(1106),(1110)式の第2種ベッセル関数Yを用いた解を、(1209)式に従い、
【数535】
JP0005886309B2_000537t.gif
【数536】
JP0005886309B2_000538t.gif
と表す。これらの関数の主変数はr、θであり、n、ωは助変数である。これら変数を全て記すと式が長くなるので、以下では、注目する変数だけを記す。例えば、等式の両辺でnは変化しても、r、θ、ωは変化しない場合には,引数を(n)と表現する。このとき、次のように導関数
【数537】
JP0005886309B2_000539t.gif
が得られる.さらに,(1270)~(1273)式からは,関係式
【数538】
JP0005886309B2_000540t.gif
が得られ,(1274)~(1277)式からは,関係式
【数539】
JP0005886309B2_000541t.gif
が得られる。

【0345】
ここで,(1278)と(1283)式,(1279)と(1282)式,(1281)と(1284)式,(1280)と(1285)式の組合せに注目して、
【数540】
JP0005886309B2_000542t.gif
が成立する。よって,(1286)~(1289)式でωを(1099)式のωに差し替えた組合せは、主連立固有値問題(1092),(1093)式を満足することがわかる。即ち、φ、φは主固有関数である。さらに、(1151),(1152)式から、φ、φの符号を反転させたものは双対固有関数φ*、φ*となる。以上により、主連立微分方程式(1087),(1088)を満たす関数セットφ、φ、φ*、φ*が4組得られたことになる。

【0346】
また,(74)~(76)式を用いて、(1286)~(1289)式のそれぞれに応力が計算できて、
【数541】
JP0005886309B2_000543t.gif
となる。よって、ωを(1099)式のωに差し替えた組合せは、主固有関数に基づく応力σx、σy、τxyを表す。なお、(82)~(84)式及び、(1151),(1152)式から、σx、σy、τxyの符号を反転させれば、双対固有関数に基づく応力σ*σ*τxy*を表す。以上により、応力を表す関数セットσx、σy、τxy、σ*σ*τxy*が4組得られたことになる。

【0347】
改めて定義すると、
【数542】
JP0005886309B2_000544t.gif
【数543】
JP0005886309B2_000545t.gif
の4組である。まだ境界条件は反映されていないが,任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0348】
同様に、(1282)と(1279)式,(1283)と(1278)式, (1285)と(1280)式,(1284)と(1281)式の組合せに注目して、
【数544】
JP0005886309B2_000546t.gif
が成立する。よって、(1298)~(1301)式でωを(1100)式のωに差し替えた組合せは、主連立固有値問題(1092),(1093)式を満足することがわかる。即ち、φ、φは主固有関数である。さらに、(1167),(1168)式から、φ、φの符号を反転させたものは双対固有関数φ*、φ*となる。以上により、主連立微分方程式(1087),(1088)を満たす関数セット φ、φ、φ*、φ*が4組得られたことになる。

【0349】
また,(74)~(76)式を用いて,(1298)~(1301)式のそれぞれに応力が計算できて、
【数545】
JP0005886309B2_000547t.gif
となる。よって、ωを(1100)式のωに差し替えた組合せは、主固有関数に基づく応力σx、σy、τxyを表す。なお、(82)~(84)式及び,(1167),(1168)式から、σx、σy、τxyの符号を反転させれば,双対固有関数に基づく応力σ*σ*τxy*を表す。以上により、応力を表す関数セットσx、σy、τxy、σ*σ*τxy*が4組得られたことになる。

【0350】
改めて定義すると、
【数546】
JP0005886309B2_000548t.gif
の4組である。まだ境界条件は反映されていないが、任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0351】
11.4.23 第1種変形ベッセル関数を用いた解関数

【0352】
固有値ω、ωを総称してωと表し,(1108),(1112)式の第1種変形ベッセル関数Iを用いた解を、(1209)式に従い、
【数547】
JP0005886309B2_000549t.gif
と表す。これらの関数の主変数はr、θであり、n、ωは助変数である。これら変数を全て記すと式が長くなるので、以下では、注目する変数だけを記す。例えば、等式の両辺でnは変化しても、r、θ、ωは変化しない場合には、引数を(n)と表現する。このとき、次のように導関数
【数548】
JP0005886309B2_000550t.gif
【数549】
JP0005886309B2_000551t.gif
が得られる。

【0353】
ここで,(1322)と(1327)式,(1323)と(1326)式,(1325)と(1328)式,(1324)と(1329)式の組合せに注目して、
【数550】
JP0005886309B2_000552t.gif
が成立する。よって,(1330)~(1333)式でωを(1099)式のωに差し替えた組合せは,主連立固有値問題(1092),(1093)式を満足することがわかる。即ち、φ、φは主固有関数である。さらに,(1157),(1158)式から、φ、φと同じものが双対固有関数φ、φとなる。以上により,主連立微分方程式(1087),(1088)を満たす関数セットφ、φ、φ、φが4組得られたことになる。

【0354】
また,(74)~(76)式を用いて,(1330)~(1333)式のそれぞれに応力が計算できて、
【数551】
JP0005886309B2_000553t.gif
となる。よって、ωを(1099)式のωに差し替えた組合せは、主固有関数に基づく応力σx、σy、τxyを表す。なお、(82)~(84)式及び、(1157),(1158)式から、σx、σy、τxyと同じものが双対固有関数に基づく応力σ*、σ*、τxy*を表す。以上により,応力を表す関数セットσx、σy、τxy、σ*、σ*、τxy*が4組得られたことになる。

【0355】
改めて定義すると、
【数552】
JP0005886309B2_000554t.gif
の4組である。まだ境界条件は反映されていないが、任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0356】
同様に、(1326)と(1323)式,(1327)と(1322)式,(1329)と(1324)式,(1328)と(1325)式の組合せに注目して、
【数553】
JP0005886309B2_000555t.gif
が成立する。よって,(1342)~(1345)式でωを(1100)式のωに差し替えた組合せは、主連立固有値問題(1092),(1093)式を満足することがわかる。即ち、φ、φは主固有関数である。さらに,(1173),(1174)式から、φ、φと同じものが双対固有関数φ*、φ*となる。以上により、主連立微分方程式(1087),(1088)を満たす関数セットφ、φy、φ*、φ*が4組得られたことになる。

【0357】
また,(74)~(76)式を用いて,(1342)~(1345)式のそれぞれに応力が計算できて、
【数554】
JP0005886309B2_000556t.gif
となる。よって、ωを(1100)式のωに差し替えた組合せは、主固有関数に基づく応力σx、σy、τxyを表す。なお、(82)~(84)式及び、(1173),(1174)式から、σx、σy、τxyと同じものが双対固有関数に基づく応力σ*σ*τxy*を表す。以上により、応力を表す関数セットσx、σy、τxy、σ*σ*τxy*が4組得られたことになる。

【0358】
改めて定義すると、
【数555】
JP0005886309B2_000557t.gif
の4組である。まだ境界条件は反映されていないが、任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0359】
11.4.24 第2種変形ベッセル関数を用いた解関数

【0360】
固有値ω1、ωを総称してωと表し、(1108),(1112)式の第2種変形ベッセル関数Kを用いた解を、(1209)式に従い、
【数556】
JP0005886309B2_000558t.gif
と表す。これらの関数の主変数はr、θであり、n、ωは助変数である。これら変数を全て記すと式が長くなるので、以下では、注目する変数だけを記す。例えば、等式の両辺でnは変化しても、r、θ、ωは変化しない場合には、引数を(n)と表現する。このとき、次のように導関数
【数557】
JP0005886309B2_000559t.gif
が得られる。さらに、(1358)~(1361)式からは、関係式
【数558】
JP0005886309B2_000560t.gif
が得られる。

【0361】
ここで、(1366)と(1371)式,(1367)と(1370)式,(1369)と(1372)式,(1368)と(1373)式の組合せに注目して、
【数559】
JP0005886309B2_000561t.gif
の4組を定義すると、それぞれ、(1120)式
【数560】
JP0005886309B2_000562t.gif
が成立する。よって、(1374)~(1377)式でωを(1099)式のωに差し替えた組合せは、主連立固有値問題(1092),(1093)式を満足することがわかる。即ち、φ、φは主固有関数である。さらに、(1157),(1158)式から、φ、φと同じものが双対固有関数φ*、φ*となる。以上により、主連立微分方程式(1087),(1088)を満たす関数セットφ、φy、φ*、φ*が4組得られたことになる。

【0362】
また,(74)~(76)式を用いて,(1374)~(1377)式のそれぞれに応力が計算できて、
【数561】
JP0005886309B2_000563t.gif
となる。よって、ωを(1099)式のωに差し替えた組合せは、主固有関数に基づく応力σx、σy、τxyを表す。なお、(82)~(84)式及び、(1157),(1158)式から、σx、σy、τxyと同じものが双対固有関数に基づく応力σ*σ*τxy*を表す.以上により,応力を表す関数セットσx、σy、τxy、σ*σ*τxy*が4組得られたことになる。

【0363】
改めて定義すると、
【数562】
JP0005886309B2_000564t.gif
【数563】
JP0005886309B2_000565t.gif
の4組である。まだ境界条件は反映されていないが、任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0364】
同様に,(1370)と(1367)式,(1371)と(1366)式,(1373)と(1368)式,(1372)と(1369)式の組合せに注目して、
【数564】
JP0005886309B2_000566t.gif
が成立する。よって、(1386)~(1389)式でωを(1100)式のωに差し替えた組合せは、主連立固有値問題(1092),(1093)式を満足することがわかる。即ち、φ、φは主固有関数である。さらに、(1173),(1174)式から、φ、φと同じものが双対固有関数φ*、φ*となる。以上により、主連立微分方程式(1087),(1088)を満たす関数セットφ、φy、φ*、φ*が4組得られたことになる。

【0365】
また、(74)~(76)式を用いて、(1386)~(1389)式のそれぞれに応力が計算できて、
【数565】
JP0005886309B2_000567t.gif
となる。よって、ωを(1100)式のωに差し替えた組合せは、主固有関数に基づく応力σx、σy、τxyを表す。なお、(82)~(84)式及び、(1173),(1174)式から、σx、σy、τxyと同じものが双対固有関数に基づく応力σ*σ*τxy*を表す。以上により、応力を表す関数セットσx、σy、τxy、σ*σ*τxy*が4組得られたことになる。

【0366】
改めて定義すると、
【数566】
JP0005886309B2_000568t.gif
の4組である。まだ境界条件は反映されていないが、任意のnに対して主連立微分方程式(1087),(1088)を満たす。

【0367】
11.4.25 主固有関数と双対固有関数の関数セット

【0368】
11.4.21節から11.4.24節に示した関数セットを組み合わせたものに境界条件を課すことで,固有値とそれに属する固有関数が定まる。モードSA,AS,SS,AAの4つの変形それぞれの、主固有関数φ、φ及びこれに基づく応力σx、σy、τxy,そして,双対固有関数 φ*、φ*及びこれに基づく応力σ*σ*τxy*の具体的な数式を示す。

【0369】
[SA]
モードSAには,(1250)式のJ1SA,(1262)式のJ2SA,(1294)式のY1SA,(1306)式のY2SA,(1338)式のI1SA,(1350)式のI2SA,(1382)式のK1SA,(1394)式のK2SA,の合計8セットが属する。任意の整数nに対しても解関数となるから,総数は8nとなる。

【0370】
(1)J1SAについては、
【数567】
JP0005886309B2_000569t.gif
となる.(1398)式を(969),(970)式により変換して、
【数568】
JP0005886309B2_000570t.gif
を得る。

【0371】
(2)J2SAについては、
【数569】
JP0005886309B2_000571t.gif
となる。(1402)式を(969),(970)式により変換して、
【数570】
JP0005886309B2_000572t.gif
を得る。

【0372】
(3)Y1SAについては、
【数571】
JP0005886309B2_000573t.gif
となる.(1406)式を(969),(970)式により変換して、
【数572】
JP0005886309B2_000574t.gif
を得る.
(4)Y2SAについては
【数573】
JP0005886309B2_000575t.gif
【数574】
JP0005886309B2_000576t.gif
を得る.(1411)式を(977),(978),(979)式により変換して,
【数575】
JP0005886309B2_000577t.gif
を得る.
(5) I1SAについては,
【数576】
JP0005886309B2_000578t.gif
を得る.(1415)式を(977),(978),(979)式により変換して,
【数577】
JP0005886309B2_000579t.gif
を得る.
(6) I2SAについては,
【数578】
JP0005886309B2_000580t.gif
となる.(1418)式を(969),(970)式により変換して、
【数579】
JP0005886309B2_000581t.gif
を得る.
(7)K1SAについては,
【数580】
JP0005886309B2_000582t.gif
となる.(1422)式を(969),(970)により変換して,
【数581】
JP0005886309B2_000583t.gif
を得る.
(8)K2SAについては,
【数582】
JP0005886309B2_000584t.gif
【数583】
JP0005886309B2_000585t.gif
を得る.

【0373】
[AS]
モードASには,(1251)式のJ1AS,(1263)式のJ2AS ,(1295)式のY1AS ,(1307)式のY2AS ,(1339)式のI1AS,(1351)式のI2AS,(1383)式のK1AS,(1395)式のK2AS ,の合計8セットが属する.任意の整数nに対しても解関数となるから,総数は 8nとなる.
(1) J1ASについては,
【数584】
JP0005886309B2_000586t.gif
【数585】
JP0005886309B2_000587t.gif
を得る.(1431)式を(977),(978),(979)式により変換して,
【数586】
JP0005886309B2_000588t.gif
を得る.
(2) J2ASについては,
【数587】
JP0005886309B2_000589t.gif
となる.(1434)式を(969),(970)式により変換して,
【数588】
JP0005886309B2_000590t.gif
を得る.
(3) Y1ASについては,
【数589】
JP0005886309B2_000591t.gif
となる.(1438)式を(969),(970)式により変換して,
【数590】
JP0005886309B2_000592t.gif
を得る.
(4) Y2ASについては,
【数591】
JP0005886309B2_000593t.gif
となる.(1442)式を(969),(970)式により変換して,
【数592】
JP0005886309B2_000594t.gif
を得る.
(5) I1ASについては,
【数593】
JP0005886309B2_000595t.gif
となる.(1446)式を(969),(970)式により変換して,
【数594】
JP0005886309B2_000596t.gif
を得る.
(6) I2ASについては,
【数595】
JP0005886309B2_000597t.gif
となる.(1450)式を(969),(970)式により変換して,
【数596】
JP0005886309B2_000598t.gif
を得る.(1451)式を(977),(978),(979)式により変換して,
【数597】
JP0005886309B2_000599t.gif
を得る.
(7) K1ASについては,
【数598】
JP0005886309B2_000600t.gif
となる.(1454)式を(969),(970)式により変換して,
【数599】
JP0005886309B2_000601t.gif
を得る.
(8) K2ASについては,
【数600】
JP0005886309B2_000602t.gif
となる.(1458)式を(969),(970)式により変換して,
【数601】
JP0005886309B2_000603t.gif
を得る.

【0374】
[SS]
モードSSには,(1252)式のJ1SS,(1264)式のJ2SS ,(1296)式のY1SS,(1308)式のY2SS,(1340)式のI1SS,(1352)式のI2SS,(1384)式のK1SS,(1396)式の K2SS,の合計8セットが属する.任意の整数nに対しても解関数となるから,総数は 8nとなる.
(1) J1SSについては,
【数602】
JP0005886309B2_000604t.gif
となる.(1462)式を(969),(970)式により変換して,
【数603】
JP0005886309B2_000605t.gif
を得る.
(2) J2SSについては,
【数604】
JP0005886309B2_000606t.gif
となる.(1466)式を(969),(970)式により変換して,
【数605】
JP0005886309B2_000607t.gif
を得る.
(3) Y1SSについては,
【数606】
JP0005886309B2_000608t.gif
となる.(1470)式を(969),(970)式により変換して,
【数607】
JP0005886309B2_000609t.gif
を得る.(1471)式を(977),(978),(979)式により変換して,
【数608】
JP0005886309B2_000610t.gif
を得る.
(4) Y2SSについては,
【数609】
JP0005886309B2_000611t.gif
となる.(1474)式を(969),(970)式により変換して,
【数610】
JP0005886309B2_000612t.gif
を得る.
(5) I1SSについては,
【数611】
JP0005886309B2_000613t.gif
となる.(1478)式を(969),(970)式により変換して,
【数612】
JP0005886309B2_000614t.gif
を得る.
(6) I2SSについては,
【数613】
JP0005886309B2_000615t.gif
となる.(1482)式を(969),(970)式により変換して,
【数614】
JP0005886309B2_000616t.gif
を得る.
(7) K1SSについては,
【数615】
JP0005886309B2_000617t.gif
となる.(1486)式を(969),(970)式により変換して,
【数616】
JP0005886309B2_000618t.gif
を得る.
(8) K2SSについては,
【数617】
JP0005886309B2_000619t.gif
となる.(1490)式を(969),(970)式により変換して,
【数618】
JP0005886309B2_000620t.gif
を得る.

【0375】
[AA]
モードAAには,(1253)式のJ1AA,(1265)式のJ2AA,(1297)式のY1AA,(1309)式のY2AA,(1341)式のI1AA,(1353)式のI2AA,(1385)式のK1AA,(1397)式のK2AA,の合計8セットが属する.任意の整数nに対しても解関数となるから,総数は8nとなる.
(1) J1AAについては,
【数619】
JP0005886309B2_000621t.gif
となる.(1494)式を(969),(970)式により変換して,
【数620】
JP0005886309B2_000622t.gif
を得る.
(2) J2AAについては,
【数621】
JP0005886309B2_000623t.gif
となる.(1498)式を(969),(970)式により変換して,
【数622】
JP0005886309B2_000624t.gif
を得る.
(3) Y1AAについては,
【数623】
JP0005886309B2_000625t.gif
【数624】
JP0005886309B2_000626t.gif
を得る.
(4) Y2AAについては,
【数625】
JP0005886309B2_000627t.gif
となる.(1506)式を(969),(970)式により変換して,
【数626】
JP0005886309B2_000628t.gif
を得る.
(5) I1AAについては,
【数627】
JP0005886309B2_000629t.gif
となる.(1510)式を(969),(970)式により変換して,
【数628】
JP0005886309B2_000630t.gif
を得る.
(6) I2AAについては,
【数629】
JP0005886309B2_000631t.gif
となる.(1514)式を(969),(970)式により変換して,
【数630】
JP0005886309B2_000632t.gif
を得る.
(7) K1AAについては,
【数631】
JP0005886309B2_000633t.gif
となる.(1518)式を(969),(970)式により変換して,
【数632】
JP0005886309B2_000634t.gif
を得る.
(8) K2AAについては,
【数633】
JP0005886309B2_000635t.gif
となる.(1522)式を(969),(970)式により変換して,
【数634】
JP0005886309B2_000636t.gif
を得る.

【0376】
11.4.26 境界条件と随伴境界条件および固有関数セット
11.4.10節の解析解と比較するために,11.4.9節で与えた境界条件で問題を解く.即ち,y方向の一様重力下の円環で,内辺の境界条件を,表面力がゼロ,変位もゼロとするには,外辺をどのように支持すれば良いか求める.従来は求解不能とされた問題である.
固有関数に課される境界条件は,
【数635】
JP0005886309B2_000637t.gif
となる.添字Eは4.2節と同じで,固有関数によって生成される応力や表面力を表す.上式は,条件式(1018),(1019)や,(1069)式と同等である.(1085)式を得るのと同様にして,(91)式の境界項R は,極座標成分で
【数636】
JP0005886309B2_000638t.gif
となる.これは,(1070)式と同等で,外辺(r=R)で表面力がゼロかつ変位もゼロとする条件である.(1526)式と(1528)式では条件が異なるので,非自己随伴境界条件である.

【0377】
y方向の一様重力の解を表わせるのは,モードASの関数セットだけであり,他の3つのモード,SA,SS,AAは寄与しない.モードASには,(1251)式のJ1AS,(1263)式のJ2AS,(1295)式のY1AS,(1307)式のY2AS,(1339)式のI1AS,(1351)式のI2AS,(1383)式のK1AS,(1395)式のK2AS,の合計8セットが属し,任意の整数nに対しても解関数となるから,総数は8nとなる.しかし,一様重力の解を表し得るのはモードASの中でもn=0の関数セットだけで,n0の関数セットは寄与しない.

【0378】
主固有関数と双対固有関数および,これらに基づく応力関数をひとまとめにしてFとおく.そして,(1251)式のJ1ASに係数cJ1を,(1263)式のJ2ASに係数cJ2を,(1295)式のY1ASに係数cY1を,(1307)式のY2ASに係数cY2を,(1339)式のI1ASに係数cI1を,(1351)式のI2ASに係数cI2を,(1383)式のK1ASに係数cK1を,(1395)式のK2ASに係数cK2を,それぞれ乗じて加え合わせてFと等置する.即ち,
【数637】
JP0005886309B2_000639t.gif
を作成し,課された境界条件をFが満足よう係数cを決定する.境界条件は(1526),(1528)式の如く極座標成分で与えるのが便利で,J1ASを(1432)と(1433)式で,J2ASを(1436)と(1437)式で,Y1ASを(1440)と(1441)式で,Y2ASを(1444)と(1445)式で, I1ASを(1448)と(1449)式で,I2ASを(1452)と(1453)式で,K1ASを(1456)と(1457)式で,K2ASを(1460)と(1461)式で,それぞれ表現する.境界条件は合計8個,係数cも合計8個であるから,連立方程式
【数638】
JP0005886309B2_000640t.gif
【数639】
JP0005886309B2_000641t.gif
とおき,方程式を整理すれば,行列Aの成分が次のように簡約できる.ここで,n=0の関数セットのみならず,n0の関数セットも境界条件を満たせることに注意する. y方向の一様重力荷重では n=0の関数セットが解を表すが,その他の荷重状態ではn0の関数セットが解を表し得る.従って,Aの成分はnを用いて表しておく.Aの第1,2列は,
【数640】
JP0005886309B2_000642t.gif
とおいて,
【数641】
JP0005886309B2_000643t.gif
である.Aの第3,4列は,
【数642】
JP0005886309B2_000644t.gif
とおいて,
【数643】
JP0005886309B2_000645t.gif
である.Aの第5,6列は,
【数644】
JP0005886309B2_000646t.gif
とおいて,
【数645】
JP0005886309B2_000647t.gif
である.Aの第7,8列は,
【数646】
JP0005886309B2_000648t.gif
とおいて,
【数647】
JP0005886309B2_000649t.gif
である。

【0379】
Aの行列式はχ、χの関数となるが,(1533),(1534)式から,結局は固有値λの関数となる。行列式をゼロとする固有方程式を解けば,固有値λを得る。さらに,(1531)式より固有値λに応じた固有ベクトルcが得られる。固有方程式および固有ベクトルは非常に長大で,残念ながら,ここにその数式を記すことはできない.しかし,解いた結果を数値的に記すことはできる。

【0380】
11.4.27 主固有関数と双対固有関数

【0381】
11.4.10節の解析解と同様に,平面応力状態で,材料定数を,
【数648】
JP0005886309B2_000650t.gif
とする。

【0382】
n=0とした(1531)式から得られる固有方程式を解いて,低次から4つの固有値λは、
【数649】
JP0005886309B2_000651t.gif
付近にある。それぞれの固有値λに応じた固有ベクトルcを求めて(1530)式に代入すればFが定まる。Fには,主固有関数と双対固有関数および,これらに基づく応力関数が含まれている。

【0383】
第1モード(λ=5.317×10)での,変形状態と応力分布の様子を図47に示す。上段は主固有関数による様子で,内辺にて表面力がゼロかつ変位もゼロとなっている.下段は双対固有関数による様子で,外辺にて表面力がゼロかつ変位もゼロとなっている.上段の変形の様子は,図42の解析解の変形によく似ているが,応力分布は相当異なる.例えば,解析解ではσに応力の分布がみられたのに対し,主固有関数ではσθに同様の応力の分布がみられる。

【0384】
第2モード(λ=1.882×10)での,変形状態と応力分布の様子を図48に示す。上段は主固有関数による様子で,内辺にて表面力がゼロかつ変位もゼロとなっている.下段は双対固有関数による様子で,外辺にて表面力がゼロかつ変位もゼロとなっている.主固有関数によるσθに注目すると,第1モードでのσθの様子とよく似た分布であることがわかる.また、σの分布は解析解に似ていることもわかる。従って,第1モードと第2モードの差により,解析解によく似た様子が得られることが期待できる。

【0385】
第3モード(λ=5.080×10)での,変形状態と応力分布の様子を図49に示す。変形の様子は幾分複雑になっている。

【0386】
第4モード(λ=1.157×10)での,変形状態と応力分布の様子を図50に示す。変形の様子はさらに複雑になっている。

【0387】
11.4.28 固有関数法

【0388】
5.1節の固有関数法により解が求められる.11.4.10節の解析解と同様に,荷重係数c、cを、
【数650】
JP0005886309B2_000652t.gif
とおいて,変位を外径Rで無次元化する。

【0389】
第2モードまで用いた計算結果を図51に示す。応力σθはキャンセルされ、σの分布は解析解に似ていることがわかる。

【0390】

第3モードまで用いた計算結果を図52に示す.応力σの分布が解析解に近づいていることがわかる.


【0391】
第10モードまで用いた計算結果を図53に示す.応力分布が解析解に近づいていることがわかる.

【0392】
第30モードまで用いた計算結果を図54に示す.第10モードまで用いた結果とほとんど変わらず,応力分布も変形の様子も,全般に解析解に近いことがわかる.


【0393】
11.5 湧き出しと渦度
本節は,10章に含める方が良い.
11.5.1 微分方程式
2次元渦なし流れの問題は,直交座標系で,連続の式
【数651】
JP0005886309B2_000653t.gif
を満たすものと定義すれば,上式を(19)式に代入して,常に等号が成立する.従って,速度ポテンシャルを用いる事は,それだけで渦なし条件(19)式を満足する.また,(1544)式を(18)式に代入すれば,ラプラス方程式[Laplace’s equation]
【数652】
JP0005886309B2_000654t.gif
を満たすものと定義すれば,上式を(18)式に代入して,常に等号が成立する.従って,流れ関数を用いる事は,それだけで連続の式(18)を満足する.また,(1546)式を(19)式に代入すれば,ラプラス方程式
【数653】
JP0005886309B2_000655t.gif
が得られ,渦なし条件を表す.
湧き出しを注意深く表現した場合,(18)式および(1545)式の右辺はディラックのδ関数[Dirac’s delta function]となる.同様に,渦度を注意深く表現した場合,(19)式および(1547)式の右辺はディラックのδ関数となる.湧き出しや渦度の特性によって,方程式を切り替える必要を無くすためには,(18)式を
【数654】
JP0005886309B2_000656t.gif
と捉えて,fを渦度分布と認識すれば良い.また,(1548),(1549)式の連立微分方程式で扱えば,境界関数を速度u,uで指定できる利点も生じる.
主解となる速度u,uを主速度と呼びu,uと書く.方程式の右辺をf,fと書けば,連立偏微分方程式(1548),(1549)式は,(23)式の形となり,
【数655】
JP0005886309B2_000657t.gif
と表せる.上式は,流体力学の問題にも固有関数法を適用する基礎となる.

【0394】
11.5.2 調和方程式の極座標解
11.4.8節では重調和方程式の解を得たが,その幾つかは調和方程式(ラプラス方程式)を満足する.微分作用素の極座標形式は(992)式
【数656】
JP0005886309B2_000658t.gif
となる.ここで,νは任意の定数である.
(1552)式から得られる微分方程式
【数657】
JP0005886309B2_000659t.gif
である.同様に,(1552)式から得られる微分方程式
【数658】
JP0005886309B2_000660t.gif
となる.速度ポテンシャルの導関数は速度を表すため,θ方向に周期性の解を得る必要から,νを整数nとして,
【数659】
JP0005886309B2_000661t.gif
となる.
特殊な場合として,(1552)式でν=0とすれば,微分方程式
【数660】
JP0005886309B2_000662t.gif
となる.速度ポテンシャルの導関数は速度を表すため,θ方向に周期性の解を得る必要から, θ・lnrは解として不適切である.また,l・lは速度ゼロしか表し得ない.よって,
【数661】
JP0005886309B2_000663t.gif
の2つが解となる.上式第1の解lnrは湧き出しを表し,第2の解θは渦度を表す.

【0395】
11.5.3 調和方程式の直交座標解
調和方程式(ラプラス方程式)
【数662】
JP0005886309B2_000664t.gif
となる.ここで,νは任意の定数である.
(1566)式から得られる微分方程式
【数663】
JP0005886309B2_000665t.gif
と理解すれば,(1572)式から得られる微分方程式
【数664】
JP0005886309B2_000666t.gif
となる.
特殊な場合として,(1572)式でν=0とすれば,微分方程式
【数665】
JP0005886309B2_000667t.gif
となる.

【0396】
11.5.4 調和固有方程式の直交座標解
調和固有方程式
【数666】
JP0005886309B2_000668t.gif
とおくと,
【数667】
JP0005886309B2_000669t.gif
となる.

【0397】
11.5.5 湧き出しの解
座標原点にある湧き出しの強さをQ[m/s]を として,速度ポテンシャルφは,(1564)式より,
【数668】
JP0005886309B2_000670t.gif
となり,確かに,(18),(19)式を満足する.しかし,ガウスの発散定理[Gauss’s divergence theorem]
JP0005886309B2_000671t.gif
からみた場合には,少し様子が異なる.即ち,左辺は(1596)式によりゼロとなるのに対し,右辺はQとなる.左辺もQとするためには,div(v)=Q・δ(x,y)であれば良い.従って,微分方程式を,
【数669】
JP0005886309B2_000672t.gif
と認識する方が,座標原点にある湧き出しを表現できていることになる.この形式は,(1548),(1549)式に含まれる.なお,(1598)式で,nは境界Cでの外向き単位法線ベクトルで,領域を円にとると,
【数670】
JP0005886309B2_000673t.gif
である.

【0398】
11.5.6 渦度の解
座標原点にある渦度の強さ(循環)をΓ[m/s]として,速度ポテンシャルφは,(1564)式より,
【数671】
JP0005886309B2_000674t.gif
となり,確かに,(18),(19)式を満足する.しかし,平面内のストークスの定理[Stokes’s theorem]
【数672】
JP0005886309B2_000675t.gif
【数673】
JP0005886309B2_000676t.gif
と認識する方が,座標原点にある湧き出しを表現できていることになる.この形式は,(1548),(1549)式に含まれる.なお,(1606)式で,tは境界Cでの単位接線ベクトルで,領域を円にとると,
【数674】
JP0005886309B2_000677t.gif
である.

【0399】
11.5.7 矩形領域での湧き出しと渦度の様子
【数675】
JP0005886309B2_000678t.gif
となる.

【0400】
湧き出しによる代表速度をQ/(2πR)と定義して,(1595)式の速度u,uを無次元化すると,
【数676】
JP0005886309B2_000679t.gif
となる.観測領域を正方形として,無次元速度の様子を図55に示す.

【0401】
渦度による代表速度をΓ/(2πR)と定義して,(1603)式の速度u,uを無次元化すると,
【数677】
JP0005886309B2_000680t.gif
となる.観測領域を正方形として,無次元速度の様子を図56に示す.

【0402】
11.5.8 随伴微分作用素と随伴境界条件
【数678】
JP0005886309B2_000681t.gif
【数679】
JP0005886309B2_000682t.gif
【数680】
JP0005886309B2_000683t.gif
【数681】
JP0005886309B2_000684t.gif
以上により,(59)式の境界項 は(1629)式と(1638)式の2通りで理解でき,使い易い方を利用すれば良いことがわかる.双対問題の変関数u,uをどのような物理量に関連づけるかは,(1634)式のように,本解法を利用する者が決定する事項であり,それがどのような物理量と関連づけられていようとも,主問題の解に影響することはない.

【0403】
11.5.9 境界条件の同次化と境界項
非同次境界条件を満たす項に添字Bを付けてuBjとおき,同次境界条件を満たす項に添字Hを付けてuHjとおく.主速度 をその和で,
【数682】
JP0005886309B2_000685t.gif
【数683】
JP0005886309B2_000686t.gif
【数684】
JP0005886309B2_000687t.gif
【数685】
JP0005886309B2_000688t.gif
と置き換えると,(1654)式は(1549)式の渦条件を表しており,(1655)式は(1548)式の連続条件を表している.そこで,双対問題の物理的な解釈として,ここでは,-uH2 を x方向の速度,-uH1をy方向の速度と理解する.この視点を持てば,(1647)式を変形して
【数686】
JP0005886309B2_000689t.gif
以上により,(62)式の境界項Rは(1651)式と(1660)式の2通りで理解でき,使い易い方を利用すれば良いことがわかる.双対問題の変関数uH1,uH2をどのような物理量に関連づけるかは,(1656)式のように,本解法を利用する者が決定する事項であり,それがどのような物理量と関連づけられていようとも,主問題の解に影響することはない.

【0404】
11.5.10 連立固有値問題の連立解除
代表長さをRとして,重み定数を
【数687】
JP0005886309B2_000690t.gif
【数688】
JP0005886309B2_000691t.gif

【0405】
11.5.11 解関数
(1669)式と(1673)式は同じ微分方程式であるから,結局,φ,φと、φ,φ は,同じ関数セットから構成されることがわかる.よって,(1669)式を解く事に照準を合わせて,
【数689】
JP0005886309B2_000692t.gif

となる.これら,主・双対連立微分方程式を同時に満たす解は,(1593)式の組み合せから,以下の8つ得られる.
【数690】
JP0005886309B2_000693t.gif
【数691】
JP0005886309B2_000694t.gif
主連立微分方程式(1675),(1676)の固有値ωがゼロとなる解は,同次解である.調和方程式の解
【数692】
JP0005886309B2_000695t.gif
【数693】
JP0005886309B2_000696t.gif
以上の解には,まだ境界条件は反映されていないが,主連立微分方程式(1675),(1676)を満たす.

【0406】
11.5.12 境界条件の同次化
11.5.1節では,流体力学の問題も,連立偏微分方程式(1548),(1549)式は,(23)式の形となり,
【数694】
JP0005886309B2_000697t.gif

【0407】
[湧き出し]
座標原点にある湧き出しによる流れ場を固有関数法で再現できるかどうか確認しようとするとき,連立偏微分方程式(23)として,(1599),(1600)式
【数695】
JP0005886309B2_000698t.gif
【数696】
JP0005886309B2_000699t.gif
【数697】
JP0005886309B2_000700t.gif

【0408】
11.5.13 固有関数名の規則
次節に示す固有関数は、固有値α,βが
【数698】
JP0005886309B2_000701t.gif

【0409】
11.5.14 主固有関数と双対固有関数
【数699】
JP0005886309B2_000702t.gif
【数700】
JP0005886309B2_000703t.gif
【数701】
JP0005886309B2_000704t.gif
(1)から(8)の主固有関数と双対固有関数の第1モードの様子を、図57~図64に示す。

【0410】
11.5.15 固有関数法
5.1節の固有関数法により解が求められる.得られた解に,11.5.7節と同様の無次元化を施して,流れの様子を示す.

【0411】
[湧き出し]
観測領域を正方形として,固有関数のモード番号m×nを最大30×30とした計算を行った.外力項fが(1698)式であることから,(1704)式の固有関数f11SSのみで特解部分uHjを表現できる.無次元速度の様子を図65に示す.この結果は,(1599),(1600)式の特解である.
図65は、解析解による様子(図55)と良く似ており,原点から流体が湧き出す様子が表現できている.しかし、特解部分のみ表現しているため,解析解との差がどのようなものか調べておく必要がある.(1595)式の解析解から,固有関数解を差し引いた様子を図66に示す.
図66より,境界周辺部での差が大きいことがわかる.これを残差流と呼ぶ.残差流を表現するのは,同次解部分のuBjである.特解部は(1704)式の固有関数F11SSのみで表現されており、固有値は、

【数702】
JP0005886309B2_000705t.gif

【数703】
JP0005886309B2_000706t.gif
【数704】
JP0005886309B2_000707t.gif
【数705】
JP0005886309B2_000708t.gif
【数706】
JP0005886309B2_000709t.gif
【数707】
JP0005886309B2_000710t.gif
【数708】
JP0005886309B2_000711t.gif
【数709】
JP0005886309B2_000712t.gif
【数710】
JP0005886309B2_000713t.gif
【数711】
JP0005886309B2_000714t.gif
【数712】
JP0005886309B2_000715t.gif
【数713】
JP0005886309B2_000716t.gif
【数714】
JP0005886309B2_000717t.gif
【数715】
JP0005886309B2_000718t.gif
【数716】
JP0005886309B2_000719t.gif
【数717】
JP0005886309B2_000720t.gif
【数718】
JP0005886309B2_000721t.gif
【数719】
JP0005886309B2_000722t.gif
【数720】
JP0005886309B2_000723t.gif
【数721】
JP0005886309B2_000724t.gif
【数722】
JP0005886309B2_000725t.gif
【数723】
JP0005886309B2_000726t.gif
【数724】
JP0005886309B2_000727t.gif
【数725】
JP0005886309B2_000728t.gif
【数726】
JP0005886309B2_000729t.gif
【数727】
JP0005886309B2_000730t.gif
【数728】
JP0005886309B2_000731t.gif
【数729】
JP0005886309B2_000732t.gif
【数730】
JP0005886309B2_000733t.gif
【数731】
JP0005886309B2_000734t.gif
【数732】
JP0005886309B2_000735t.gif
【数733】
JP0005886309B2_000736t.gif
【数734】
JP0005886309B2_000737t.gif
【数735】
JP0005886309B2_000738t.gif
【数736】
JP0005886309B2_000739t.gif
【数737】
JP0005886309B2_000740t.gif
【数738】
JP0005886309B2_000741t.gif
【数739】
JP0005886309B2_000742t.gif
【数740】
JP0005886309B2_000743t.gif
【数741】
JP0005886309B2_000744t.gif
【数742】
JP0005886309B2_000745t.gif
【数743】
JP0005886309B2_000746t.gif
【数744】
JP0005886309B2_000747t.gif
【数745】
JP0005886309B2_000748t.gif
【数746】
JP0005886309B2_000749t.gif
【数747】
JP0005886309B2_000750t.gif
【数748】
JP0005886309B2_000751t.gif
【数749】
JP0005886309B2_000752t.gif
【数750】
JP0005886309B2_000753t.gif
【数751】
JP0005886309B2_000754t.gif
【数752】
JP0005886309B2_000755t.gif
図面
【図1】
0
【図2】
1
【図3】
2
【図4】
3
【図5】
4
【図6】
5
【図7】
6
【図8】
7
【図13A】
8
【図13B】
9
【図13C】
10
【図14】
11
【図15】
12
【図16】
13
【図17】
14
【図18】
15
【図19】
16
【図20】
17
【図21】
18
【図22】
19
【図23】
20
【図24】
21
【図25】
22
【図26】
23
【図27】
24
【図28】
25
【図29】
26
【図30】
27
【図31】
28
【図32】
29
【図33】
30
【図34】
31
【図35】
32
【図36】
33
【図37】
34
【図38】
35
【図39】
36
【図40A】
37
【図40B】
38
【図41】
39
【図43】
40
【図44】
41
【図45】
42
【図46】
43
【図55】
44
【図56】
45
【図57】
46
【図58】
47
【図59】
48
【図60】
49
【図61】
50
【図62】
51
【図63】
52
【図64】
53
【図65】
54
【図66】
55
【図67】
56
【図68】
57
【図69】
58
【図70】
59
【図71】
60
【図72】
61
【図9】
62
【図10】
63
【図11】
64
【図12】
65
【図42】
66
【図47】
67
【図48】
68
【図49】
69
【図50】
70
【図51】
71
【図52】
72
【図53】
73
【図54】
74