並列有限要素計算システム
【課題】ＣＧＣＧ法を発展させることで，より短時間で解が求まる並列有限要素計算システムを提供する。
【解決手段】解析対象領域Ωと，そのメッシュ分割を与える。並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する。並列プロセスにわたる第０階層領域分割を生成する。ν個の階層の領域分割を生成する。全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る。全世代の２組の剛性行列と２組の外力ベクトルを作る。ν個のコース空間の方程式を直接法で解く。第ν世代の共役空間の方程式を反復法(前処理付きＣＧ法)で解く。ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える。この方法は，最初に「領域分割」を全部与えておき直列に計算する方法である。領域分割を「共役空間」における反復法の収束の様子を見たうえで順次与え，つぎつぎと計算していく方法も本発明に含まれる。
【解決手段】解析対象領域Ωと，そのメッシュ分割を与える。並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する。並列プロセスにわたる第０階層領域分割を生成する。ν個の階層の領域分割を生成する。全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る。全世代の２組の剛性行列と２組の外力ベクトルを作る。ν個のコース空間の方程式を直接法で解く。第ν世代の共役空間の方程式を反復法(前処理付きＣＧ法)で解く。ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える。この方法は，最初に「領域分割」を全部与えておき直列に計算する方法である。領域分割を「共役空間」における反復法の収束の様子を見たうえで順次与え，つぎつぎと計算していく方法も本発明に含まれる。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は，自由度１００万以上の大規模構造問題を効果的に解く，大規模並列有限要素法ソルバアルゴリズムに関するもので，特に，並列ＣＧアルゴリズムを基盤とし，領域分割（Domain Decomposition）に基づくアルゴリズムである。
【０００２】
本明細書において，該アルゴリズムを「ＣＧＣＧ法」（Coarse Grid CG method,CGCG method）と称する。Ｋを当該構造問題の剛性行列とし，領域分割による部分領域単位の自由度の空間をコース（coarse）空間，それに共役直交な空間をファイン空間と称する。ＣＧＣＧ法は，コース空間とファイン空間を分離し，それぞれの解を重ね合わせるものである。この際，領域分割は，コース空間の生成及びファイン空間上のＣＧ法の並列化に用いられる。
【背景技術】
【０００３】
このようなＣＧＣＧ法を利用した並列有限要素計算システムとして，特許文献１に記載のものが知られている。この発明は，本出願人が提案したものである。この先行発明は，ＣＧＣＧ法を使用することで，自由度が１００万以上の大規模構造問題を解く際，解が発散することなく求まり，反復計算の回数が少なく，短時間の計算時間で解が求まる，大規模構造問題を解くシステムを提供するものである。
【０００４】
このＣＧＣＧ法では，解析対象の領域分割を行ない，それに基づいて部分領域の内部と境界の区別をすることなしにコース空間を設定し，前処理として簡便な対角スケーリングを採用した前処理つき射影ＣＧ法を用いる。すなわち，並列処理用に分割された領域にわたる“コース運動”を考え，解空間をコース空間とその共役直交空間に分解する。コース空間は直接法で解き，共役直交空間には反復法（共役射影並列ＣＧ法）を用いる。コース運動として剛体運動をとることで，従来の解法に比較して高速に解を得ることができる。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】特許第３８２４５８２号公報
【０００６】
【非特許文献１】J. Mandel, "Balancing Domain Decomposition". Communications on Numerical Methods in Engineering, 1993, 9, pp.233-241
【発明の概要】
【発明が解決しようとする課題】
【０００７】
ＣＧＣＧ法は大規模な有限要素法問題を解くための非常に高速なアルゴリズムであるが，薄板構造への適用性に困難があり，また，よりローカルな問題が全体の収束性に与える場合など，問題依存による収束性の向上の必要がある場合がある。本発明アルゴリズムは，直接法を適用する小さなコース空間と反復法を適用する大きなファイン空間を階層的に多段に取っていくものであるが，コース空間を多段にすることによって，直接法の負荷が大きくなり，薄板構造がより解きやすい方向に改善される。また，構造のよりローカルな運動が表現可能になり，構造依存による収束性の向上がはかれる。
【０００８】
また，ＣＧＣＧ法は，反復法を基本にするため，固有値問題，複数荷重問題など，剛性行列は同一で，荷重条件など，右辺のみが異なる問題などが不得意であるが，直接法の負荷を大きくすることによってこれらがより解きやすい方向に改善される。
【０００９】
本発明は前記のような要望に基づいて提案されたものであって，その目的は，ＣＧＣＧ法を発展させることで，より短時間で解が求まる並列有限要素計算システムを提供することにある。
【課題を解決するための手段】
【００１０】
前記の目的を達成するために，本発明の並列有限要素計算システムは，次の構成を有することを特徴とする。
(1) 解析対象領域Ωと，そのメッシュ分割を与える。
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する。
(3) 並列プロセスにわたる第０階層領域分割を生成する。
(4) ν個の階層の領域分割を生成する。
【００１１】
(5) 全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る。
(6) 全世代の２組の剛性行列と２組の外力ベクトルを作る。
(7) ν個のコース空間の方程式を直接法で解く。
(8) 第ν世代の共役空間の方程式を反復法(前処理付きＣＧ法)で解く。
(9) ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える。
【００１２】
また，次の構成も，本発明の一態様である。
(1) 解析対象領域Ωと，そのメッシュ分割を与える。
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する。
(3) 並列プロセスにわたる第０階層領域分割を生成する。
(4) 第ｋ階層領域分割をもとに第ｋ世代コース空間Ｗ(k)を作る。
【００１３】
(5) 直前世代の共役空間Ｖ(k-1)をコース空間とそのＫ(k-1)共役直交空間に分割する。
(6) 第ｋ世代の２つの剛性行列と外力ベクトルを作る。
(7) 第ｋ世代コース空間の方程式を直接法で解く。
(8) 第ｋ世代の共役空間の方程式について，反復法の収束性をチェックする。
(9) 十分な収束性のチェックを行う。収束性が得られなければ，ｋをｋ＋１として(5)に戻る。
(10) 十分な収束性が得られた場合には，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を与える。
【発明の効果】
【００１４】
本発明によれば，直接法を適用する小さなコース空間と反復法を適用する大きなファイン空間を階層的に多段に取っていくものであるが，コース空間を多段に取ることによって，直接法の負荷が大きくなり，薄板構造がより解きやすい方向に改善される。また，構造のよりローカルな運動が表現可能になり，構造依存による収束性の向上がはかれる。また，固有値問題，複数荷重問題など剛性行列は同一で，荷重条件など，右辺のみが異なる問題が，直接法の負荷を大きくすることによってこれらがより解きやすい方向に改善される。
【図面の簡単な説明】
【００１５】
【図１】本発明の領域分割のイメージを示す模式図。
【図２】本発明の実施例１の処理を示すフローチャート。
【図３】本発明の実施例２の処理を示すフローチャート。
【図４】本発明のシステムを実現するコンピュータのハードウェアを示すブロック図。
【図５】本発明の実施例１の各手段を示すブロック図。
【図６】本発明の実施例２のの各手段を示すブロック図。
【発明を実施するための形態】
【００１６】
本発明の前提となる事項について，説明する。
（１）概要
本発明において，Ωを解析対象領域とする。この上の構造問題を有限要素法にしたがって離散化する。目標はここから帰着される線形方程式
Ｋｕ＝ｆ……（１．１）
を解くことである。
【００１７】
(共役直交射影子)
ベクトル空間Ｖ上の正定値対称変換ＡとＶ上の射影Ｐが与えられたとき，Ｖの部分空間ＰＶと（１−Ｐ）Ｖに直和分解され，つぎの２式がなりたつ。
Ｐ2=Ｐ……(５．４)
ＡＰ=ＰTＡ……(５．５)
あるいはつぎのようにもいえる。ベクトル空間Ｖ上の正定値対称変換ＡとＶの部分空間Ｗが与えられたとき，ＶからＷへの射影Ｐは(５．４)式をみたし，ＡとＰの間には関係(５．５)がなりたつ。
逆に，(５．４)式と(５．５)式をみたすＶ上の線形変換Ｐと，正定値対称変換Ａが与えられたとき，Ｖは部分空間ＰＶと（１−Ｐ）Ｖに直和分解され，これら２つの部分空間はＡに関して共役直交である。
【００１８】
（部分空間の共役直交基底）
ＷをＶのｒ次元部分空間とするとき，Ｗの基底ｆ1，…，ｆrを互いにＡ共役直交であるようにえらべる。すなわち，Ｗの基底ｆ1，…，ｆrを適当にとることにより，つぎの関係がなりたつようにできる。
(ｆj，Ａｆi)=０ (ｉ≠ｊ)……(５．６)
【００１９】
（ベクトル空間の共役直交直和分解）
ＷをＶの，Ｖに一致しない部分空間とするとき，Ｖは，ＷとＡに関して共役直交な部分空間Ｖ′に直和分解できる。
【数１】
【００２０】
(共役直交射影の表現)
Ｖをベクトル空間，Ｗをｒ次元部分空間，ｆ1，…，ｆrをＷの共役直交基底，ＡをＶ上の正定値対称変換，Ｐを共役直交射影子とするとき，Ｐは
【数２】
で与えられる。
【００２１】
（２）ＣＧＣＧ法
（領域分割）
Ωは有限要素分割されている。その分割に基づいてΩを，オーバラップしないが全部の和集合はΩに一致するような部分領域に分ける。その第１の分割は並列プロセスにわたる分割である。これをさらに分割して，ΩI，I=１，…，ｍとする。ΩIのある面または稜線は別の部分領域と接していて，これらの面または稜線は複数の部分領域に共有される。Ω内部での部分領域どうしの境界面を内部境界と呼ぶ。この面はあるΩIのある面のグループは外部境界面(Ωの境界面)の一部であることもある。
【００２２】
以下では２つの空間ＶIおよび~ＶIを定義するが，これは，以下で述べる部分空間から全空間Ｖへの延長とその逆の制限を定義し，これをもとにコースグリッドの項で述べた方法でコース空間を定義するのに用いるためである。
【００２３】
Ω上で定義された有限要素基底関数を｛ψj｝とする。ψjをΩIに制限した関数のうち，恒等的にゼロでないものを~ψIｊと書く。これらの｛~ψIj｝はΩIの有限要素基底関数をなす。｛~ψIj｝で離散化されたΩI上の変位場の空間を~ＶIと書き，切断された局所空間空間と呼ぶ。
【００２４】
一方，ΩI(およびその境界)に含まれる節点上に定義される許容変位場がなすベクトル空間をＶIと書く。ΩIの基底関数~ψjはもとの基底関数ψjの外側を切り取った(ただし~ψjの定義域には境界面も含む)ものなので，~ＶIはＶIとは異なる。ＶIはＶの部分空間だが，~ＶIは部分空間ではないことに注意せよ。ΩIはオーバラップしないが，ＶIは内部境界のまわりでオーバラップする。
【００２５】
~ＶIの各節点における基底関数に，その節点に対応するＶIの(実はＶの)基底関数の対応させる写像は１対１である。ＶIはＶの部分空間なので，この写像
【数３】
を~ＶIからＶへの，あるいは~ＶIからＶIへの延長(または延長写像)(prolongation)と呼ぶ。~ＶIはＶの部分空間ではないが，ＮIによって部分空間に埋め込まれる。つまり，ＮI~ＶIは部分空間である。一方，ＮIの転置
【数４】
をＶから~ＶIへの制限(または制限写像)(restriction)と呼ぶ。
【００２６】
グローバルな節点を任意にとる。その節点は，あるIに対する部分領域ΩI上にあるとする。対応する自由度を制限(ＮI)Tによってローカルな自由度に写し，その後，延長ＮIによってグローバルな自由度に写す。内部境界上の節点の自由度は，(ＮI)Tによって~ＶIを含むいくつかの切断された局所空間のローカルな自由度に写され，延長によって元の節点を共有する複数の部分領域に対応する自由度にもどる。この重複を解消するため，つぎのような対角行列ＤIを定義する。領域分割の際，内部境界上の重複する節点の数を求めておき，その自由度に対応する成分を(１/重複節点数)，重複のない成分を１とする。この作用は
【数５】
と表わせる。
【００２７】
Ｖは部分空間ＶI全部の和からなり，ＶIはＮIによる~ＶIの延長で得られるが，重複を考慮すると，これらの関係は
【数６】
と表わせる。
【００２８】
（コースグリッド）
部分領域ΩI，Ｉ=１，…，ｍの剛体運動を考える。要素分割したΩの全節点数をｎとし，グローバルな節点番号をｉ(ｉ=１，…，ｎ)で表わす。初期の節点座標を
【数７】
とする。肩に付けた添え字Ｉは，任意の節点はいずれかの部分領域ΩIにあるので，これを明示するために付けたものである。
【００２９】
剛体変位後の座標を
【数８】
とする。
【００３０】
剛体変位は平行移動３自由度，回転３自由度をもち，
【数９】
と表わせる。
【００３１】
ただし
【数１０】
であり，Ｕ1，Ｕ2，Ｕ3は各軸に沿った平行移動，
【数１１】
は各軸の周りの，ベクトル
【数１２】
の右ねじ方向の周りの，大きさ｜θ｜の回転を表わす。
【００３２】
回転は微小とすると，(６．３)は
【数１３】
と書ける。したがって，対応する節点変位をｕIｉとすると
【数１４】
となる。
【００３３】
以上をもとに，部分領域ΩIの剛体変位場(節点座標を変数にとり，剛体変位を与える~ＶIのベクトル)の基底をつぎのように定義する。ただし，ΩIの一部が全体領域の外表面に接しており，この部分に変位拘束が課されている場合を除く。
【数１５】
【００３４】
このように，剛体変位場ベクトルの自由度は１節点あたり６，つまり，部分領域ΩIの剛体変位を与える空間の基底(６．４)の次元は６である。上で除外した，ΩIの表面が全体空間の外表面に一部をなし，ここに変位拘束が課されている場合には，その次元が５以下になることもありうるが，独立な基底をもとにこれを拡大して６にそろえる。
【００３５】
このようにして得られた空間を~ＷIと書くことにする。ΩIの剛体変位はもとよりΩIの許容変位場なので，ベクトル(６．４)はΩIの有限要素基底関数｛~ψIj｝で表現できる。すなわち，~ＷIは~ＶIの部分空間である。ＮIによって~ＶIがＶIに延長されたように，~ＷIもＮIによってＶIに延長される。
【００３６】
ただし，自由度の重複の解消は行わなければならない。(６．２)と同じように，内部境界自由度の重複を，１の分割ＤIを用いて解消し，
【数１６】
とおく。これを剛体運動によるＶのコース空間と呼ぶ。コース空間ＷはＶの部分空間である。
【００３７】
（グローバル空間の共役直交分解）
解析対象モデルΩはメッシュ分割によりｎ個の節点をもち，したがって自由度は３ｎである。グローバルな解空間Ｖの次元および(１．１)の左辺のグローバルな剛性行列Ｋの次元も３ｎである。
【００３８】
コースグリッドの項で説明した通り，コース空間Ｗとその共役空間Ｖ′に分解する。これは(５．７)式に対応する。Ｗの基底は(６．４)式を部分領域の数ｍ個分を並べたものであり，したがってＷの次元は６ｍである。Ｖ′は，(５．７)式のように作られ，ＷとＶ′はＫに関して共役直交である。Ｖ′の次元は３ｎ−６ｍであり，これをあらためて(ｆ1，…，ｆ6m)とする。６ｍは３ｎに比べて，直接法が適用できるほど十分小さいものとする。Ｖの基底は，Ｗの基底に，これと共役直交なＶ′の基底を並べて作られる。
【００３９】
Ｖ′の基底を(ｇ6m+1，…，ｇ3n)とする。Ｋの表現はこの基底の並び方に合わせ，(４．１)のように作られているものとする。つまり，基底の並べ替えに応じてＫはつぎのような表現をもつ。
【数１７】
【００４０】
射影Ｐの表現は，（５．８）式により６ｍ次元単位行列である。より正確には，Ｐは，次元を合わせて
【数１８】
で定義される。
【００４１】
Ｐは６ｍ×３ｎの横長の行列であり，３ｎ次元の解空間Ｖを６ｍ次元のコース空間Ｗに射影する。
Ｗ=ＰＶ
【００４２】
ここで，Ｐはべき等行列である。
Ｐの転置行列
【数１９】
は，ＷをＶの中のＷと同型な部分空間に射影する３ｎ×６ｍの縦長の行列である。
【数２０】
【００４３】
コース空間Ｗ上のＫの作用は，(４．１５)と同じように
~Ｋ=ＰＫＰT
と書ける。
【００４４】
共役空間Ｖ′への射影は１−Ｐだが，(６．７)と同様，次元を合わせて書くと，
【数２１】
となる。
【００４５】
このように，１−Ｐは(３ｎ−６ｍ)×３ｎの横長の行列であり，３ｎ次元の解空間Ｖを３ｎ−６ｍ次元の共役空間Ｖ′に射影する。１−Ｐの転置行列
【数２２】
は，Ｖ′をＶの中のＶ′と同型な部分空間に射影する３ｎ×(３ｎ−６ｍ)の縦長の行列である。すなわち，
【数２３】
【００４６】
共役空間Ｖ′上のＫの作用は
Ｋ′=(１−Ｐ)Ｋ(１−ＰT)
と書け，これは３ｎ−６ｍ次元の正定値対称行列である。
【００４７】
（５．６）式により，Ｋは~ＫとＫ′により，
【数２４】
と書ける。
【００４８】
（ＣＧＣＧ法）
前記の通り，解空間Ｖをコース空間Ｗとその共役空間とに直和分解する。
【数２５】
【００４９】
この上でグローバルな線形方程式，すなわち
Ｋｖ=ｆ……(６．１３)
を解く。
【００５０】
まず，(６．１３)をＷおよびＶ′にそれぞれＰおよび１−Ｐによって部分空間に写す。
【数２６】
【００５１】
この上で，ｖ∈ＶをＰTおよび１−ＰTによってそれぞれＷおよびＶ′の元に限定すると，
【数２７】
が得られる。
【００５２】
これを(６．１４)に入れると
【数２８】
が得られる。ここで
【数２９】
とおく。
【００５３】
一方，ＰＶ=Ｗの双対空間Ｗ*はＷと同一視される。ｆはＫの像空間ＫＶに含まれるので，ＰｆはＰＫＰTによるＷの像空間(ＰＫＰT)Ｗに含まれる。同様に，(１−Ｐ)Ｖ=Ｖ′の双対空間Ｖ′*はＶ′と同一視され，(１−Ｐ)ｆは((１−Ｐ)Ｋ(１−ＰT))Ｖ=Ｖ′に含まれる。そこで，
【数３０】
とおくと，(６．１４)は
【数３１】
に帰着する。
【００５４】
すでに述べたとおり，(６．１９)は直接法で解く。(６．２０)は，これと同等な方程式
ＧＫ′ｕ=Ｇｆ′……(６．２１)
について，つぎの前処理付きＣＧ法で解く。ここで，ＧはＶ′上の適当な正定値対称行列である。
【数３２】
【００５５】
ＷとＶ′の共役直交性により，数値誤差を考えなければ，残差ｒkは反復過程でつねにＶ′(正確にいえばＶ′の双対空間)にとどまり，その結果，勾配ベクトルｐk，解の近似ベクトルｕkもＶ′(こちらは双対空間ではない)にとどまる。収束解は剛体運動と重ねあわされ，幅広い運動スペクトルを表現するものと期待できる。
【００５６】
（３）本発明のアルゴリズム
本発明においては，以上述べたＣＧＣＧ法のアルゴリズムを拡張する。これをここでは，本発明のアルゴリズムを，仮に，強化版ＣＧＣＧ法(enhanced CGCG method)と呼んでおく。領域分割はν(ν≧１)階層に分割される。ν＝１の場合が，前節までに述べた通常のＣＧＣＧ法であり，部分領域の数はｍである。
【００５７】
(0) 解析対象領域Ωの並列プロセスにわたるｍ０個の部分領域Ω(0,1)，…，Ω(0,m0)。この分割を親領域分割(parent domain decomposition)と呼ぶ。
【００５８】
(1) 部分領域Ω(0,1)，…，Ω(0,m0)を含む，より細かい領域分割Ω(1,1)，…，Ω(1,m1)。これを第１世代領域分割(first generation domain decomposition)と呼ぶ。これがＣＧＣＧ法の説明に用いた部分領域Ω1，…，Ωmにあたり，第１世代部分領域の数はｍ1＝ｍである。分割Ω(1,1)，…，Ω(1,m1)は分割Ω(0,1)，…，Ω(0,m0)を含み，Ω(1,1)，…，Ω(1,m1)の全体が全領域Ωに一致する。したがってｍ1≧ｍ0である。ｍ1＝ｍ0の場合はＣＧＣＧ法に一致しするので，強化版では意味ある分割，すなわちｍ1＞ｍ0とする必要がある。
【００５９】
(2) 部分領域Ω(1,1),…,Ω(1,m1)を含む，より細かい領域分割Ω(2,1),…,Ω(2,m2)。ただし，Ω(2,1),…,Ω(2,m2)のうちの最初のある領域のグループはΩ(1,1)に含まれ，つぎのグループはΩ(1,2)に含まれ，最後のグループはΩ(1,m1)に含まれる。これを第２世代領域分割(second generation domain decomposition)と呼ぶ。
【００６０】
(ν) 部分領域Ω(ν-1,1),…,Ω(ν-1,mν-1)を含む，より細かい領域分割Ω(ν,1),…,Ω(ν,mν)。
ただし，Ω(ν,1),…,Ω(ν,mν)のうちの最初のある領域のグループはΩ(ν-1,1)に含まれ，つぎのグループはΩ(ν-1,2)に含まれ，最後のグループはΩ(ν-1,mν-1)に含まれる。これを第ν世代領域分割(ν-th generation domain decomposition) と呼ぶ。
【００６１】
なお、前記(ν)のΩ(ν-1,mν-1)及びΩ(ν,mν)の標記中、上付き文字のmν-1及びmνは、正確には、ｍν-1及びｍνである。同様に前記(0)から(2)や本文以降等の他の式においても、上付き文字中の下付き文字については、便宜的に下付き文字表示を省略している。
【００６２】
上の階層的領域分割に従って，コース空間をつぎのように階層的に作る。
【００６３】
(1) 第０世代を含む第１世代のコース空間を，コースグリッドの項で述べた方法で構成する。これをＷ(1)とする。
【００６４】
(2) Ω(1,1),…, Ω(1,m1)に対応するそれぞれのコース空間を，コースグリッドの項で述べた方法で各グループ単位の領域分割でコース空間を作る。ただし，各第２世代部分領域は第１世代部分領域の内部境界上の節点は含まず，内部境界面で囲まれている節点のみを用る。
これを第２世代のコース空間と呼び，Ｗ(2)と書く。Ｗ(2)はｍ1個のコース空間からなる，互いに他から切り離されたローカルな空間である。ｍ1個の各コース空間は，そこにグループ化されただけの数の第２世代部分領域を持ち，第２世代部分領域数は全体でｍ2なので，Ｗ(2)の次元は６ｍ2である。第１世代のコース運動の中には，Ｗ(2)上のコース運動では表現できないものがある。
【００６５】
(ν) Ω(ν-1,1),…,Ω(ν-1,mν-1)に対応するそれぞれのコース空間を，コースグリッドの項で述べた方法で各グループ単位の領域分割でコース空間を作る。ただし，各第ν世代部分領域は第ν−１世代部分領域の内部境界上の節点は含まず，内部境界面で囲まれている節点のみを用る。これを第ν世代のコース空間と呼び，Ｗ(ν)と書く。Ｗ(ν)はｍν-1個のコース空間からなる，互いに他から切り離されたローカルな空間である。ｍν-1個の各コース空間は，そこにグループ化されただけ数の第ν世代部分領域を持ち，27第ν世代部分領域数は全体でｍνなので，Ｗ(ν)の次元は６ｍνである。第ν−１世代のコース運動の中には，Ｗ(ν)上のコース運動では表現できないものがある。
【００６６】
本発明の強化版ＣＧＣＧ法がマルチグリッド法と異なる点は，マルチグリッド法が，細かいメッシュ分割をもとに，より粗い分割を作っていくのに対して，強化版ＣＧＣＧ法では，コースグリッドをより細かく作っていくことにある。すなわち，剛性方程式(６．１３)を解く際には，前処理行列Ｇとしては負荷の軽い対角スケーリングなどを使うのが効率的である。そのため，コース運動を解く際に，それが粗すぎると，十分な精度が得られないことが経験上わかっている。一方，非特許文献１に記載のＢＤＤ法では，分割領域ΩＩの境界面のみに剛体運動を与え，ΩＩ内部は，剛体運動で決まった境界面の変位拘束条件をもとに，部分構造としての剛性方程式を解いているため，ＢＤＤ法のコースグリッドにはＣＧＣＧ法のそれよりもより詳細な表現が与えられている。この意味で，ＣＧＣＧ法のコース運動はＢＤＤ法と比べ粗い。コースグリッドを細かく取っていくというアイデアは，この粗さを補完することに対応している。
【００６７】
Ｖは第１世代子分割に基づいて，共役直交分解される。これは前節のＣＧＣＧ法の手順と同じものである。以下で記号を改める。前節を参照して，
【数３３】
とおく。
【００６８】
第１世代子アルゴリズム(通常のＣＧＣＧ法アルゴリズム)はつぎのように書ける。
【数３４】
【００６９】
前記（６．２２）式における下から３行目のｗ(1)に関する方程式は直接法で解く。また，Ｖ(1)上の剛性方程式は第２世代分割に基づいてコース空間Ｗ(2)とそのＫ(1)共役直交空間に分解され，この上にＣＧＣＧ法が適用できる。
【数３５】
【００７０】
Ｖ(1)はＷ(1)と共役直交なので，もとの空間Ｖ＝Ｖ(0)は，
【数３６】
と書ける。
【００７１】
第２世代アルゴリズムはつぎのようになる。
【数３７】
【００７２】
前記（６．２４）式における下から４行目のｗ(2)に関する方程式は直接法で解く。
【００７３】
以下同様にして，第ν−１世代アルゴリズムが与えられたとして，第ν世代アルゴリズムはつぎのようになる。
【数３８】
【００７４】
前記（６．２５）式の下から４行目のｗ(ν)に関する方程式は直接法で解く。
【００７５】
第１世代アルゴリズム(６．２２)，第２世代アルゴリズム(６．２４)から第ν−１世代アルゴリズムに現れるｕ(1)，ｕ(2)，…，ｕ(ν-1)についての方程式は実際には解く必要はなく，それぞれ，つぎの世代のアルゴリズムに吸収される。結局，強化版ＣＧＣＧ法アルゴリズムはつぎのようにまとめられる。ν個のコース問題
【数３９】
をそれぞれ直接法で解き，最後の第ν世代アルゴリズムにおける３nν次元の方程式
Ｋνｕν＝ｆ′ν， ｕν∈Ｖν……（６．２７）
に前処理付きＣＧ法を適用して，これらで得られた解をすべて足しあわせてグローバルな解を
ｕ＝ｕ(0)＝ｗ(1)+ｗ(2)+…+ｗ(ν)+ｕ(ν)……(６．２８)
で得る。
【００７６】
Ｖ(ν)の次元３nνはつぎで与えられる。
３nν＝３n0−６m1−６m2−…−６mν……(６．２９)
【００７７】
また，解空間の共役直交分解による階層的な関係はつぎのように書ける。
【数４０】
【００７８】
本発明の強化版ＣＧＣＧ法においては，コース空間を多段にとることで，より細かいローカルな運動の表現が可能になって，ＣＧＣＧ法に比べて収束性が向上することが期待できる。また，階層を深くとれば直接法の負荷が大きくなるため，反復法では解きにくい問題も，より解きやすい方向に改善されるものと期待される。
【実施例１】
【００７９】
以下，本発明の実施例１を図面に従って説明する。
まず，図４を参照して，実施例１の並列有限要素計算システムを構成するコンピュータのハードウェアを説明する。このハードウェアは，バス３０を介して相互に接続されるＣＰＵ３１と，ＲＡＭ３２と，記憶装置３３と，入力装置３４と，出力装置３５とを具備する。記憶装置３３はハードディスクやメモリ等の外部記憶装置である。入力装置３４は，キーボードやマウス等のユーザによって操作されることで，各種データをＣＰＵ３１や記憶装置３３に出力する。出力装置３５は，モニタやプリンタに例示され，ＣＰＵ３１から出力される求解の結果をユーザに対し視認可能に出力する。
【００８０】
ＣＰＵ３１は，入力装置３４からの入力に応答して，記憶装置３３内の強化版ＣＧＣＧ法のアルゴリズムを利用した求解プログラムを実行する。この際，記憶装置３３からの各種データやプログラムはＲＡＭ３２に一時格納され，ＣＰＵ３１はＲＡＭ３２内のデータを用いて各種処理を実行する。
【００８１】
前記強化版ＣＧＣＧ法を実現するために，実施例１の並列有限要素計算システムにおいては，前記のようなハードウェアを有するコンピュータを，問題入力手段１，条件設定手段２，第０階層領域分割手段３，各階層領域分割手段４，共役空間生成手段５，全世代剛性行列・外力ベクトル作成手段６，直接法計算手段７，反復法計算手段８及び解出力手段９として，機能させる。
【００８２】
すなわち，コンピュータの有するキーボード，マウス，ネットワーク，記憶媒体などの入力装置３４から入力される問題を受け付ける問題入力手段１，同じく設定すべき条件を受け付ける条件設定手段２を備える。
【００８３】
前記第０階層領域分割手段３，各階層領域分割手段４，共役空間生成手段５，全世代剛性行列・外力ベクトル作成手段６，直接法計算手段７，反復法計算手段８及び解出力手段９は，コンピュータのメモリ（ＲＡＭ３２）やＣＰＵ３１上に展開されるプログラムによって，コンピュータ上に構築される手段で，これら各手段が前記強化版ＣＧＣＧ法のアルゴリズムに従って，与えられた問題に対する解を求める。
【００８４】
解出力手段９は，コンピュータに設けられたディスプレイ，プリンタ，記憶媒体，或いはネットワークを通じて接続された他のコンピュータなどの出力装置３５に接続され，これらのハードウェア上に求められた解を出力する。
【００８５】
次に，各手段が実行する処理を図２のフローチャートに従って説明する。
問題入力手段１は，ユーザからの問題として，解析対象領域Ωと，そのメッシュ分割を受け付ける（ステップ１）。条件設定手段２は，問題を解決する場合に要求される並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを含む設定条件を，ユーザから受け付ける（ステップ２）。この処理は，前記アルゴリズムの階層構造の深さを決めることに相当する。
【００８６】
前記第０階層領域分割手段３は，前記のようにして受け付けた解析対象領域Ωに対して，Ｖ(0)＝Ｖ，Ｋ(0)＝Ｋ，ｆ(0)＝ｆ，Ｐ(0)＝Ｐ，Ｗ(1)＝Ｗ，Ｖ(1)＝Ｖ′とおき，これに基づいて，並列プロセスにわたる第０階層領域分割を生成する（ステップ３）。
【００８７】
各階層領域分割手段４は，設定条件として受け付けた階層数νに基づいてν個の階層の領域分割を生成する（ステップ４）。この処理は，前記アルゴリズムの（６．２５）式のＷ(1)，…，Ｗ(ν)を決めることに相当する。
【００８８】
共役空間生成手段５は，前記各階層領域分割手段４によって生成された全ν階層の領域分割をもとに，第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作成する（ステップ５）。この処理は，前記アルゴリズムの（６．２５）式の第２式を作成することに相当する。
【００８９】
全世代剛性行列・外力ベクトル作成手段６は，全世代の２組の剛性行列と２組の外力ベクトルを作成する（ステップ６）。この処理は，前記アルゴリズムの（６．２６）式を作成することに相当する。
【００９０】
直接法計算手段７は，前記のようにして作成された剛性行列と外力ベクトルに基づき，ν個のコース空間の方程式を直接法で解く（ステップ７）。この処理は，前記アルゴリズムの（６．２６）式を解くことに相当する。同様に，反復法計算手段８は，第ν世代の共役空間の方程式を反復法で解く（ステップ８）。この処理は，前記アルゴリズムの（６．２７）式を解くことに相当する。
【００９１】
解出力手段９は，直接法計算手段７によって得られた解と，反復法計算手段８によって得られた解を加算し，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を求め，これをコンピュータの出力装置に出力する（ステップ９）。この処理は，前記アルゴリズムの（６．２８）式を実行することに相当する。
【００９２】
以上の通り，実施例１によれば，コース空間Ｗに対しては直接法を適用し，共役空間Ｖ′には反復法を適用して，両者の解を足すことで解を求めることができる。その場合に，予め条件設定手段２で受け付けた階層数νによって解析対象領域Ωを階層的に分割することにより，階層数νに応じた数のコース空間を解空間Ｖから分離することができ，その分，直接法により求解する領域が増大し，反復法により求解する領域が少なくなり，求解の速度が向上する。また，予め階層数νを設定しておくため，求解の途中で反復法における収束性のチェックが不要であり，その点からも計算速度の向上が見込まれる。
【実施例２】
【００９３】
本発明の実施例２を図３のフローチャート及び図６のブロック図によって説明する。なお，実施例１と同様な構成並びに作用については，説明を省略する。
【００９４】
実施例２においても，前記実施例１と同様に，図４のブロック図に示したコンピュータのハードウェアを使用し，コンピュータを，問題入力手段１１，条件設定手段１２，第０階層その領域分割手段１３，第ｋ階層領域分割手段１４，共役空間生成手段１５，第ｋ世代剛性行列・外力ベクトル作成手段１６，直接法計算手段１７，反復法計算手段１８及び解出力手段１９として，機能させる。
【００９５】
このような構成を有する実施例２の作用は，図３のフローチャートに示すとおりである。すなわち，問題入力手段１は，ユーザからの問題として，解析対象領域Ωと，そのメッシュ分割を受け付ける（ステップ１１）。条件設定手段２は，問題を解決する場合に要求される並列プロセス数ｍ0，最大階層数νおよび各階層の領域数νmνを含む設定条件を，ユーザから受け付ける（ステップ１２）。
【００９６】
前記第０階層領域分割手段３は，前記のようにして受け付けた解析対象領域Ωに対して，Ｖ(0)＝Ｖ，Ｋ(0)＝Ｋ，ｆ(0)＝ｆ，Ｐ(0)＝Ｐ，Ｗ(1)＝Ｗ，Ｖ(1)＝Ｖ′とおき，これに基づいて，並列プロセスにわたる第０階層領域分割を生成する（ステップ１３）。以上のステップ１１から１３は，実施例１のステップ１から３と同じ処理である。
【００９７】
第ｋ階層領域分割手段１４は，第ｋ階層の領域分割を生成する。この場合，第１階層は第０階層の領域のすべてを含み，各階層は直前の階層の領域をすべて含む。（ステップ１４）。共役空間生成手段１５は，前記第ｋ階階層領域分割手段１４によって生成された領域分割をもとに，第ｋ世代コース空間Ｗ(k)を作る。この場合，直前世代の共役空間はＶ(k-1)である。Ｖ(k-1)をコース空間とそのｋ(k-1)共役直交空間に分割する（ステップ１５）。
【００９８】
全世代剛性行列・外力ベクトル作成手段１６は，第ｋ世代の２つの剛性行列と外力ベクトルを作成する（ステップ１６）。直接法計算手段１７は，前記のようにして作成された剛性行列と外力ベクトルに基づき，第ｋ世代コース空間の方程式を直接法で解く（ステップ１７）。同様に，反復法計算手段１８は，第ｋ世代の共役空間の方程式を反復法（前処理付きＣＧ法）で解き，その収束性をチェックする（ステップ１８）。反復法計算手段１８による収束性のチェックの結果，十分な収束性が得られていれば，次の解出力ステップ２０に進み，収束性が得られなければｋをｋ＋１として前記ステップ１４に戻る。また，ｋ＝νの場合（予め設定した階層数νに分割数が達した場合）は反復を完結させた上で，ステップ２０へ進む（ステップ１９）。
【００９９】
解出力手段１９は，直接法計算手段１７によって得られた解と，反復法計算手段１８によって得られた解を加算し，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を求め，これをコンピュータの出力装置に出力する（ステップ２０）。
【０１００】
この実施例２においては，階層的に領域分割を行う都度，共役空間に対して反復法による求解を行い，その収束性をチェックするので，比較的早めの段階で共役空間について解が求まる場合にはそれ以上の分割処理が不要となり，計算速度の向上が可能となる。また，問題の規模によっては，予め適切な階層数νを決定（予想）できない場合があるが，実施例２では適切な階層数が収束性のチェックによって定まるので，不要な階層的分割を行うこともない。
【産業上の利用可能性】
【０１０１】
本発明は，大規模な問題を有効に解く有限要素法ソルバアルゴリズムであり，有限要素法に基づく解析一般に適用可能である。たとえば，有限要素法によって定式化される連続体力学全般，固体力学，構造力学，流体力学，熱伝導問題，電磁場問題に利用できる。また，有限要素法は偏微分方程式を支配方程式とする境界値問題の一般的解法とみなされるので，その意味で，ＣＧＣＧ法は，偏微分方程式を支配方程式とする境界値問題全般の解法として適用可能である。
【符号の説明】
【０１０２】
１…問題入力手段
２…条件設定手段
３…前記第０階層領域分割手段
４…各階層領域分割手段
５…共役空間生成手段
６…全世代剛性行列・外力ベクトル作成手段
７…直接法計算手段
８…反復法計算手段
９…解出力手段
１１…問題入力手段
１２…条件設定手段
１３…第０階層その領域分割手段
１４…第ｋ階層領域分割手段
１５…共役空間生成手段
１６…第ｋ世代剛性行列・外力ベクトル作成手段
１７…直接法計算手段
１８…反復法計算手段
１９…解出力手段
【技術分野】
【０００１】
本発明は，自由度１００万以上の大規模構造問題を効果的に解く，大規模並列有限要素法ソルバアルゴリズムに関するもので，特に，並列ＣＧアルゴリズムを基盤とし，領域分割（Domain Decomposition）に基づくアルゴリズムである。
【０００２】
本明細書において，該アルゴリズムを「ＣＧＣＧ法」（Coarse Grid CG method,CGCG method）と称する。Ｋを当該構造問題の剛性行列とし，領域分割による部分領域単位の自由度の空間をコース（coarse）空間，それに共役直交な空間をファイン空間と称する。ＣＧＣＧ法は，コース空間とファイン空間を分離し，それぞれの解を重ね合わせるものである。この際，領域分割は，コース空間の生成及びファイン空間上のＣＧ法の並列化に用いられる。
【背景技術】
【０００３】
このようなＣＧＣＧ法を利用した並列有限要素計算システムとして，特許文献１に記載のものが知られている。この発明は，本出願人が提案したものである。この先行発明は，ＣＧＣＧ法を使用することで，自由度が１００万以上の大規模構造問題を解く際，解が発散することなく求まり，反復計算の回数が少なく，短時間の計算時間で解が求まる，大規模構造問題を解くシステムを提供するものである。
【０００４】
このＣＧＣＧ法では，解析対象の領域分割を行ない，それに基づいて部分領域の内部と境界の区別をすることなしにコース空間を設定し，前処理として簡便な対角スケーリングを採用した前処理つき射影ＣＧ法を用いる。すなわち，並列処理用に分割された領域にわたる“コース運動”を考え，解空間をコース空間とその共役直交空間に分解する。コース空間は直接法で解き，共役直交空間には反復法（共役射影並列ＣＧ法）を用いる。コース運動として剛体運動をとることで，従来の解法に比較して高速に解を得ることができる。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】特許第３８２４５８２号公報
【０００６】
【非特許文献１】J. Mandel, "Balancing Domain Decomposition". Communications on Numerical Methods in Engineering, 1993, 9, pp.233-241
【発明の概要】
【発明が解決しようとする課題】
【０００７】
ＣＧＣＧ法は大規模な有限要素法問題を解くための非常に高速なアルゴリズムであるが，薄板構造への適用性に困難があり，また，よりローカルな問題が全体の収束性に与える場合など，問題依存による収束性の向上の必要がある場合がある。本発明アルゴリズムは，直接法を適用する小さなコース空間と反復法を適用する大きなファイン空間を階層的に多段に取っていくものであるが，コース空間を多段にすることによって，直接法の負荷が大きくなり，薄板構造がより解きやすい方向に改善される。また，構造のよりローカルな運動が表現可能になり，構造依存による収束性の向上がはかれる。
【０００８】
また，ＣＧＣＧ法は，反復法を基本にするため，固有値問題，複数荷重問題など，剛性行列は同一で，荷重条件など，右辺のみが異なる問題などが不得意であるが，直接法の負荷を大きくすることによってこれらがより解きやすい方向に改善される。
【０００９】
本発明は前記のような要望に基づいて提案されたものであって，その目的は，ＣＧＣＧ法を発展させることで，より短時間で解が求まる並列有限要素計算システムを提供することにある。
【課題を解決するための手段】
【００１０】
前記の目的を達成するために，本発明の並列有限要素計算システムは，次の構成を有することを特徴とする。
(1) 解析対象領域Ωと，そのメッシュ分割を与える。
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する。
(3) 並列プロセスにわたる第０階層領域分割を生成する。
(4) ν個の階層の領域分割を生成する。
【００１１】
(5) 全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る。
(6) 全世代の２組の剛性行列と２組の外力ベクトルを作る。
(7) ν個のコース空間の方程式を直接法で解く。
(8) 第ν世代の共役空間の方程式を反復法(前処理付きＣＧ法)で解く。
(9) ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える。
【００１２】
また，次の構成も，本発明の一態様である。
(1) 解析対象領域Ωと，そのメッシュ分割を与える。
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する。
(3) 並列プロセスにわたる第０階層領域分割を生成する。
(4) 第ｋ階層領域分割をもとに第ｋ世代コース空間Ｗ(k)を作る。
【００１３】
(5) 直前世代の共役空間Ｖ(k-1)をコース空間とそのＫ(k-1)共役直交空間に分割する。
(6) 第ｋ世代の２つの剛性行列と外力ベクトルを作る。
(7) 第ｋ世代コース空間の方程式を直接法で解く。
(8) 第ｋ世代の共役空間の方程式について，反復法の収束性をチェックする。
(9) 十分な収束性のチェックを行う。収束性が得られなければ，ｋをｋ＋１として(5)に戻る。
(10) 十分な収束性が得られた場合には，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を与える。
【発明の効果】
【００１４】
本発明によれば，直接法を適用する小さなコース空間と反復法を適用する大きなファイン空間を階層的に多段に取っていくものであるが，コース空間を多段に取ることによって，直接法の負荷が大きくなり，薄板構造がより解きやすい方向に改善される。また，構造のよりローカルな運動が表現可能になり，構造依存による収束性の向上がはかれる。また，固有値問題，複数荷重問題など剛性行列は同一で，荷重条件など，右辺のみが異なる問題が，直接法の負荷を大きくすることによってこれらがより解きやすい方向に改善される。
【図面の簡単な説明】
【００１５】
【図１】本発明の領域分割のイメージを示す模式図。
【図２】本発明の実施例１の処理を示すフローチャート。
【図３】本発明の実施例２の処理を示すフローチャート。
【図４】本発明のシステムを実現するコンピュータのハードウェアを示すブロック図。
【図５】本発明の実施例１の各手段を示すブロック図。
【図６】本発明の実施例２のの各手段を示すブロック図。
【発明を実施するための形態】
【００１６】
本発明の前提となる事項について，説明する。
（１）概要
本発明において，Ωを解析対象領域とする。この上の構造問題を有限要素法にしたがって離散化する。目標はここから帰着される線形方程式
Ｋｕ＝ｆ……（１．１）
を解くことである。
【００１７】
(共役直交射影子)
ベクトル空間Ｖ上の正定値対称変換ＡとＶ上の射影Ｐが与えられたとき，Ｖの部分空間ＰＶと（１−Ｐ）Ｖに直和分解され，つぎの２式がなりたつ。
Ｐ2=Ｐ……(５．４)
ＡＰ=ＰTＡ……(５．５)
あるいはつぎのようにもいえる。ベクトル空間Ｖ上の正定値対称変換ＡとＶの部分空間Ｗが与えられたとき，ＶからＷへの射影Ｐは(５．４)式をみたし，ＡとＰの間には関係(５．５)がなりたつ。
逆に，(５．４)式と(５．５)式をみたすＶ上の線形変換Ｐと，正定値対称変換Ａが与えられたとき，Ｖは部分空間ＰＶと（１−Ｐ）Ｖに直和分解され，これら２つの部分空間はＡに関して共役直交である。
【００１８】
（部分空間の共役直交基底）
ＷをＶのｒ次元部分空間とするとき，Ｗの基底ｆ1，…，ｆrを互いにＡ共役直交であるようにえらべる。すなわち，Ｗの基底ｆ1，…，ｆrを適当にとることにより，つぎの関係がなりたつようにできる。
(ｆj，Ａｆi)=０ (ｉ≠ｊ)……(５．６)
【００１９】
（ベクトル空間の共役直交直和分解）
ＷをＶの，Ｖに一致しない部分空間とするとき，Ｖは，ＷとＡに関して共役直交な部分空間Ｖ′に直和分解できる。
【数１】
【００２０】
(共役直交射影の表現)
Ｖをベクトル空間，Ｗをｒ次元部分空間，ｆ1，…，ｆrをＷの共役直交基底，ＡをＶ上の正定値対称変換，Ｐを共役直交射影子とするとき，Ｐは
【数２】
で与えられる。
【００２１】
（２）ＣＧＣＧ法
（領域分割）
Ωは有限要素分割されている。その分割に基づいてΩを，オーバラップしないが全部の和集合はΩに一致するような部分領域に分ける。その第１の分割は並列プロセスにわたる分割である。これをさらに分割して，ΩI，I=１，…，ｍとする。ΩIのある面または稜線は別の部分領域と接していて，これらの面または稜線は複数の部分領域に共有される。Ω内部での部分領域どうしの境界面を内部境界と呼ぶ。この面はあるΩIのある面のグループは外部境界面(Ωの境界面)の一部であることもある。
【００２２】
以下では２つの空間ＶIおよび~ＶIを定義するが，これは，以下で述べる部分空間から全空間Ｖへの延長とその逆の制限を定義し，これをもとにコースグリッドの項で述べた方法でコース空間を定義するのに用いるためである。
【００２３】
Ω上で定義された有限要素基底関数を｛ψj｝とする。ψjをΩIに制限した関数のうち，恒等的にゼロでないものを~ψIｊと書く。これらの｛~ψIj｝はΩIの有限要素基底関数をなす。｛~ψIj｝で離散化されたΩI上の変位場の空間を~ＶIと書き，切断された局所空間空間と呼ぶ。
【００２４】
一方，ΩI(およびその境界)に含まれる節点上に定義される許容変位場がなすベクトル空間をＶIと書く。ΩIの基底関数~ψjはもとの基底関数ψjの外側を切り取った(ただし~ψjの定義域には境界面も含む)ものなので，~ＶIはＶIとは異なる。ＶIはＶの部分空間だが，~ＶIは部分空間ではないことに注意せよ。ΩIはオーバラップしないが，ＶIは内部境界のまわりでオーバラップする。
【００２５】
~ＶIの各節点における基底関数に，その節点に対応するＶIの(実はＶの)基底関数の対応させる写像は１対１である。ＶIはＶの部分空間なので，この写像
【数３】
を~ＶIからＶへの，あるいは~ＶIからＶIへの延長(または延長写像)(prolongation)と呼ぶ。~ＶIはＶの部分空間ではないが，ＮIによって部分空間に埋め込まれる。つまり，ＮI~ＶIは部分空間である。一方，ＮIの転置
【数４】
をＶから~ＶIへの制限(または制限写像)(restriction)と呼ぶ。
【００２６】
グローバルな節点を任意にとる。その節点は，あるIに対する部分領域ΩI上にあるとする。対応する自由度を制限(ＮI)Tによってローカルな自由度に写し，その後，延長ＮIによってグローバルな自由度に写す。内部境界上の節点の自由度は，(ＮI)Tによって~ＶIを含むいくつかの切断された局所空間のローカルな自由度に写され，延長によって元の節点を共有する複数の部分領域に対応する自由度にもどる。この重複を解消するため，つぎのような対角行列ＤIを定義する。領域分割の際，内部境界上の重複する節点の数を求めておき，その自由度に対応する成分を(１/重複節点数)，重複のない成分を１とする。この作用は
【数５】
と表わせる。
【００２７】
Ｖは部分空間ＶI全部の和からなり，ＶIはＮIによる~ＶIの延長で得られるが，重複を考慮すると，これらの関係は
【数６】
と表わせる。
【００２８】
（コースグリッド）
部分領域ΩI，Ｉ=１，…，ｍの剛体運動を考える。要素分割したΩの全節点数をｎとし，グローバルな節点番号をｉ(ｉ=１，…，ｎ)で表わす。初期の節点座標を
【数７】
とする。肩に付けた添え字Ｉは，任意の節点はいずれかの部分領域ΩIにあるので，これを明示するために付けたものである。
【００２９】
剛体変位後の座標を
【数８】
とする。
【００３０】
剛体変位は平行移動３自由度，回転３自由度をもち，
【数９】
と表わせる。
【００３１】
ただし
【数１０】
であり，Ｕ1，Ｕ2，Ｕ3は各軸に沿った平行移動，
【数１１】
は各軸の周りの，ベクトル
【数１２】
の右ねじ方向の周りの，大きさ｜θ｜の回転を表わす。
【００３２】
回転は微小とすると，(６．３)は
【数１３】
と書ける。したがって，対応する節点変位をｕIｉとすると
【数１４】
となる。
【００３３】
以上をもとに，部分領域ΩIの剛体変位場(節点座標を変数にとり，剛体変位を与える~ＶIのベクトル)の基底をつぎのように定義する。ただし，ΩIの一部が全体領域の外表面に接しており，この部分に変位拘束が課されている場合を除く。
【数１５】
【００３４】
このように，剛体変位場ベクトルの自由度は１節点あたり６，つまり，部分領域ΩIの剛体変位を与える空間の基底(６．４)の次元は６である。上で除外した，ΩIの表面が全体空間の外表面に一部をなし，ここに変位拘束が課されている場合には，その次元が５以下になることもありうるが，独立な基底をもとにこれを拡大して６にそろえる。
【００３５】
このようにして得られた空間を~ＷIと書くことにする。ΩIの剛体変位はもとよりΩIの許容変位場なので，ベクトル(６．４)はΩIの有限要素基底関数｛~ψIj｝で表現できる。すなわち，~ＷIは~ＶIの部分空間である。ＮIによって~ＶIがＶIに延長されたように，~ＷIもＮIによってＶIに延長される。
【００３６】
ただし，自由度の重複の解消は行わなければならない。(６．２)と同じように，内部境界自由度の重複を，１の分割ＤIを用いて解消し，
【数１６】
とおく。これを剛体運動によるＶのコース空間と呼ぶ。コース空間ＷはＶの部分空間である。
【００３７】
（グローバル空間の共役直交分解）
解析対象モデルΩはメッシュ分割によりｎ個の節点をもち，したがって自由度は３ｎである。グローバルな解空間Ｖの次元および(１．１)の左辺のグローバルな剛性行列Ｋの次元も３ｎである。
【００３８】
コースグリッドの項で説明した通り，コース空間Ｗとその共役空間Ｖ′に分解する。これは(５．７)式に対応する。Ｗの基底は(６．４)式を部分領域の数ｍ個分を並べたものであり，したがってＷの次元は６ｍである。Ｖ′は，(５．７)式のように作られ，ＷとＶ′はＫに関して共役直交である。Ｖ′の次元は３ｎ−６ｍであり，これをあらためて(ｆ1，…，ｆ6m)とする。６ｍは３ｎに比べて，直接法が適用できるほど十分小さいものとする。Ｖの基底は，Ｗの基底に，これと共役直交なＶ′の基底を並べて作られる。
【００３９】
Ｖ′の基底を(ｇ6m+1，…，ｇ3n)とする。Ｋの表現はこの基底の並び方に合わせ，(４．１)のように作られているものとする。つまり，基底の並べ替えに応じてＫはつぎのような表現をもつ。
【数１７】
【００４０】
射影Ｐの表現は，（５．８）式により６ｍ次元単位行列である。より正確には，Ｐは，次元を合わせて
【数１８】
で定義される。
【００４１】
Ｐは６ｍ×３ｎの横長の行列であり，３ｎ次元の解空間Ｖを６ｍ次元のコース空間Ｗに射影する。
Ｗ=ＰＶ
【００４２】
ここで，Ｐはべき等行列である。
Ｐの転置行列
【数１９】
は，ＷをＶの中のＷと同型な部分空間に射影する３ｎ×６ｍの縦長の行列である。
【数２０】
【００４３】
コース空間Ｗ上のＫの作用は，(４．１５)と同じように
~Ｋ=ＰＫＰT
と書ける。
【００４４】
共役空間Ｖ′への射影は１−Ｐだが，(６．７)と同様，次元を合わせて書くと，
【数２１】
となる。
【００４５】
このように，１−Ｐは(３ｎ−６ｍ)×３ｎの横長の行列であり，３ｎ次元の解空間Ｖを３ｎ−６ｍ次元の共役空間Ｖ′に射影する。１−Ｐの転置行列
【数２２】
は，Ｖ′をＶの中のＶ′と同型な部分空間に射影する３ｎ×(３ｎ−６ｍ)の縦長の行列である。すなわち，
【数２３】
【００４６】
共役空間Ｖ′上のＫの作用は
Ｋ′=(１−Ｐ)Ｋ(１−ＰT)
と書け，これは３ｎ−６ｍ次元の正定値対称行列である。
【００４７】
（５．６）式により，Ｋは~ＫとＫ′により，
【数２４】
と書ける。
【００４８】
（ＣＧＣＧ法）
前記の通り，解空間Ｖをコース空間Ｗとその共役空間とに直和分解する。
【数２５】
【００４９】
この上でグローバルな線形方程式，すなわち
Ｋｖ=ｆ……(６．１３)
を解く。
【００５０】
まず，(６．１３)をＷおよびＶ′にそれぞれＰおよび１−Ｐによって部分空間に写す。
【数２６】
【００５１】
この上で，ｖ∈ＶをＰTおよび１−ＰTによってそれぞれＷおよびＶ′の元に限定すると，
【数２７】
が得られる。
【００５２】
これを(６．１４)に入れると
【数２８】
が得られる。ここで
【数２９】
とおく。
【００５３】
一方，ＰＶ=Ｗの双対空間Ｗ*はＷと同一視される。ｆはＫの像空間ＫＶに含まれるので，ＰｆはＰＫＰTによるＷの像空間(ＰＫＰT)Ｗに含まれる。同様に，(１−Ｐ)Ｖ=Ｖ′の双対空間Ｖ′*はＶ′と同一視され，(１−Ｐ)ｆは((１−Ｐ)Ｋ(１−ＰT))Ｖ=Ｖ′に含まれる。そこで，
【数３０】
とおくと，(６．１４)は
【数３１】
に帰着する。
【００５４】
すでに述べたとおり，(６．１９)は直接法で解く。(６．２０)は，これと同等な方程式
ＧＫ′ｕ=Ｇｆ′……(６．２１)
について，つぎの前処理付きＣＧ法で解く。ここで，ＧはＶ′上の適当な正定値対称行列である。
【数３２】
【００５５】
ＷとＶ′の共役直交性により，数値誤差を考えなければ，残差ｒkは反復過程でつねにＶ′(正確にいえばＶ′の双対空間)にとどまり，その結果，勾配ベクトルｐk，解の近似ベクトルｕkもＶ′(こちらは双対空間ではない)にとどまる。収束解は剛体運動と重ねあわされ，幅広い運動スペクトルを表現するものと期待できる。
【００５６】
（３）本発明のアルゴリズム
本発明においては，以上述べたＣＧＣＧ法のアルゴリズムを拡張する。これをここでは，本発明のアルゴリズムを，仮に，強化版ＣＧＣＧ法(enhanced CGCG method)と呼んでおく。領域分割はν(ν≧１)階層に分割される。ν＝１の場合が，前節までに述べた通常のＣＧＣＧ法であり，部分領域の数はｍである。
【００５７】
(0) 解析対象領域Ωの並列プロセスにわたるｍ０個の部分領域Ω(0,1)，…，Ω(0,m0)。この分割を親領域分割(parent domain decomposition)と呼ぶ。
【００５８】
(1) 部分領域Ω(0,1)，…，Ω(0,m0)を含む，より細かい領域分割Ω(1,1)，…，Ω(1,m1)。これを第１世代領域分割(first generation domain decomposition)と呼ぶ。これがＣＧＣＧ法の説明に用いた部分領域Ω1，…，Ωmにあたり，第１世代部分領域の数はｍ1＝ｍである。分割Ω(1,1)，…，Ω(1,m1)は分割Ω(0,1)，…，Ω(0,m0)を含み，Ω(1,1)，…，Ω(1,m1)の全体が全領域Ωに一致する。したがってｍ1≧ｍ0である。ｍ1＝ｍ0の場合はＣＧＣＧ法に一致しするので，強化版では意味ある分割，すなわちｍ1＞ｍ0とする必要がある。
【００５９】
(2) 部分領域Ω(1,1),…,Ω(1,m1)を含む，より細かい領域分割Ω(2,1),…,Ω(2,m2)。ただし，Ω(2,1),…,Ω(2,m2)のうちの最初のある領域のグループはΩ(1,1)に含まれ，つぎのグループはΩ(1,2)に含まれ，最後のグループはΩ(1,m1)に含まれる。これを第２世代領域分割(second generation domain decomposition)と呼ぶ。
【００６０】
(ν) 部分領域Ω(ν-1,1),…,Ω(ν-1,mν-1)を含む，より細かい領域分割Ω(ν,1),…,Ω(ν,mν)。
ただし，Ω(ν,1),…,Ω(ν,mν)のうちの最初のある領域のグループはΩ(ν-1,1)に含まれ，つぎのグループはΩ(ν-1,2)に含まれ，最後のグループはΩ(ν-1,mν-1)に含まれる。これを第ν世代領域分割(ν-th generation domain decomposition) と呼ぶ。
【００６１】
なお、前記(ν)のΩ(ν-1,mν-1)及びΩ(ν,mν)の標記中、上付き文字のmν-1及びmνは、正確には、ｍν-1及びｍνである。同様に前記(0)から(2)や本文以降等の他の式においても、上付き文字中の下付き文字については、便宜的に下付き文字表示を省略している。
【００６２】
上の階層的領域分割に従って，コース空間をつぎのように階層的に作る。
【００６３】
(1) 第０世代を含む第１世代のコース空間を，コースグリッドの項で述べた方法で構成する。これをＷ(1)とする。
【００６４】
(2) Ω(1,1),…, Ω(1,m1)に対応するそれぞれのコース空間を，コースグリッドの項で述べた方法で各グループ単位の領域分割でコース空間を作る。ただし，各第２世代部分領域は第１世代部分領域の内部境界上の節点は含まず，内部境界面で囲まれている節点のみを用る。
これを第２世代のコース空間と呼び，Ｗ(2)と書く。Ｗ(2)はｍ1個のコース空間からなる，互いに他から切り離されたローカルな空間である。ｍ1個の各コース空間は，そこにグループ化されただけの数の第２世代部分領域を持ち，第２世代部分領域数は全体でｍ2なので，Ｗ(2)の次元は６ｍ2である。第１世代のコース運動の中には，Ｗ(2)上のコース運動では表現できないものがある。
【００６５】
(ν) Ω(ν-1,1),…,Ω(ν-1,mν-1)に対応するそれぞれのコース空間を，コースグリッドの項で述べた方法で各グループ単位の領域分割でコース空間を作る。ただし，各第ν世代部分領域は第ν−１世代部分領域の内部境界上の節点は含まず，内部境界面で囲まれている節点のみを用る。これを第ν世代のコース空間と呼び，Ｗ(ν)と書く。Ｗ(ν)はｍν-1個のコース空間からなる，互いに他から切り離されたローカルな空間である。ｍν-1個の各コース空間は，そこにグループ化されただけ数の第ν世代部分領域を持ち，27第ν世代部分領域数は全体でｍνなので，Ｗ(ν)の次元は６ｍνである。第ν−１世代のコース運動の中には，Ｗ(ν)上のコース運動では表現できないものがある。
【００６６】
本発明の強化版ＣＧＣＧ法がマルチグリッド法と異なる点は，マルチグリッド法が，細かいメッシュ分割をもとに，より粗い分割を作っていくのに対して，強化版ＣＧＣＧ法では，コースグリッドをより細かく作っていくことにある。すなわち，剛性方程式(６．１３)を解く際には，前処理行列Ｇとしては負荷の軽い対角スケーリングなどを使うのが効率的である。そのため，コース運動を解く際に，それが粗すぎると，十分な精度が得られないことが経験上わかっている。一方，非特許文献１に記載のＢＤＤ法では，分割領域ΩＩの境界面のみに剛体運動を与え，ΩＩ内部は，剛体運動で決まった境界面の変位拘束条件をもとに，部分構造としての剛性方程式を解いているため，ＢＤＤ法のコースグリッドにはＣＧＣＧ法のそれよりもより詳細な表現が与えられている。この意味で，ＣＧＣＧ法のコース運動はＢＤＤ法と比べ粗い。コースグリッドを細かく取っていくというアイデアは，この粗さを補完することに対応している。
【００６７】
Ｖは第１世代子分割に基づいて，共役直交分解される。これは前節のＣＧＣＧ法の手順と同じものである。以下で記号を改める。前節を参照して，
【数３３】
とおく。
【００６８】
第１世代子アルゴリズム(通常のＣＧＣＧ法アルゴリズム)はつぎのように書ける。
【数３４】
【００６９】
前記（６．２２）式における下から３行目のｗ(1)に関する方程式は直接法で解く。また，Ｖ(1)上の剛性方程式は第２世代分割に基づいてコース空間Ｗ(2)とそのＫ(1)共役直交空間に分解され，この上にＣＧＣＧ法が適用できる。
【数３５】
【００７０】
Ｖ(1)はＷ(1)と共役直交なので，もとの空間Ｖ＝Ｖ(0)は，
【数３６】
と書ける。
【００７１】
第２世代アルゴリズムはつぎのようになる。
【数３７】
【００７２】
前記（６．２４）式における下から４行目のｗ(2)に関する方程式は直接法で解く。
【００７３】
以下同様にして，第ν−１世代アルゴリズムが与えられたとして，第ν世代アルゴリズムはつぎのようになる。
【数３８】
【００７４】
前記（６．２５）式の下から４行目のｗ(ν)に関する方程式は直接法で解く。
【００７５】
第１世代アルゴリズム(６．２２)，第２世代アルゴリズム(６．２４)から第ν−１世代アルゴリズムに現れるｕ(1)，ｕ(2)，…，ｕ(ν-1)についての方程式は実際には解く必要はなく，それぞれ，つぎの世代のアルゴリズムに吸収される。結局，強化版ＣＧＣＧ法アルゴリズムはつぎのようにまとめられる。ν個のコース問題
【数３９】
をそれぞれ直接法で解き，最後の第ν世代アルゴリズムにおける３nν次元の方程式
Ｋνｕν＝ｆ′ν， ｕν∈Ｖν……（６．２７）
に前処理付きＣＧ法を適用して，これらで得られた解をすべて足しあわせてグローバルな解を
ｕ＝ｕ(0)＝ｗ(1)+ｗ(2)+…+ｗ(ν)+ｕ(ν)……(６．２８)
で得る。
【００７６】
Ｖ(ν)の次元３nνはつぎで与えられる。
３nν＝３n0−６m1−６m2−…−６mν……(６．２９)
【００７７】
また，解空間の共役直交分解による階層的な関係はつぎのように書ける。
【数４０】
【００７８】
本発明の強化版ＣＧＣＧ法においては，コース空間を多段にとることで，より細かいローカルな運動の表現が可能になって，ＣＧＣＧ法に比べて収束性が向上することが期待できる。また，階層を深くとれば直接法の負荷が大きくなるため，反復法では解きにくい問題も，より解きやすい方向に改善されるものと期待される。
【実施例１】
【００７９】
以下，本発明の実施例１を図面に従って説明する。
まず，図４を参照して，実施例１の並列有限要素計算システムを構成するコンピュータのハードウェアを説明する。このハードウェアは，バス３０を介して相互に接続されるＣＰＵ３１と，ＲＡＭ３２と，記憶装置３３と，入力装置３４と，出力装置３５とを具備する。記憶装置３３はハードディスクやメモリ等の外部記憶装置である。入力装置３４は，キーボードやマウス等のユーザによって操作されることで，各種データをＣＰＵ３１や記憶装置３３に出力する。出力装置３５は，モニタやプリンタに例示され，ＣＰＵ３１から出力される求解の結果をユーザに対し視認可能に出力する。
【００８０】
ＣＰＵ３１は，入力装置３４からの入力に応答して，記憶装置３３内の強化版ＣＧＣＧ法のアルゴリズムを利用した求解プログラムを実行する。この際，記憶装置３３からの各種データやプログラムはＲＡＭ３２に一時格納され，ＣＰＵ３１はＲＡＭ３２内のデータを用いて各種処理を実行する。
【００８１】
前記強化版ＣＧＣＧ法を実現するために，実施例１の並列有限要素計算システムにおいては，前記のようなハードウェアを有するコンピュータを，問題入力手段１，条件設定手段２，第０階層領域分割手段３，各階層領域分割手段４，共役空間生成手段５，全世代剛性行列・外力ベクトル作成手段６，直接法計算手段７，反復法計算手段８及び解出力手段９として，機能させる。
【００８２】
すなわち，コンピュータの有するキーボード，マウス，ネットワーク，記憶媒体などの入力装置３４から入力される問題を受け付ける問題入力手段１，同じく設定すべき条件を受け付ける条件設定手段２を備える。
【００８３】
前記第０階層領域分割手段３，各階層領域分割手段４，共役空間生成手段５，全世代剛性行列・外力ベクトル作成手段６，直接法計算手段７，反復法計算手段８及び解出力手段９は，コンピュータのメモリ（ＲＡＭ３２）やＣＰＵ３１上に展開されるプログラムによって，コンピュータ上に構築される手段で，これら各手段が前記強化版ＣＧＣＧ法のアルゴリズムに従って，与えられた問題に対する解を求める。
【００８４】
解出力手段９は，コンピュータに設けられたディスプレイ，プリンタ，記憶媒体，或いはネットワークを通じて接続された他のコンピュータなどの出力装置３５に接続され，これらのハードウェア上に求められた解を出力する。
【００８５】
次に，各手段が実行する処理を図２のフローチャートに従って説明する。
問題入力手段１は，ユーザからの問題として，解析対象領域Ωと，そのメッシュ分割を受け付ける（ステップ１）。条件設定手段２は，問題を解決する場合に要求される並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを含む設定条件を，ユーザから受け付ける（ステップ２）。この処理は，前記アルゴリズムの階層構造の深さを決めることに相当する。
【００８６】
前記第０階層領域分割手段３は，前記のようにして受け付けた解析対象領域Ωに対して，Ｖ(0)＝Ｖ，Ｋ(0)＝Ｋ，ｆ(0)＝ｆ，Ｐ(0)＝Ｐ，Ｗ(1)＝Ｗ，Ｖ(1)＝Ｖ′とおき，これに基づいて，並列プロセスにわたる第０階層領域分割を生成する（ステップ３）。
【００８７】
各階層領域分割手段４は，設定条件として受け付けた階層数νに基づいてν個の階層の領域分割を生成する（ステップ４）。この処理は，前記アルゴリズムの（６．２５）式のＷ(1)，…，Ｗ(ν)を決めることに相当する。
【００８８】
共役空間生成手段５は，前記各階層領域分割手段４によって生成された全ν階層の領域分割をもとに，第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作成する（ステップ５）。この処理は，前記アルゴリズムの（６．２５）式の第２式を作成することに相当する。
【００８９】
全世代剛性行列・外力ベクトル作成手段６は，全世代の２組の剛性行列と２組の外力ベクトルを作成する（ステップ６）。この処理は，前記アルゴリズムの（６．２６）式を作成することに相当する。
【００９０】
直接法計算手段７は，前記のようにして作成された剛性行列と外力ベクトルに基づき，ν個のコース空間の方程式を直接法で解く（ステップ７）。この処理は，前記アルゴリズムの（６．２６）式を解くことに相当する。同様に，反復法計算手段８は，第ν世代の共役空間の方程式を反復法で解く（ステップ８）。この処理は，前記アルゴリズムの（６．２７）式を解くことに相当する。
【００９１】
解出力手段９は，直接法計算手段７によって得られた解と，反復法計算手段８によって得られた解を加算し，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を求め，これをコンピュータの出力装置に出力する（ステップ９）。この処理は，前記アルゴリズムの（６．２８）式を実行することに相当する。
【００９２】
以上の通り，実施例１によれば，コース空間Ｗに対しては直接法を適用し，共役空間Ｖ′には反復法を適用して，両者の解を足すことで解を求めることができる。その場合に，予め条件設定手段２で受け付けた階層数νによって解析対象領域Ωを階層的に分割することにより，階層数νに応じた数のコース空間を解空間Ｖから分離することができ，その分，直接法により求解する領域が増大し，反復法により求解する領域が少なくなり，求解の速度が向上する。また，予め階層数νを設定しておくため，求解の途中で反復法における収束性のチェックが不要であり，その点からも計算速度の向上が見込まれる。
【実施例２】
【００９３】
本発明の実施例２を図３のフローチャート及び図６のブロック図によって説明する。なお，実施例１と同様な構成並びに作用については，説明を省略する。
【００９４】
実施例２においても，前記実施例１と同様に，図４のブロック図に示したコンピュータのハードウェアを使用し，コンピュータを，問題入力手段１１，条件設定手段１２，第０階層その領域分割手段１３，第ｋ階層領域分割手段１４，共役空間生成手段１５，第ｋ世代剛性行列・外力ベクトル作成手段１６，直接法計算手段１７，反復法計算手段１８及び解出力手段１９として，機能させる。
【００９５】
このような構成を有する実施例２の作用は，図３のフローチャートに示すとおりである。すなわち，問題入力手段１は，ユーザからの問題として，解析対象領域Ωと，そのメッシュ分割を受け付ける（ステップ１１）。条件設定手段２は，問題を解決する場合に要求される並列プロセス数ｍ0，最大階層数νおよび各階層の領域数νmνを含む設定条件を，ユーザから受け付ける（ステップ１２）。
【００９６】
前記第０階層領域分割手段３は，前記のようにして受け付けた解析対象領域Ωに対して，Ｖ(0)＝Ｖ，Ｋ(0)＝Ｋ，ｆ(0)＝ｆ，Ｐ(0)＝Ｐ，Ｗ(1)＝Ｗ，Ｖ(1)＝Ｖ′とおき，これに基づいて，並列プロセスにわたる第０階層領域分割を生成する（ステップ１３）。以上のステップ１１から１３は，実施例１のステップ１から３と同じ処理である。
【００９７】
第ｋ階層領域分割手段１４は，第ｋ階層の領域分割を生成する。この場合，第１階層は第０階層の領域のすべてを含み，各階層は直前の階層の領域をすべて含む。（ステップ１４）。共役空間生成手段１５は，前記第ｋ階階層領域分割手段１４によって生成された領域分割をもとに，第ｋ世代コース空間Ｗ(k)を作る。この場合，直前世代の共役空間はＶ(k-1)である。Ｖ(k-1)をコース空間とそのｋ(k-1)共役直交空間に分割する（ステップ１５）。
【００９８】
全世代剛性行列・外力ベクトル作成手段１６は，第ｋ世代の２つの剛性行列と外力ベクトルを作成する（ステップ１６）。直接法計算手段１７は，前記のようにして作成された剛性行列と外力ベクトルに基づき，第ｋ世代コース空間の方程式を直接法で解く（ステップ１７）。同様に，反復法計算手段１８は，第ｋ世代の共役空間の方程式を反復法（前処理付きＣＧ法）で解き，その収束性をチェックする（ステップ１８）。反復法計算手段１８による収束性のチェックの結果，十分な収束性が得られていれば，次の解出力ステップ２０に進み，収束性が得られなければｋをｋ＋１として前記ステップ１４に戻る。また，ｋ＝νの場合（予め設定した階層数νに分割数が達した場合）は反復を完結させた上で，ステップ２０へ進む（ステップ１９）。
【００９９】
解出力手段１９は，直接法計算手段１７によって得られた解と，反復法計算手段１８によって得られた解を加算し，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を求め，これをコンピュータの出力装置に出力する（ステップ２０）。
【０１００】
この実施例２においては，階層的に領域分割を行う都度，共役空間に対して反復法による求解を行い，その収束性をチェックするので，比較的早めの段階で共役空間について解が求まる場合にはそれ以上の分割処理が不要となり，計算速度の向上が可能となる。また，問題の規模によっては，予め適切な階層数νを決定（予想）できない場合があるが，実施例２では適切な階層数が収束性のチェックによって定まるので，不要な階層的分割を行うこともない。
【産業上の利用可能性】
【０１０１】
本発明は，大規模な問題を有効に解く有限要素法ソルバアルゴリズムであり，有限要素法に基づく解析一般に適用可能である。たとえば，有限要素法によって定式化される連続体力学全般，固体力学，構造力学，流体力学，熱伝導問題，電磁場問題に利用できる。また，有限要素法は偏微分方程式を支配方程式とする境界値問題の一般的解法とみなされるので，その意味で，ＣＧＣＧ法は，偏微分方程式を支配方程式とする境界値問題全般の解法として適用可能である。
【符号の説明】
【０１０２】
１…問題入力手段
２…条件設定手段
３…前記第０階層領域分割手段
４…各階層領域分割手段
５…共役空間生成手段
６…全世代剛性行列・外力ベクトル作成手段
７…直接法計算手段
８…反復法計算手段
９…解出力手段
１１…問題入力手段
１２…条件設定手段
１３…第０階層その領域分割手段
１４…第ｋ階層領域分割手段
１５…共役空間生成手段
１６…第ｋ世代剛性行列・外力ベクトル作成手段
１７…直接法計算手段
１８…反復法計算手段
１９…解出力手段
【特許請求の範囲】
【請求項１】
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０階層領域分割手段と，
(4) ν個の階層の領域分割を生成する各階層領域分割手段と，
(5) 全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る共役空間生成手段と，
(6) 全世代の２組の剛性行列と２組の外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) ν個のコース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ν世代の共役空間の方程式を反復法で解く反復法計算手段と，
(9) ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える解出力手段と，
を有することを特徴とする並列有限要素計算システム。
【請求項２】
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０層領域分割手段と，
(4) 第ｋ階層領域分割をもとに第ｋ世代コース空間Ｗ(k)を作る第ｋ世代コース空間生成手段と，
(5) 直前世代の共役空間Ｖ(k-1)をコース空間とそのＫ(k-1)共役直交空間に分割する共役空間分割手段と，
(6) 第ｋ世代の２つの剛性行列と外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) 第ｋ世代コース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ｋ世代の共役空間の方程式について，反復法の収束性をチェックし，十分な収束性収束性が得られなければ，ｋをｋ＋１として，前記共役空間分割手段による処理を再度行わせる収束性判定手段と，
(9) 前記収束性判定手段において十分な収束性が得られた場合には，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を与える解出力手段と，
を有することを特徴とする並列有限要素計算システム。
【請求項３】
コンピュータを，
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０階層領域分割手段と，
(4) ν個の階層の領域分割を生成する各階層領域分割手段と，
(5) 全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る共役空間生成手段と，
(6) 全世代の２組の剛性行列と２組の外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) ν個のコース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ν世代の共役空間の方程式を反復法で解く反復法計算手段と，
(9) ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える解出力手段，
として機能させることを特徴とする並列有限要素計算プログラム。
【請求項４】
コンピュータを，
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０層領域分割手段と，
(4) 第ｋ階層領域分割をもとに第ｋ世代コース空間Ｗ(k)を作る第ｋ世代コース空間生成手段と，
(5) 直前世代の共役空間Ｖ(k-1)をコース空間とそのＫ(k-1)共役直交空間に分割する共役空間分割手段と，
(6) 第ｋ世代の２つの剛性行列と外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) 第ｋ世代コース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ｋ世代の共役空間の方程式について，反復法の収束性をチェックし，十分な収束性収束性が得られなければ，ｋをｋ＋１として，前記共役空間分割手段による処理を再度行わせる収束性判定手段と，
(9) 前記収束性判定手段において十分な収束性が得られた場合には，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を与える解出力手段，
として機能させることを特徴とする並列有限要素計算プログラム。
【請求項１】
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０階層領域分割手段と，
(4) ν個の階層の領域分割を生成する各階層領域分割手段と，
(5) 全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る共役空間生成手段と，
(6) 全世代の２組の剛性行列と２組の外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) ν個のコース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ν世代の共役空間の方程式を反復法で解く反復法計算手段と，
(9) ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える解出力手段と，
を有することを特徴とする並列有限要素計算システム。
【請求項２】
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０層領域分割手段と，
(4) 第ｋ階層領域分割をもとに第ｋ世代コース空間Ｗ(k)を作る第ｋ世代コース空間生成手段と，
(5) 直前世代の共役空間Ｖ(k-1)をコース空間とそのＫ(k-1)共役直交空間に分割する共役空間分割手段と，
(6) 第ｋ世代の２つの剛性行列と外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) 第ｋ世代コース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ｋ世代の共役空間の方程式について，反復法の収束性をチェックし，十分な収束性収束性が得られなければ，ｋをｋ＋１として，前記共役空間分割手段による処理を再度行わせる収束性判定手段と，
(9) 前記収束性判定手段において十分な収束性が得られた場合には，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を与える解出力手段と，
を有することを特徴とする並列有限要素計算システム。
【請求項３】
コンピュータを，
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０階層領域分割手段と，
(4) ν個の階層の領域分割を生成する各階層領域分割手段と，
(5) 全ν階層の領域分割をもとに第１，…，ν世代コース空間Ｗ(1)，…，Ｗ(ν)と，Ｗ(ν)の共役空間Ｖ(ν)を作る共役空間生成手段と，
(6) 全世代の２組の剛性行列と２組の外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) ν個のコース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ν世代の共役空間の方程式を反復法で解く反復法計算手段と，
(9) ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(ν)＋ｕ(ν)によって解を与える解出力手段，
として機能させることを特徴とする並列有限要素計算プログラム。
【請求項４】
コンピュータを，
(1) 解析対象領域Ωと，そのメッシュ分割を与える問題入力手段と，
(2) 並列プロセス数ｍ0，階層数νおよび各階層の領域数νmνを入力する条件設定手段と，
(3) 並列プロセスにわたる第０階層領域分割を生成する第０層領域分割手段と，
(4) 第ｋ階層領域分割をもとに第ｋ世代コース空間Ｗ(k)を作る第ｋ世代コース空間生成手段と，
(5) 直前世代の共役空間Ｖ(k-1)をコース空間とそのＫ(k-1)共役直交空間に分割する共役空間分割手段と，
(6) 第ｋ世代の２つの剛性行列と外力ベクトルを作る全世代剛性行列・外力ベクトル作成手段と，
(7) 第ｋ世代コース空間の方程式を直接法で解く直接法計算手段と，
(8) 第ｋ世代の共役空間の方程式について，反復法の収束性をチェックし，十分な収束性収束性が得られなければ，ｋをｋ＋１として，前記共役空間分割手段による処理を再度行わせる収束性判定手段と，
(9) 前記収束性判定手段において十分な収束性が得られた場合には，ｕ＝ｕ(0)＝ｗ(1)＋…＋ｗ(k)＋ｕ(k)によって解を与える解出力手段，
として機能させることを特徴とする並列有限要素計算プログラム。
【図２】
【図３】
【図４】
【図５】
【図６】
【図１】
【図３】
【図４】
【図５】
【図６】
【図１】
【公開番号】特開２０１１−２４２８１８（Ｐ２０１１−２４２８１８Ａ）
【公開日】平成２３年１２月１日（２０１１．１２．１）
【国際特許分類】
【出願番号】特願２０１０−１０６３２２（Ｐ２０１０−１０６３２２）
【出願日】平成２２年５月６日（２０１０．５．６）
【出願人】（５９９１７２７８２）株式会社アライドエンジニアリング (3)
【Ｆターム（参考）】
【公開日】平成２３年１２月１日（２０１１．１２．１）
【国際特許分類】
【出願日】平成２２年５月６日（２０１０．５．６）
【出願人】（５９９１７２７８２）株式会社アライドエンジニアリング (3)
【Ｆターム（参考）】
[ Back to top ]