(57)【要約】【課題】メッシュ間のデータ授受情報を求める計算を効率的に行う。【解決手段】端部が重なるメッシュを物体毎に作成する（Ｓ１１）。２以上の物体が相対移動した時に、移動方向に位置するセルを物体の移動につれて伸縮させ、物体とセルの位置関係を維持してメッシュを変形するとともに、前記メッシュは、メッシュ間の重畳部において、隣接する境界線に垂直な方向におけるメッシュの間隔が変化しないように維持して、変形する（Ｓ１４）。２以上の物体が相対移動する空間をメッシュによりセルに分割し、セル毎に状態を特定することで物体周りの流体解析を行う（Ｓ１３）。【選択図】図１
2009211555.tif 000002
【特許請求の範囲】
【請求項１】
２以上の物体が相対移動する空間をメッシュによりセルに分割し、セル毎に状態を特定することで物体周りの流体の解析を行う装置であって、
端部が重なるメッシュを物体毎に作成し、２以上の物体が相対移動した時に、移動方向に位置するセルを物体の移動につれて伸縮させ、物体とセルの位置関係を維持してメッシュを変形するとともに、
前記メッシュは、メッシュ間の重畳部において、隣接する境界線に垂直な方向におけるメッシュの間隔が変化しないように維持して、変形することを特徴とする解析装置。
【請求項２】
請求項１に記載の解析装置において、
前記メッシュは、水平な方向の両端を移動させずに変形することを特徴とする解析装置。
【請求項３】
２以上の物体が相対移動する空間をメッシュによりセルに分割し、セル毎に状態を特定することで物体周りの流体の解析をコンピュータに実行させるプログラムであって、
コンピュータに、
端部が重なるメッシュを物体毎に作成させ、２以上の物体が相対移動した時に、移動方向に位置するセルを物体の移動につれて伸縮させ、物体とセルの位置関係を維持してメッシュを変形させるとともに、
前記メッシュを、メッシュ間の重畳部において、隣接する境界線に垂直な方向におけるメッシュの間隔が変化しないように維持して、変形させることを特徴とする解析プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
２以上の物体が相対移動する空間をメッシュによりセルに分割し、セル毎に状態を特定することで物体周りの流体の解析を行う解析装置および解析プログラムに関する。
【背景技術】
【０００２】
従来より、流体の解析を行うシミュレータが知られており、各種の流体解析に利用されている。例えば、自動車の走行時の空気流れ解析があり、ボディーのデザインの際などに広く利用されている。なお、シミュレータは、通常汎用コンピュータにシミュレーションプログラムをインストールして行われる。
【０００３】
ここで、自動車がトラックなど他車から追い抜かれた場合、流体力によって車体が揺すられ、運転手は危険を感じる場合がある。そこで、このような追い抜き時の現象を扱うための物体が相対運動する際の流れ解析が望まれる。
【０００４】
自動車まわりの流れ計算において、例えば、空力抵抗係数（Ｃｄ）のような時間平均量を算出する場合には、時間に依存しない方程式を用いる定常計算がよく用いられる。一方、物体が相対運動する計算では、現象が定常ではないため、非定常計算をする必要がある。非定常計算では、一般的に微少な時間を時間ステップとした多くの反復計算により、時間を進行させる。このため、定常計算に比べて計算時間が大きい。計算時間の増大量は一概には示すことはできないが、定常計算と比較すると、おおよそ１０倍から１００倍であると考えられる。
【０００５】
また、多くの流体計算手法では、物体まわりにメッシュ（格子）を作成する。２つ以上の物体が相対運動する場合、メッシュを変形、または、移動させる必要が生じる。このため、非定常計算により計算時間が増大するのに加えて、このメッシュの変形・移動に関わる計算時間が発生する。計算方法によっては、メッシュ変形・移動に関わる計算時間の方が大きくなる場合もある。
【０００６】
従って、相対運動する物体まわりの流体を計算する際の、メッシュの変形・移動に関わる計算時間をできるだけ短縮したいという要求がある。
【０００７】
【非特許文献１】K. Nakahashi and F. Togashi, "Intergrid-Boundary Definition Method for Overset Unstructured Grid Approach", AIAA J. Vol.38、 No. 11, 2000
【非特許文献２】高村他、「自動車の追い越し空力シミュレーション」、自技会学術講演前刷集、No. 8-01, 20015049, 2001
【発明の開示】
【発明が解決しようとする課題】
【０００８】
ここで、上述のようなメッシュの変形、移動については、以下の３つの技術が挙げられる。
【０００９】
＜メッシュ変形法＞
メッシュの構成を変化させずに、単にメッシュの各セルの変形のみで物体の運動を扱う方法である。この方法には、メッシュが大きく歪んだ場合に計算精度が悪化する課題がある。さらに、例えば、車両が遥か後方から前方へ他車を追い抜く場合には、メッシュの変形のみではほとんど対応できないという課題がある。
【００１０】
＜移動する重合格子（メッシュ）法＞
複数のメッシュを重ね合わせて、セルを変形させずに、それらを相対移動させる方法である（非特許文献１参照）。この方法では、お互いのメッシュの間でデータを授受する必要があり、この授受するセル位置と空間補間するため係数を計算するための時間が大きい点が課題である。
【００１１】
＜リメッシュ法＞
メッシュは全領域で１つであり、時間を追ってメッシュを作り直す方法である。メッシュの変形では、大きく相対位置が変わる場合には対応できないのに対し、メッシュを作り直せば、ほとんど任意の運動に対応できる。ただし、この方法では、メッシュを再作成するための時間が大きい点が課題である。
【課題を解決するための手段】
【００１２】
本発明は、２以上の物体が相対移動する空間をメッシュによりセルに分割し、セル毎に状態を特定することで物体周りの流体の解析を行う装置であって、端部が重なるメッシュを物体毎に作成し、２以上の物体が相対移動した時に、移動方向に位置するセルを物体の移動につれて伸縮させ、物体とセルの位置関係を維持してメッシュを変形するとともに、前記メッシュは、メッシュ間の重畳部において、隣接する境界線に垂直な方向におけるメッシュの間隔が変化しないように維持して、変形することを特徴とする。
【００１３】
また、本発明は、２以上の物体が相対移動する空間をメッシュによりセルに分割し、セル毎に状態を特定することで物体周りの流体の解析をコンピュータに実行させるプログラムであって、コンピュータに、端部が重なるメッシュを物体毎に作成させ、２以上の物体が相対移動した時に、移動方向に位置するセルを物体の移動につれて伸縮させ、物体とセルの位置関係を維持してメッシュを変形させるとともに、前記メッシュを、メッシュ間の重畳部において、隣接する境界線に垂直な方向におけるメッシュの間隔が変化しないように維持して、変形させることを特徴とする。
【発明の効果】
【００１４】
本発明は、重合格子法とメッシュ変形法を合成した相対運動する物体まわりの流れ計算を行う。重合格子法では、移動物体に適用すると計算時間が大きくなる問題点があるが、本発明では、重ね合わせるメッシュにある制限を付けることによって、計算時間を小さくする。
【発明を実施するための最良の形態】
【００１５】
以下、本発明の実施形態について、図面に基づいて説明する。
【００１６】
図１に、本実施形態の全体処理の手順を示す。まず、初期メッシュを作成する。この初期メッシュは１つの物体についてそれを取り囲む１つのメッシュとする（Ｓ１１）。各メッシュは、その端部同士で重畳するように形成する。次にメッシュ間のデータ授受情報を計算する（Ｓ１２）。これは、重畳している部分におけるセル間におけるデータを授受する相手セルの位置と、その相手との補間係数の計算である。Ｓ１２の計算が終了した場合には、主計算を行う（Ｓ１３）。これは、メッシュを利用した物体周りの流体の演算であり、複数メッシュ間では、Ｓ１２において計算した相手セルの位置と補間係数を利用して各種物理量の授受を行い、全体としての流体計算を行う。
【００１７】
Ｓ１３の主計算が終了した場合には、次の時間の位置に物体を移動するとともにメッシュを変形させて、Ｓ１２に戻り、これを終了まで繰り返す。また、Ｓ１３の各時点での計算結果は、それぞれ出力され（Ｓ１５）、後処理を経て可視化される（Ｓ１６）。このようにして、例えば、車両のすれ違い時における周辺流体の解析が行える。
【００１８】
＜初期メッシュの作成（Ｓ１１）＞
２つの物体（直方体）が相対運動する場合を例に挙げ、初期メッシュの作成について説明する。なお、このメッシュの構成は、物体が３つ以上の場合においても応用できる。
【００１９】
図２に初期メッシュの例を示す。下側のメッシュをメッシュ１、上側のメッシュをメッシュ２とする。ここでは、メッシュ１は変形せず、メッシュ２のみ変形するものとする。なお、本実施形態は、メッシュ１が変形する場合にも、適用可能である。本実施形態に係る計算は、この２つのメッシュの関係に特徴がある。
【００２０】
メッシュ構成方法を以下に示す。（ξ，η）を計算空間座標とし、図２のようにξ方向線とη方向線を定義する。メッシュのξ方向のノード番号（メッシュ点番号）をｉ＝０，１，２，３，４，．．．，ｎｘとする。このとき、例えば、ξ＝１．５は、物理空間では、ｉ＝１とｉ＝２のノード点の中点を指す。
【００２１】
ξ方向の全ての行に対し、η方向にお互いのメッシュを複数の層にわたって重畳させる。何層重畳するかは、流体計算アルゴリズムに依存するが、本実施形態においては、５層重ねるものとする。η方向をｊ＝０，１，２，．．．，ｎｙとすると、メッシュ１のｊ＝ｎｙ−５，．．．，ｎｙまでのメッシュがメッシュ２のｊ＝０，１，２，．．．，５までのメッシュと重なる。この際、重ねた層においては、ξ方向は、端の点０と点ｎｘでそれぞれ相手のメッシュのノード点を一致させておく。ただし、点０と点ｎｘ以外では、一致させる必要はない。また、各メッシュのｎｘの値が異なっていてもよい。さらに、この層内のη方向のメッシュ間隔をξ方向に渡って一定にすることと、各η線がそれぞれξ方向に対し、直角方向であるようにメッシュを作成する。
【００２２】
＜メッシュの変形（Ｓ１４）＞
ここで、初回以外はＳ１４において、変形したメッシュについてＳ１２の計算を行うため、メッシュの変形について先に説明する。
【００２３】
メッシュを上記初期メッシュの条件を満たしながら、時間進行に伴い、徐々にメッシュを変形する。
【００２４】
本実施形態では、２つの物体が相対移動した時に、移動方向に位置するセルを物体の移動につれて伸縮させ、物体とセルの位置関係を維持してメッシュを変形する。本実施形態では、各メッシュの左右両端の位置は変更せず、内部のメッシュ線のみを移動する。従って、データ授受情報を計算する際に、対応メッシュが限定される。
【００２５】
ここで、前記隣接する２つのメッシュは、隣接する境界線に垂直な方向におけるメッシュ線の間隔が変化しないように変形する。すなわち、重畳部で、セル列方向の間隔は変化させない。
【００２６】
＜データ授受情報の算出（Ｓ１２）＞
データ授受情報の計算は、一般的に、３次元空間で実施される。ただし、ここでは説明のため、単純化して、奥行き方向（紙面垂直方向）のメッシュは一定であるとし、２次元空間において考える。
【００２７】
メッシュ１とメッシュ２とは、重合格子法により接続されている。このため、両メッシュ１，２の端部で、一方のメッシュから計算物理量を他のメッシュに与える必要がある。図５には、実線でメッシュ１、点線でメッシュ２を示している。この図で、メッシュ２の点Ｐについてメッシュ１における対応セルおよび補間係数を算出する場合について説明する。なお、メッシュのデータ補間を行う点Ｐは、図６に示すように、メッシュの交点であるノードとする場合と、セルの中心とする場合があるが、以下の説明では、セルの中心とする。なお、ノード毎に物理量を補間する場合もあるが、これはセル毎に計算したものと方法は同様であり、本出願においてはセル毎に計算したものに含まれる。
【００２８】
計算物理量を受けるメッシュ（この場合はメッシュ２）の点Ｐに対して、計算物理量を渡すメッシュ（この場合はメッシュ１）の（ξ，η）を算出し、（ξ，η）の値からデータを内挿する。（ξ，η）は、点Ｐが他方のメッシュにおけるセルの位置と、その補間係数を表す。この（ξ，η）の算出は以下のように行う。
【００２９】
点Ｐ（Ｘ，Ｙ）が一方のメッシュのどの位置のメッシュセルに入っているかは、あらかじめ判らないので、一旦、仮に（ｉ，ｊ）、（ｉ＋１，ｊ）、（ｉ，ｊ＋１）、（ｉ＋１，ｊ＋１）の４頂点のセルの中に入っているものとする（各位置の物理座標をｘｉ，ｊのように表す）。すると、線形補間式より、以下のように（Ｘ，Ｙ）が表せる。
【００３０】
【数１】
2009211555.tif 000003
【００３１】
【数２】
2009211555.tif 000004
ここで、ξ’＝ξ−ｉ，η’＝η−ｊである。
【００３２】
既知の量（Ｘ，Ｙ）と（ｘｉ，ｊ，ｙｉ，ｊ），．．．，（ｘｉ＋１，ｊ＋１，ｙｉ＋１，ｊ＋１）により、ξ’，η’が定まり、ξ，ηが求まる。ただし、この方程式は解析的には解けないため、反復法にて計算しなければならず、計算時間は大きい。
【００３３】
点Ｐが（ｘｉ，ｊ，ｙｉ，ｊ），．．．，（ｘｉ＋１，ｊ＋１，ｙｉ＋１，ｊ＋１）のセル内（境界上を含む）に入っている条件は、
０≦ξ’≦１かつ０≦η’≦１          （３）
である。すなわち、図７に示すように、メッシュ１のセルにおいて、メッシュ２のＰ点の位置は、ξ’，η’により決定され、セル幅を１とした場合に、点Ｐが存在するセルを式（３）により検出することができる。
【００３４】
この計算を各メッシュセルに対して実行し、式（３）の条件を満たすセルを探索する必要がある。反復計算にかかる時間をｓとし、探索対象となるセル数がＮとすると、計算時間は１つの点Ｐに対し、最大ｓ×Ｎである。
【００３５】
しかし、本実施形態では計算の構成１、２のメッシュの条件の下では、η方向には変形しないため、η，η’は既知である。すると、未知量はξのみとなり、式（１）と（２）を連立させる必要はなくなり、以下の解析式で算出でき、反復計算が不要になる。
【００３６】
【数３】
2009211555.tif 000005
【００３７】
従って、図８に示すように、あるセルに対して、ξ’を求め（Ｓ２１）、求めたξ’が０≦ξ’≦１の条件を満たすかを判定する（Ｓ２２）。そして、条件を満たさない場合には、セルを移動し（Ｓ２３）、Ｓ２１に戻る。そして、条件を満たすまでこれを繰り返し、Ｓ２２において条件を満たすことで、セルが特定され、補間係数ξ’が特定される（Ｓ２４）。
【００３８】
このように、本実施形態では、次のような特徴がある。
（ｉ）２つのメッシュを一定の間隔の層で重ね、その層内にある、層に垂直な方向のメッシュの間隔を一定とし、重なりあった部分の垂直方向の計算空間方向線が互いのメッシュで等間隔である。
（ｉｉ）以上のメッシュ生成方法を保ちながら計算メッシュを変形させることにより、重合メッシュ間のデータ授受情報を短時間で算出して、相対運動する物体まわりの流れを高速に計算できる。
【００３９】
＜実施形態の効果＞
実際の計算は、３次元であるが、前述のように奥行き方向のメッシュは変化しないと考え、まず２次元空間での計算時間を見積もる。計算時間のオーダーを見積もるため、ｉ，ｊの各方向のメッシュ数を代表させてそれぞれｎｘ点とする。
【００４０】
・本実施形態による計算時間：
探索する点Ｐの数ＮＰは、メッシュ１つに対して最大ｎｘ点である（重合格子法のアルゴリズムにより２層の情報を授受する場合もあるが、流体計算ソフト固有の問題であるため、ここでは１層のみと考え、ｎｘ点としておく）。そして、２つのメッシュが重合しているため、これらの間でお互いに情報を授受することから、ＮＰ＝２×ｎｘ点になる。次に、探索しなければならないセルの数は、ｉ方向（ξ方向）のみ考慮すればよいことから、Ｎ＝ｎｘである。反復計算は１０回で収束するものと考えれば、本実施形態の計算ではせいぜい１０分の１の計算時間でξ、η（ηは既知）が求まるため、１点あたりが１セルに対して求める計算時間は、Ｓ＝０．１×ｓとなる。ここで、Ｓは前述のように反復計算にかかる時間である。結局、本計算方法のメッシュの移動に伴う計算時間は、
ＮＰ×Ｎ×Ｓ＝（２×ｎｘ）×（ｎｘ）×（０．１×ｓ）      （５）
と見積もられる。
【００４１】
・移動する重合格子法（従来技術）の計算時間：
図３のように、メッシュ２がメッシュ１の内部に含まれるようなメッシュの構成となる。メッシュ２は全周囲でメッシュ１と重合しているため、点Ｐの数は、本計算方法の４倍となり（ただし、メッシュ２の点Ｐは、層状な界面ではなく、メッシュ１の内部にある点全てを探索する必要があり、さらに点数が多くなると考えられる。ただし、簡単のため４倍と見積もる）、ＮＰ＝２×４×ｎｘとなる。また、探索するセルは、２次元面上を探索しなければならず、最大、Ｎ＝ｎｘ×ｎｘとなる。よって、計算時間は、
ＮＰ×Ｎ×Ｓ＝（２×４×ｎｘ）×（ｎｘ×ｎｘ）×（ｓ）        （６）
と見積もられる。
【００４２】
式（６）と式（５）と比較すると、移動する重合格子法では（４×ｎｘ×１０）倍の計算時間を要することがわかる。ｎｘ＝１００〜２００程度とすると、数千倍の計算時間を要することになる。
【００４３】
実際の流れの計算は、全メッシュ点数に依存した大きさのマトリックスからなる方程式を反復して計算する。移動する重合格子法では、全メッシュ点数に依存してセルの探索時間が増大することから、実際の流れの計算時間に対し、データ授受情報を算出するための時間は、大まかな見積として同等であると考えられる。
【００４４】
非特許文献１では、計算上の工夫により、この情報の授受に関する計算時間の低減が図られている。ただし、この非特許文献１によると、メッシュの移動のための計算時間は、全体の約３１％を占めている。本実施形態の計算手法では、この計算時間がほとんど無視できる程度になる。
【００４５】
・リメッシュ法の計算時間、計算精度：
リメッシュ法の場合は、そのアルゴリズムに大きく依存するため、前述のような計算時間の見積は難しい。ただし、メッシュを再構築する場合には、例えば、楕円方程式を用いて反復計算する方法がしばしば用いられる。よって、メッシュ点数に応じた反復計算が必要であると考えられ、流れ計算と同等の計算時間が必要になる。つまり、リメッシュ法においても移動する重合格子法と同程度の計算時間が必要であると考えられる。
【００４６】
非特許文献２では、リメッシュの計算時間を削減するため、メッシュ変形とリメッシュを組み合わせている。計算例では、リメッシュ回数をできる限り減らす工夫をしており、計算実行ステップの２００回に１回のみリメッシュをしている。しかし、この場合でも、リメッシュに要する計算時間は小さくなく、全体の１９％を占めている。
【００４７】
・３次元の場合
さらに、３次元の場合には、本実施形態の効果がさらに大きくなる。メッシュ点数を一方向ｎｘとすると、３次元の場合総点数はｎｘ×ｎｘ×ｎｘとなる。この場合の計算時間の見積りは図９に示すようになる。すなわち、補間される点Ｐの個数は、従来の場合２×６×ｎｘ、本実施形態の場合２×ｎｘ、セル移動回数は、従来の場合ｎｘ×ｎｘ×ｎｘ、本実施形態の場合ｎｘ×ｎｘ、式（４）を計算する時間は、従来の場合ｓ、本実施形態の場合０．１ｓとなる。実際の計算時間は、（補間される点Ｐの個数）×（セル移動回数）×（式（４）を計算する時間）となる。このため、３次元の場合において、本実施形態の効果は２次元の場合に比べて各段に大きくなる。
【００４８】
＜適用例＞
図２の構成のメッシュを実際に変形させた例を図４に示す。ここでは、メッシュの変形は、ξ方向の疎密の変化のみで与えている。この変形により、メッシュ２の物体がメッシュ１の物体を追い抜く現象を模擬している。このように、この例では、メッシュ１とメッシュ２の上下方向（η方向）は変形せず、左右（ξ方向）について、両端固定した状態で、内部の縦方向のメッシュを移動する。
【００４９】
この例では、相対運動を扱うための重合・変形に関わる時間はほとんど無視できる程度である。
【図面の簡単な説明】
【００５０】
【図１】実施形態の処理の手順を示すフローチャートである。
【図２】２つのメッシュの関係を示す図である。
【図３】メッシュ２がメッシュ１の内部に含まれるようなメッシュの構成を示す図である。
【図４】物体が移動した場合のメッシュの変形状態を示す図である。
【図５】２つのメッシュが重なった部分の対応を示す図である。
【図６】セルとノードの関係を示す図である。
【図７】補間係数の定義を示す図である。
【図８】補間係数の特定手順を示すフローチャートである。
【図９】計算時間を示す図である。
【図１】
2009211555.tif 000006
【図２】
2009211555.tif 000007
【図３】
2009211555.tif 000008
【図５】
2009211555.tif 000009
【図６】
2009211555.tif 000010
【図７】
2009211555.tif 000011
【図８】
2009211555.tif 000012
【図９】
2009211555.tif 000013
【図４】
2009211555.tif 000014