WO2012026383A1

WO2012026383A1 - 計算用データ生成装置、計算用データ生成方法及び計算用データ生成プログラム

Info

Publication number: WO2012026383A1
Application number: PCT/JP2011/068662
Authority: WO
Inventors: 植村　健; 一貴柳原; 恒洋齊藤
Original assignee: 旭硝子株式会社
Priority date: 2010-08-24
Filing date: 2011-08-18
Publication date: 2012-03-01
Also published as: CN103080941A; EP2610772A4; US20130173239A1; KR20130098306A; RU2013112926A; EP2610772A1; JP5045853B2; JPWO2012026383A1

Abstract

　数値解析を行う際の計算用データを生成する。　パラメータ記憶手段に記憶されているパラメータに基づき、対象物体を含む解析領域内を複数の直方体に分割したボクセルデータを定義し、各ボクセルに対して、ボクセル属性を付与して、ボクセルデータ記憶手段に記憶する手段と、ボクセルの中心点を使用して、ボクセルの数より少ない初期点データを生成して、初期点データ記憶手段に記憶する手段と、ボクセルデータ記憶手段に記憶されている各ボクセルの属性と、初期点データ記憶手段に記憶されている初期点データとに基づき、対象物体内を、複数のボクセルで構成する分割領域を定義し、定義した分割領域データを分割領域データ記憶手段に記憶する手段と、分割領域データ記憶手段に記憶されている分割領域データに基づき、各分割領域の境界面データを生成し、該境界面データを計算用データとして計算用データ記憶手段に記憶する手段とを備える。

Description

計算用データ生成装置、計算用データ生成方法及び計算用データ生成プログラム

　本発明は、数値解析を行うための計算用データを生成する計算用データ生成装置、計算用データ生成方法及び計算用データ生成プログラムに関する。

　従来から、流速分布、応力分布及び熱分布等を数値解析によって求めるための数値解析手法として、例えば、有限要素法、有限体積法、ボクセル法及び粒子法が知られている。このような数値解析手法は、一般的にプリ処理とソルバ処理とポスト処理とから構成されている。そして、プリ処理の中で計算用データモデルを作成し、ソルバ処理にて当該計算用データモデル及び離散化された支配方程式（以下、離散化支配方程式と呼ぶ）を用いて上記物理量の計算を行っている。

　従来の有限体積法は、例えば、解析領域を複数の領域に分割し、各分割領域の体積、隣り合う分割領域の境界面の面積及び当該境界面の法線ベクトルを用いて各分割領域における物理量を計算する。有限体積法では、プリ処理にて各分割領域の頂点の座標（Vertex）を含む計算用データモデル（通常、メッシュと呼ばれる）を作成し、ソルバ処理にて当該計算用データモデルに含まれるVertex等を用いて前述の分割領域の体積、境界面の面積及び境界面の法線ベクトルを算出し、これらの値を用いて物理量の計算を行う。Vertexは分割領域の幾何学的形状を規定するための量である。よって、有限体積法においては、ソルバ処理にて分割領域の幾何学的形状を用いて、分割領域の体積、境界面の面積及び境界面の法線ベクトルの算出を行っていると言える。

　さらに、有限体積法では、隣り合う分割領域での頂点共有の条件が部分的に満足されていない部分を有してもよい。このため、有限体積法では、分割領域に対する制約が若干緩和される場合があるが、利用する解析要素タイプは例えばテトラ要素、ヘキサ要素、プリズム要素、ピラミッド要素などに限定される。

　なお、特許文献１に示すように、解析要素タイプを限定しない有限体積法も提案されている。ただし、このような解析要素タイプを限定しない有限体積法であっても、前述した従来の有限体積法と同様に、プリ処理にて各分割領域の頂点の座標（Vertex）を含む計算用データモデルを作成し、ソルバ処理にて当該計算用データモデルに含まれるVertex等を用いて物理量の計算を行っている。

　また、有限要素法は、周知のように補間関数を用いて各分割領域における物理量を算出する方法であるが、有限体積法と同様に、ソルバ処理にてVertex等によって規定される分割領域の幾何学的形状を用いている。

　ボクセル法及び粒子法は、有限要素法や有限体積法と比較して、容易に計算用データモデルが作成可能な数値解析手法である。ボクセル法は、基本的に同一サイズの直方体形状の複数のボクセル（直交格子）によって解析領域を定義するボクセルデータを計算用データモデルとして作成し、このボクセルデータを用いた物理量計算を行うことによって数値解析を行う方法である。ボクセル法としては、重み付き残差積分法による支配方程式を用いる重み付き残差積分法型と、例えばセル－オートマトンモデルや格子－ボルツマンモデル等を用いる非積分法型に大別される。そして、このボクセル法によれば、ボクセルデータとして、Vertex等は必要ない。このようなボクセル法によれば、解析領域をボクセルによって分割することによって解析領域を容易に定義することができ、短時間で計算用データモデルを作成することができる。

　一方、粒子法は、解析領域を複数の粒子によって定義する粒子データを計算用データモデルとして作成し、この粒子データを用いた物理量計算を行うことによって数値解析を行う方法である。粒子法は、支配方程式として、非積分法型で粒子間相互作用モデルを利用している。粒子法では、分割領域がないため、Vertex等を必要としない。このような粒子法によれば、解析領域に例えば粒子を均一に配置することによって解析領域を容易に定義することができ、計算用データモデルを短時間で作成することができる。

米国特許出願公開第２００８／００２１６８４号明細書

　従来の有限要素法や有限体積法等の数値解析手法のように、ソルバ処理において分割領域の幾何学的形状を用いる場合には、当然ながら計算用データモデルに対して、分割領域の幾何学的形状を示すデータを持たせることが必須となる。
　分割領域の幾何学的形状を定義するためには、Vertexに加えて、頂点の連結情報（Connectivity of Vertex、以下ではConnectivityと略記する）が必要となる。したがって、有限要素法や有限体積法では、計算用データモデルにVertexとConnectivityと持たせる必要がある。

　なお、具体的には、Connectivityは、全分割領域の頂点に対して順番に付された全体ノード番号と、１つの分割領域内において頂点に対して順番に付された局所ノード番号との対応情報によって定義される。

　このようなVertexとConnectivityとを有する計算用データモデルは、周知のように作成に非常に膨大な作業が必要となる。例えば、有限要素法で用いられる計算用データモデルでは、図１に示すように、隣接する分割領域が必ずVertexを共有するという条件を満足するように計算用データモデルを作成する必要があり、全ての分割領域がこの条件を満足させるためには、非常に膨大な時間を必要とする。

　一方、有限体積法で用いられる計算用データモデルは、図２に示すように、隣接する分割領域に共有されないVertexの存在を許容できる分、有限要素法と比較すればメッシュ生成の自由度が増す。しかしながら、有限体積法においても、共有されないVertexが少なくとも隣接する分割領域の辺上に存在すること、また一般的に分割領域の形状を予め設定された解析要素タイプに合わせることといった条件の中で計算用データモデルを作成する必要があり、メッシュ生成の自由度が高いとは言い難い。

　また、近年においては、３次元ＣＡＤ（Computer Aided Design）データ等の３次元形状データから抽出される解析領域に対する数値解析が行われている。しかしながら、３次元形状データは、数値解析用に形成されたデータではなく、面の重なり、面の交差、面間の隙間、微小穴等を示すデータが含まれており、VertexとConnectivityとを有する計算用データモデルの作成に適さない条件を多く含んでいる。このため、これらのVertexとConnectivityとを有する計算用データモデルの作成が出来るように、３次元形状データを修正あるいは変更する必要がある。そして、VertexとConnectivityとを有する計算用データモデルが作成可能なように、３次元形状データを修正あるいは変更するためには、経験や試行錯誤を要する非常に膨大な手作業を行う必要がある。これが、有限要素法や有限体積法を実用で利用する際の大きな問題である。

　また、有限体積法のようにソルバ処理において分割領域の体積、境界面の面積及び境界面の法線ベクトルの算出を行う場合には、ソルバ処理における計算量がさらに増えることとなり、ソルバ処理における計算負荷がさらに増大する。

　ボクセル法では、短時間で計算用データモデルを作成することができるものの、以下の点で問題がある。ボクセル法は、基本的に解析領域が全て同一サイズのボクセル（直交格子）によって定義される。通常、有限要素法や有限体積法では、より高い解析精度を得たい領域の要素サイズ（分割領域のサイズ）を小さく設定することによって当該領域について正確な物理量計算を行い、また他の領域の要素サイズを大きく設定することによって当該領域についての計算負荷を低減させている。ところが、ボクセル法においては、全てのボクセルが基本的に同一サイズであるため、ボクセルを小さく設定した場合には計算負荷が非常に大きくなり、ボクセルを大きく設定した場合には解析精度が悪化することとなる。

　また、ボクセル法では、同一サイズのボクセル（直交格子）を配列することによって解析領域を定義する必要があることから、外部領域との境界付近において解析領域を滑らかにすることができずに階段状となる場合がある。つまり、実際に解析したい領域が斜面や曲面等を有している場合であっても、ボクセルデータにおいてその領域は階段状に表されることとなる。このため、ボクセル法における解析領域形状が実際に解析したい領域形状と異なることとなり、解析精度は悪化する。

　これに対して、ボクセルデータの階段状の領域を、実際に解析したい領域が有する斜面や曲面に沿って切断（境界補正）するカットセル法と呼ばれる改良方法が提案されている。しかしながら、この改良方法によると、この境界補正によって非常に小さな分割領域が生成されやすく、このような小さな分割領域が生成された場合には、解析精度の悪化を招く。また、この改良方法では、カットセルの形成のため及びソルバ処理において、Vertexを利用することとなる。

　以上のように、境界補正を行わないボクセル法においては、Vertex等が必要ないものの、ボクセルの生成、いわゆるメッシュ生成には限界がある。すなわち、十分な解析精度を得ようとすると、ボクセルの数が増え、ソルバ処理における計算負荷も増え、問題となる。さらに、境界補正を行うボクセル法の改良方法では、結果としてVertexが必要となるため、結局、分割領域の幾何学的形状の影響を受けることなり、外部領域との境界周辺での分割領域形成のための処理に、経験や試行錯誤を要する非常に膨大な手作業を伴い、形状データモデルを短時間に作成することができないことになる。

　一方、粒子法では、ある特定の粒子と他の粒子との結合関係を算出する必要がある。このため、当該特定粒子の近傍に存在する粒子を探索する必要がある。そして、この粒子の近傍探索処理は、原則として全ての粒子に対して行われる。しかしながら、粒子法においては、時刻の変化と共に各粒子が移動し、これによって常に粒子同士の結合関係が変化する。このため、解析における時刻が変化するたびに近傍探索処理を行う必要があり、計算負荷の増大を招く。近傍探索の対象となる粒子を厳選する等して近傍探索処理の計算負荷を低減させる試みがなされているものの、例えば、解析精度を向上させるために粒子数を増大させた場合には、粒子数の２乗に比例して計算負荷が増大する。

　このような粒子法において実用時間内における数値解析を実現するためには、大型の並列計算機で多くのＣＰＵ（Central Processing Unit）を使用する必要がある。例えば、百万個程度の粒子を用いた粒子法による数値解析が、VertexとConnectivityを使用した一般的な有限体積法ソルバで、１ＣＰＵで半日の計算が、粒子法では３２ＣＰＵを使用した並列計算で１週間以上かかった実例がある。また、粒子法においても、粒子を密に配置された場合には計算負荷が非常に大きくなり、粒子を粗く配置した場合には解析精度が悪化することとなる。

　さらに、粒子法では、後に詳説するが、流体、構造、熱、拡散等の物理量の保存則に基づく物理現象を解析する場合には、その保存性が十分に満足されていない。例えば、解析領域と外部領域との境界面に臨んで配置される粒子が、境界面においてどれだけの面積を有しているかの情報がない。このため、境界面から熱を入熱させる条件を与えようとした場合であっても、各粒子にどれだけの熱が入熱されるかを正確に把握することができず、精度の良い定量値が得られない。

　本発明は、前述する従来の数値解析手法である有限要素法、有限体積法、ボクセル法、ボクセル法の改良手法、及び粒子法の問題点に鑑みてなされたもので、解析精度の悪化を伴わずにソルバ処理における計算負荷の低減を図ることができる数値解析装置に入力する計算用データの作成における作業負担を軽減することができる計算用データ生成装置、計算用データ生成方法及び計算用データ生成プログラムを提供することを目的とする。

　上記課題を解決するために、本発明は、物理現象を数値的に解析する数値解析方法において物理量を計算する物理量計算方法であって、直交格子形状のみによらない複数の分割領域に分割された解析領域における物理量を計算する物理量計算工程を含み、該物理量計算工程にて、前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された支配方程式と、各前記分割領域の体積及び隣り合う前記分割領域同士の境界面の特性を示す境界面特性量を、前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する計算用データモデルと、を用いて前記物理量を計算するという構成を採用する。

　本発明で用いられる離散化支配方程式は、従来のように分割領域の幾何学的形状を規定する量（VertexとConnectivity）を含んだ形式で表現されるものではなく、分割領域の幾何学的形状を規定する量を必要としないものである。本発明で用いられる離散化支配方程式は、従来の幾何学的形状を規定する量を使用する方程式を重み付き残差積分法に基づいて導出する過程で敢えて途中にて留めることによって得ることができる。このような本発明で用いられる離散化支配方程式は、分割領域の幾何学的形状を必要としない量（すなわちVertexとConnectivityとを必要としない量）で表現され、例えば分割領域の体積と境界面特性量の２つのみに依存する形式とすることができる。

　つまり、従来の有限要素法や有限体積法では、前提として解析対象物を微小領域に分割するため、この微小領域の幾何学的形状を規定する量、すなわちVertexとConnectivityを用いることを前提にして離散化支配方程式の導出をしているが、本発明で用いられる離散化支配方程式は、従来と異なるまったく新しい発想に基づいて導出されるものである。そして、本発明は、このような新しい発想に基づいて導出された離散化支配方程式を用いることを特徴とするものであり、従来の数値解析方法と異なり、幾何学的形状に依存せず、従来の問題を解決し、種々の顕著な効果を奏する。

　ここで、分割領域の体積と境界面特性量とが、分割領域の特定の幾何学的形状を規定するVertexとConnectivityとを必要としない量であることについて説明する。なお、VertexとConnectivityとを必要としない量とは、VertexとConnectivityとを用いなくとも定義が可能な量である。例えば、分割領域の体積について考えると、分割領域の体積がある所定の値となるための分割領域の幾何学的形状は複数存在する。つまり、体積がある所定の値をとる分割領域の幾何学的形状は、立方体である場合や球である場合も考えられる。そして、例えば、分割領域の体積は、全分割領域の総和が解析領域全体の体積と一致するという制約条件の下で、例えば分割領域の体積が隣接分割領域との平均距離の３乗にできるだけ比例するような最適化計算により定義することができる。したがって、分割領域の体積は、分割領域の特定の幾何学的形状を必要としない量（VertexとConnectivityとを必要としない量）と捉えることができる。

　また、境界面特性量としては、例えば境界面の面積や、境界面の法線ベクトル、境界面の周長等が考えられるが、これらの境界面特性量がある所定の値となるための分割領域の幾何学的形状（すなわち境界面の幾何学的形状）は複数存在する。そして、例えば、境界面特性量は、各分割領域を取り囲む全境界面に対して、法線ベクトルの面積加重平均ベクトルの長さがゼロとなる制約条件の下で、境界面の法線ベクトルの方向を隣接する２つの分割領域のコントロールポイント（図５参照）を結ぶ線分に近づけ、かつ、分割領域を取り囲む全境界面面積の総和が当該分割領域の体積の２分の３乗にできるだけ比例するような最適化計算により定義することができる。したがって、境界面特性量は、分割領域の特定の幾何学的形状を必要としない量（VertexとConnectivityとを必要としない量）と捉えることができる。

　また、本発明において「直交格子形状のみによらない複数の分割領域に分割された解析領域」とは、解析領域を構成する複数の分割領域の少なくともいずれかが直交格子形状を取らないことを意味する。つまり、解析領域が直交格子形状以外の形状の分割領域を含むことを意味する。また、本発明において「VertexとConnectivityとを必要としない量のみを使用する」とは、離散化支配方程式に代入される値がVertexとConnectivityとを必要としない量のみであることを意味する。

　次に、図３の概念図を参照して、本発明を用いた数値解析手法と従来の数値解析手法とにおけるプリ処理及びソルバ処理を対比しながら、本発明の顕著な効果についてより詳細な説明を行う。

　本発明を用いる数値解析手法の場合には、図３に示すように、ソルバ処理（本発明の物理量計算工程）にて、VertexとConnectivityとを必要としない量のみを使用する離散化支配方程式を用いて分割領域における物理量の算出が行われる。このため、離散化支配方程式を解くにあたり、プリ処理にて作成される計算用データモデルにVertexとConnectivityとを含める必要がない。

　そして、本発明を用いる場合には、VertexとConnectivityとを必要としない量として、分割領域の体積と境界面特性量とが使用される。このため、プリ処理にて作成される計算用データモデルは、VertexとConnectivityとを持たず、分割領域の体積と、境界面特性量と、その他補助データ（例えば、後述する、分割領域の係合情報やコントロールポイント座標等）とを有するものとなる。

　このように本発明を用いた場合には、前述のように、分割領域の体積と上記境界面特性量、すなわち分割領域の幾何学的形状を必要としない量に基づいて各分割領域における物理量が計算可能とされる。このため、計算用データモデルに、分割領域の幾何学的形状、すなわちVertexとConnectivityとを持たせることなく物理量を算出することが可能となる。したがって、本発明を用いることにより、プリ処理において、少なくとも分割領域の体積と境界面特性量（境界面の面積及び境界面の法線ベクトル）とを有する計算用データモデルを作成すれば良くなり、VertexとConnectivityとを有する計算用データモデルを作成することなく物理量の計算を行うことが可能となる。

　VertexとConnectivityとを持たない計算用データモデルは、分割領域の幾何学的形状を必要としないため、分割領域の幾何学的形状に縛られることなく作成することができる。このため、３次元形状データの修正作業に対する規制も大幅に緩和される。よって、VertexとConnectivityとを持たない計算用データモデルは、VertexとConnectivityとを有する計算用データモデルと比較して遥かに容易に作成することが可能となる。したがって、本発明によれば、計算用データモデルの作成における作業負担を軽減することが可能となる。

　また、本発明を用いる場合であっても、プリ処理においては、VertexとConnectivityとを使用しても構わない。つまり、プリ処理においてVertexとConnectivityを用いて分割領域の体積や境界面特性値等を算出しても良い。このような場合であっても、ソルバ処理においては分割領域の体積や境界面特性値があれば物理量の計算が可能であるため、プリ処理においてVertexとConnectivityとを利用するにしても、分割領域の幾何学的形状に対する制約、例えば分割領域の歪みや捩じれ等に起因する制約がなく、計算用データモデルの作成における作業負担の軽減が可能である。

　また、本発明を用いることによって、プリ処理において、分割領域の幾何学的形状に対する制約がなくなるため、分割領域を任意の形状に変更させることができる。このため、分割領域の数を増やすことなく、解析領域を実際に解析したい領域に容易にフィッティングさせることが可能となり、計算負荷を増大させることなく解析精度を向上させることが可能となる。さらに、本発明を用いることによって分割領域の分布密度も任意に変更することができるため、必要な範囲で計算負荷の増大を許容しながらさらに解析精度を向上させることも可能である。

　また、本発明を用いることによって、従来の数値解析手法と異なり、ソルバ処理においてVertexとConnectivityとを用いて分割領域の体積及び境界面特性量を算出する必要がない。したがって、ソルバ処理における計算負荷を低減させることも可能となる。

　また、本発明では、解析領域の形状が変化しない場合には分割領域の移動を必要としないため、粒子法のように時刻が変化するごとに近傍探索処理を行う必要がなく、計算負荷が小さい。また、後に詳説するが、本発明を用いることによって、粒子法と異なり、物理量の保存則を満たしながら物理量の計算を行うことができる。

　一方、従来の数値解析手法である有限体積法は、プリ処理にて分割領域の幾何学的形状を示すVertexとConnectivityとを有する計算用データモデルを作成し、ソルバ処理にて計算用データモデルに含まれるVertexとConnectivityとを用いて分割領域の体積と境界面特性量（境界面の面積及び境界面の法線ベクトル）とを算出してから各分割領域における物理量を計算する。この場合には、幾何学的形状に対する制約、すなわちVertexとConnectivityとの関係に問題ないことが要求される。このため、分割領域の歪み、捩じれ等の制約の中で計算用データモデル（すなわちメッシュ）を作成する必要があり、前述のように計算用データモデル作成上における膨大な手作業が発生するという問題がある。

　また、有限要素法も、ソルバ処理にて計算用データモデルに含まれるVertexとConnectivityとを用いて物理量を計算するため、プリ処理にて分割領域の幾何学的形状を示すVertexとConnectivityとを有する計算用データモデルを作成する必要があり、計算用データモデル作成上における膨大な手作業が発生する。

　また、従来の数値解析手法であるボクセル法は、図３に示すように、ソルバ処理において物理量を算出するにあたり、VertexとConnectivityとが必要ないものの、分割領域の形状がボクセルに限定されることから、前述したように外部領域との境界が階段状になるという問題がある。このため、前述のように十分な解析精度を得ようとすると、ボクセルの数も増え、ソルバ処理における計算負荷が増え問題となる。また、境界補正を行うボクセル法では、結局、分割領域の体積等を算出するにあたりVertexを利用することとなり、計算用データモデルを作成するにあたり分割領域の幾何学的形状の影響を受けることになる。

　また、従来の数値解析手法である粒子法は、分割領域という概念が存在しないため、図３に示すように、ソルバ処理において物理量を算出するにあたり、VertexとConnectivityとが必要ないものの、分割領域の代わりに計算用データモデルを定義する粒子の移動に起因し、前述のように計算負荷が増大する。また、粒子法では、保存則を満たしながら物理量の計算を行うことは難しい。

　続いて、図４を参照して、本発明と従来の有限体積法とについてさらに詳しく比較を行う。従来の有限体積法では、前述のように、プリ処理において、メッシュ分割で得られる分割領域の幾何学的な形状を定義するためのVertexとConnectivityとを有する計算用データモデルを作成する。また、一般的には、ソルバ処理において分割領域の係合情報（以下linkと呼ぶ）が必要となる。このため、プリ処理においては、VertexとConnectivityとlinkとを有する計算用データモデルが作成される。

　そして、従来の有限体積法では、図４に示すように、プリ処理から、VertexとConnectivityとlinkとを有する計算用データモデルと、ソルバ処理において必要となる境界条件や初期条件等とをソルバ処理に受け渡す。ソルバ処理では、受け渡された計算用データモデルに含まれるVertexやConnectivity等を使用して離散化支配方程式を解くことによって物理量の算出を行う。

　一方、本発明では、プリ処理において、任意に配置された分割領域の体積、境界面特性量（境界面の面積、境界面の法線ベクトル）及びlinkを有する計算用データモデルを作成する。また、後に詳説するが、本発明では、必要に応じて、分割領域の内部に配置されるコントロールポイントの座標を計算用データモデルに持たせる場合もある。

　そして、本発明では、図４に示すように、プリ処理から、分割領域の体積、境界面特性量及びlink（必要に応じてコントロールポイントの座標）を有する計算用データモデルと、境界条件や初期条件等とをソルバ処理に受け渡す。ソルバ処理では、受け渡された計算用データモデルに含まれる分割領域の体積、境界面特性量等を使用して離散化支配方程式を解くことによって物理量の算出を行う。

　そして、図４から分かるように、本発明では、ソルバ処理において、VertexとConnectivityとを使用しないで物理量を計算している点が従来の有限体積法と大きく異なり、この点が本発明の大きな特徴である。このような特徴は、ソルバ処理にて、VertexとConnectivityとを必要としない量のみを使用する離散化支配方程式を用いることによって得られるものである。

　この結果、図４に示すように、本発明では、ソルバ処理にVertexとConnectivityとを受け渡す必要がなくなり、プリ処理において、VertexとConnectivityとを持たない計算用データモデルを作成すれば良いこととなる。したがって、従来の有限体積法と比較して本発明では、遥かに容易に計算用データモデルを作成することが可能となり、計算用データモデルの作成における作業負担を軽減することが可能となる。

　また、数値解析を行う解析領域の形状が時系列的に変化する場合、すなわち解析領域が移動境界を含む場合がある。このような場合、移動境界に合わせて分割領域を移動及び変形させる必要がある。

　従来の有限体積法では、移動境界の移動ごとにおけるVertexを予めストアする方法か、分割領域のいびつな変形により計算不可となった場合に領域分割を再実行する方法によって、移動境界を含む場合の物理量の計算を行っている。これに対して、本発明では、Vertexに替えて分割領域の体積や境界面特性値等を予め求めておいてストアしておく方法か、領域分割の再実行によって移動境界を含む場合の物理量の計算を行うことができる。

　従来の有限体積法または本発明において、前述のいずれの方法を採用する場合であっても、複数の計算用データモデルを作成する必要がある。しかしながら、従来の有限体積法では、１つでも膨大な作業量が必要となる計算用データモデルを複数作成するとなると、その作業量は、多くの場合において現実的に負担可能な範囲で行えるものではなくなる。

　一方、本発明では、計算用データモデルがVertexとConnectivityとを持つ必要がなく、領域分割過程でVertexやConnectivityの整合性を考慮する必要がないため、とても高速に計算用データモデルを作成することができ、移動境界を含む場合の物理量の計算を容易に行うことができる。

　ここで、前述のlinkについて補足をする。linkは、物理量のやり取りを行う分割領域同士を関連付ける情報である。そして、このlinkによって関連付けられる分割領域は、必ずしも空間的に隣接されている必要はなく、空間的に離間していても構わない。このようなlinkは、VertexやConnectivityと関連するものではなく、Vertex及びConnectivityと比較すると、極めて短時間で作成することができるものである。

　次に、本発明を用いた数値解析手法（以下、本数値解析手法と称する）の原理、すなわち重み付き残差積分法に基づいて導出された離散化支配方程式と、分割領域の体積と境界面特性量とによって物理量を算出可能となる原理について詳細に説明する。なお、以下の説明において、［］にて挟まれた文字は、図面において太字で記されたベクトルを示す。

　まず、本数値解析手法における計算用データモデルは、解析領域を分割して得られる各分割領域の体積と、隣り合う分割領域同士の境界面の特性を示す境界面特性量とを用いて定義されている。図５は、このような本数値解析手法の計算用データモデルの一例を示す概念図である。この図において、セルＲ１，Ｒ２，Ｒ３・・・は、解析領域を分割して得られる分割領域であり、各々が体積Ｖａ，Ｖｂ，Ｖｃ・・・を有している。また、境界面Ｅは、セルＲ１とセルＲ２との間において物理量の交換が行われる面であり、本発明における境界面に相当するものである。また、面積Ｓａｂは、境界面Ｅの面積を示し、本発明における境界面特性量の１つである。また、［ｎ］ａｂは、境界面Ｅの法線ベクトルを示し、本発明における境界面特性量の１つである。

　また、コントロールポイントａ，ｂ，ｃ・・・は、各セルＲ１，Ｒ２，Ｒ３の内部に配置されおり、図５においては各セルＲ１，Ｒ２，Ｒ３・・・の重心位置に配置されている。ただし、コントロールポイントａ，ｂ，ｃ・・・は、必ずしも各セルＲ１，Ｒ２，Ｒ３・・・の重心位置に配置される必要はない。また、αは、コントロールポイントａからコントロールポイントｂまでの距離を１とした場合におけるコントロールポイントａから境界面Ｅまでの距離を示し、境界面Ｅがコントロールポイントａとコントロールポイントｂとを結ぶ線分のどの内分点に存在するかを示す比率である。

　なお、境界面は、セルＲ１とセルＲ２との間のみに限らず、隣り合う全てのセル間に存在する。そして、境界面の法線ベクトル及び境界面の面積も、境界面ごとに与えられるものである。

　そして、実際の計算用データモデルは、各コントロールポイントａ，ｂ，ｃ・・・の配置データと、各コントロールポイントａ，ｂ，ｃ・・・が存在するセルＲ１，Ｒ２，Ｒ３・・・の体積Ｖａ，Ｖｂ，Ｖｃ・・・を示す体積データと、各境界面の面積を示す面積データと、各境界面の法線ベクトルを示す法線ベクトルデータとを有するデータ群として構築されている。すなわち、本数値解析手法の計算用データモデルは、セルＲ１，Ｒ２，Ｒ３・・・の体積Ｖａ，Ｖｂ，Ｖｃ・・・と、隣り合うセルＲ１，Ｒ２，Ｒ３・・・同士の境界面の特性を示す境界面特性量である境界面の面積と、隣り合うセルＲ１，Ｒ２，Ｒ３・・・同士の境界面の特性を示す境界面特性量である境界面の法線ベクトルとを有して定義されている。

　なお、各セルＲ１，Ｒ２，Ｒ３・・・は、コントロールポイントａ，ｂ，ｃ・・・を有している。このため、セルＲ１，Ｒ２，Ｒ３・・・の体積Ｖａ，Ｖｂ，Ｖｃ・・・は、コントロールポイントａ，ｂ，ｃ・・・が仮想的に占める空間（コントロールボリューム）の体積として捉えることができる。また、本数値解析手法の計算用データモデルは、必要に応じて、境界面が挟まれたコントロールポイント同士を結ぶ線分のどの内分点に存在するかの比率αを示す比率データを有する。

　以下では、前述の計算用データモデルを用いて解析領域の各セル（分割領域）における流速を求める物理量計算例について説明する。なお、ここでは、各コントロールポイントにおける流速を各セルにおける流速として求める。
　まず、本物理量計算において本数値解析手法は、流体解析の場合、下式（１）で示すナビエ・ストークスの式と、下式（２）で示す連続の式とを用いる。

　なお、式（１），（２）において、ｔが時間を示し、ｘｉ（ｉ＝１，２，３）がカーテシアン系における座標を示し、ρが流体密度を示し、ｕｉ（ｉ＝１，２，３）が流体の流速成分を示し、Ｐが圧力を示し、μが流体の粘性係数を示し、添字ｉ（ｉ＝１，２，３），ｊ（ｊ＝１，２，３）がカーテシアン座標系における各方向成分を示している。また、添字ｊに関しては総和規約に従うものとする。

　そして、式（１），（２）を、重み付き残差積分法に基づいて、コントロールボリュームの体積に対して積分して示すと、式（１）が下式（３）のように示され、式（２）が下式（４）のように示される。

　なお、式（３），（４）において、Ｖがコントロールボリュームの体積を示し、∫ＶｄＶが体積Ｖに関する積分を示し、Ｓがコントロールボリュームの面積を示し、∫ＳｄＳが面積Ｓに関する積分を示し、［ｎ］がＳの法線ベクトルを示し、ｎｉ（ｉ＝１，２，３）が法線ベクトル［ｎ］の成分を示し、∂／∂ｎが法線方向微分を示している。

　ここで、説明を簡単化するために、流体の密度ρと粘性係数μを定数とする。ただし、以下の定数化は、流体の物性値が、時間、空間、温度等によって変化する場合に対して拡張可能である。そして、図５のコントロールポイントａについて、境界面Ｅの面積Ｓａｂについて離散化し、代数方程式による近似式に変換すると、式（３）が下式（５）、式（４）が下式（６）のように示される。

　ここで、添字ａｂが付く、［ｎ］ａｂ，［ｕ］ａｂ，ｕｉａｂ，ｎｉａｂ，Ｐａｂ，（∂ｕｉ／∂ｎ）ａｂは、コントロールポイントａとコントロールポイントｂとの間の境界面Ｅ上における物理量であることを示す。また、ｎｉａｂは［ｎ］ａｂの成分である。また、ｍは、コントロールポイントａと結合関係（境界面を挟む関係）にある全てのコントロールポイントの数である。

　そして、式（５），（６）をＶａ（コントロールポイントａのコントロールボリュームの体積）で割ると、式（５）が下式（７）のように示され、式（６）が下式（８）のように示される。

ここで、下式（９）とする。

　すると、式（７）が下式（１０）のように示され、式（８）が下式（１１）のように示される。

　式（１０），（１１）において、［ｕ］ａｂ、ｕｉａｂ、Ｐａｂ、（∂ｕｉ／∂ｎ）ａｂは、コントロールポイントａとコントロールポイントｂ上の物理量の重み付け平均（移流項については、風上性を考慮した重み付け平均）により近似的に求められ、コントロールポイントａ，ｂ間の距離及び向きと、その間に存在する境界面Ｅとの位置関係（上記比率α）と、境界面Ｅの法線ベクトルの向きに依存して決定される。ただし、［ｕ］ａｂ、ｕｉａｂ、Ｐａｂ、（∂ｕｉ／∂ｎ）ａｂは、境界面Ｅの幾何学的な形状には無関係な量（すなわちセル形状を規定するVertexとConnectivityを必要としない量）である。

　また、式（９）で定義されるφａｂも（面積／体積）という量であり、コントロールボリュームの幾何学的形状には無関係な量（すなわちセル形状を規定するVertexとConnectivityを必要としない量）である。つまり、このような式（１０），（１１）は、セル形状を規定するVertexとConnectivityを必要としない量のみを使用して物理量が算出可能な、重み付き残差積分法に基づく演算式である。

　このため、物理量計算（ソルバ処理）に先立って前述の計算用データモデルを作成し、物理量計算において当該計算用データモデルと、式（１０），（１１）の離散化支配方程式とを用いることによって、物理量計算においてコントロールボリュームの幾何学的形状（すなわちセルの形状を規定するVertexとConnectivity）を全く使用せずに、流速の計算を行うことが可能である。

　このように、物理量計算においてVertexとConnectivityとを全く使用せずに流速の計算が行えることから、計算用データモデルにVertexとConnectivityとを持たせる必要がなくなる。よって、計算用データモデルの作成にあたり、セルの幾何学的形状に縛られる必要がなくなるため、セルの形状を任意に設定することができる。このため、本数値解析手法によれば、前述のように３次元形状データの修正作業に対する規制を大幅に緩和することができる。

　なお、実際に式（１０），（１１）を解くにあたり、［ｕ］ａｂやＰａｂ等の境界面Ｅ上の物理量は、通常、線形補間によって補間される。例えば、コントロールポイントａの物理量をψａ、コントロールポイントｂの物理用をψｂとすると、境界面Ｅ上の物理量ψａｂは、下式（１２）によって求めることができる。

　また、物理量ψａｂは、境界面が挟まれたコントロールポイント同士を結ぶ線分のどの内分点に存在するかの比率αを用いることによって、下式（１３）によって求めることもできる。

　したがって、計算用データモデルが比率αを示す比率データを有している場合には、式（１３）を用いて境界面Ｅ上の物理量をコントロールポイントａとコントロールポイントｂとからの離間する距離に応じた重み付け平均を用いて算出することができる。

　また、連続体モデルの方程式（ナビエ・ストークスの式等）には、式（１）に示すように、１階の偏導関数（偏微分）が含まれる。
　ここで、連続体モデルの方程式の微係数を部分積分、ガウスの発散定理、あるいは一般化されたグリーンの定理を利用して、体積分を面積分に変換し、微分の次数を下げる。これによって１次微分は０次微分（スカラー量またはベクトル量）とすることができる。例えば一般化されたグリーンの定理では、物理量をψとすると、下式（１４）という関係が成り立つ。

　なお、式（１４）において、ｎｉ（ｉ＝１，２，３）は、表面Ｓ上の単位法線ベクトル［ｎ］のｉ方向の成分である。連続体モデルの方程式の１次微分項は、体積分から面積分の変換により、境界面上ではスカラー量またはベクトル量として取り扱われる。そして、これらの値は、前述の線形補間等によって、各コントロールポイント上の物理量から補間することができる。

　また、後述するが、連続体モデルの方程式によっては、２階の偏導関数が含まれる場合もある。
　式（１４）の被積分関数をさらに１階微分した式は下式（１５）となり、連続体モデルの方程式の２次微分項は、体積分から面積分の変換により境界面Ｅ上では下式（１６）となる。

　なお、式（１５）において∂／∂ｎは法線方向微分を示し、式（１６）において∂／∂ｎａｂは、［ｎ］ａｂ方向微分を示す。つまり、連続体モデルの方程式の２次微分項は、体積分から面積分の変換により、物理量ψの法線方向微分（Ｓａｂの法線［ｎ］ａｂ方向への微分）に、［ｎ］の成分ｎｉａｂ、ｎｊａｂを乗じた形となる。

　ここで式（１６）中の∂ψ／∂ｎａｂは、下式（１７）と近似される。

　なお、コントロールポイントａとコントロールポイントｂとのコントロールポイント間ベクトル［ｒ］ａｂは、コントロールポイントａの位置ベクトル［ｒ］ａとコントロールポイントｂの位置ベクトル［ｒ］ｂから下式（１８）のように定義される。

　したがって、境界面Ｅの面積がＳａｂであるため、式（１６）は下式（１９）となり、これを利用して式（１６）を計算することができる。

　なお、式（１６）の導出にあたり、次のことがわかる。すべての線形偏微分方程式は、定数と、１次、２次、その他の偏導関数に係数を乗じた項の線形和で表わされる。式（１５）から式（１８）において、物理量ψをψの１次偏導関数に置き換えると、より高次の偏導関数の体積分を、式（１４）のように低次の偏導関数の面積分により求めることができる。この手順を、低次の偏微分から順次繰り返すと、線形偏微分方程式を構成するすべての項の偏導関数は、コントロールポイントの物理量ψと、式（１２）又は式（１３）で計算される境界面上のψであるψａｂと、式（１８）で定義されるコントロールポイント間ベクトルから求められるコントロールポイント間距離と、式（５）に示される境界面Ｅの面積Ｓab、式（１６）に示される法線ベクトルの成分ｎｉａｂとｎｊａｂから、すべて求めることができる。

　次に非線形の偏微分方程式では、例えば、式（２０）に示す非線形項であるψとψの一次偏導関数とを乗算した項や、一次偏導関数の二乗は、反復計算により数値計算することができる。すなわち、それぞれの項のψ、一次偏導関数の方を反復計算における一反復前の計算値として、近似して、反復計算をすればよい。このような方法により、偏微分方程式における非線形項はすべて数値計算が可能になる。以上から、前述では特に連続体モデルの方程式について説明したが、それ以外のいかなる偏微分方程式に対しても、VertexとConnectivityとを必要としない離散化が可能であることがわかった。ただし、保存則については別の成立条件が必要である。これについては後述する。

　さて、本数値解析手法において物理量計算にあたりVertexとConnectivityとを必要としないことは前述した。このため、計算用データモデルの作成（プリ処理）にあたり、コントロールボリュームの体積と、境界面の面積及び法線ベクトルとを、VertexとConnectivityと使用しないで求めれば、式（１０），（１１）の離散化支配方程式を用いて、コントロールボリュームの幾何学的形状（すなわちセルの幾何学的形状）を全く使用せずに、流速の計算を行うことが可能である。

　ただし、本数値解析手法においては、必ずしも、コントロールボリュームの体積と、境界面の面積及び法線ベクトルとを、コントロールボリュームの具体的な幾何学的形状を使用しないで求める必要はない。すなわち、ソルバ処理においてVertexとConnectivityを利用しないので、コントロールボリュームの具体的な幾何学的形状、具体的にはVertexとConnectivityを利用するとしても、従来の有限要素法、有限体積法のような分割領域に関わる制約、すなわち分割領域の歪みや捩じれに対する制約がないため、前述のように容易に計算用データモデルの作成が可能となる。

　さて、本数値解析手法においては、条件によっては、前述の法線ベクトルを、コントロールボリューム同士を結ぶ距離ベクトルに置き換えることも可能である。以下にその理由について説明する。図５において示す境界面Ｅの法線ベクトル［ｎ］ａｂが、コントロールポイントａとコントロールポイントｂとを結ぶ距離ベクトル［ｒ］ａｂと同じ方向を向く場合には、法線ベクトル［ｎ］は、下式（２１）のように示すことができる。

　したがって、境界面Ｅの法線ベクトル［ｎ］ａｂが距離ベクトル［ｒ］ａｂと同じ方向を向く場合には、式（１０），（１１）で示す離散化支配方程式に式（２１）を代入することによって、法線ベクトル［ｎ］ａｂがコントロールポイントａ，ｂ間方向を向いている場合の離散化支配方程式が得られる。すなわち、境界面と該境界面を挟むコントロールポイントを結ぶベクトルが直交する場合には、離散化支配方程式における法線ベクトルを距離ベクトルに置き換えることができる。このような離散化支配方程式によれば、コントロールポイントの位置座標だけから、境界面Ｅの法線ベクトル［ｎ］ａｂを決定することができる。また、境界面Ｅと距離ベクトルとを出来るだけ直交に近づけることによって物理量を計算する際の精度が向上する。したがって、法線ベクトルを、距離ベクトルに置き換えることによって、計算精度の向上を図ることができる。さらに、法線ベクトルの任意性を、コントロールポイントを結ぶ方向に固定化することができる。

　ただ、法線ベクトルの任意性がなくなることで境界面の姿勢に自由度を持たせることができなくなる場合がある。このような場合には、法線ベクトルの任意性が存在する場合と比較して、計算用データモデルを作成する際にコントロールボリュームの体積と境界面の設定に制約が生じることとなる。また、法線ベクトルを、コントロールポイントを結ぶベクトルに置き換える場合には、距離ベクトルが必要となり、該距離ベクトルを算出するためにコントロールポイントの座標が必要となるため、計算用データモデルにコントロールポイントの座標データあるいは距離ベクトルデータを持たせる必要がある。ただし、この場合でも、本数値解析手法においては、物理量計算にてVertexとConnectivityとは必要ない。

　次に、本発明と粒子法との違いを明確にするために、物理量計算にあたり、物理量の保存則が満足される条件について説明する。

　例えば、図６に示すように、Ｌ字形の流路を４つに分割するセルＲａ，Ｒｂ，Ｒｃ，Ｒｄについて考える。なお、各セルＲａ，Ｒｂ，Ｒｃ，Ｒｄの内部に配置されるコントロールポイントａ，ｂ，ｃ，ｄは、セルＲａ，Ｒｂ，Ｒｃ，Ｒｄの中心に設置されている。また、各境界面における流速ベクトルは、境界面に対して垂直であるものとする。

　なお、図６において、ＶａがセルＲａの体積（コントロールポイントａのコントロールボリュームの体積）、ＶｂがセルＲｂの体積（コントロールポイントｂのコントロールボリュームの体積）、ＶｃがセルＲｃの体積（コントロールポイントｃのコントロールボリュームの体積）、ＶｄがセルＲｄの体積（コントロールポイントｄのコントロールボリュームの体積）、ρａがセルＲａにおける密度、ρｂがセルＲｂにおける密度、ρｃがセルＲｃにおける密度、ρｄがセルＲｄにおける密度、ＳａがセルＲａと外部領域との境界面の面積、ＳｃがセルＲｃと外部領域との境界面の面積、ＳｄがセルＲｄと外部領域との境界面の面積、ＳａｂがセルＲａとセルＲｂとの境界面の面積、ＳａｃがセルＲａとセルＲｃとの境界面の面積、ＳｂｄがセルＲｂとセルＲdとの境界面の面積、ＳｃｄがセルＲｃとセルＲｄとの境界面の面積、ｕａがセルＲａと外部領域との境界面における流速、ｕｃがセルＲｃと外部領域との境界面における流速、ｕｄがセルＲｄと外部領域との境界面における流速、ｕａｂがセルＲａとセルＲｂとの境界面における流速、ｕａｃがセルＲａとセルＲｃとの境界面における流速、ｕｂｄがセルＲｂとセルＲdとの境界面における流速、ｕｃｄがセルＲｃとセルＲｄとの境界面における流速、ρａがセルＲａにおける密度、ρabがセルＲａとセルＲｂとの境界面における密度、ρacがセルＲａとセルＲcとの境界面における密度、ρbdがセルＲbとセルＲdとの境界面における密度、を示している。

　図６に示すような４つのコントロールポイントａ，ｂ，ｃ，ｄ（４つのセルＲａ，Ｒｂ，Ｒｃ，Ｒｄ）上で質量保存の式を離散化することを考える。なお、質量保存の式を離散化した離散化支配方程式は、後述の式（４３）に示す。
　コントロールポイントａにおける離散化支配方程式は、下式（２２）に示される。

　コントロールポイントｂにおける離散化支配方程式は、下式（２３）に示される。

　コントロールポイントｃにおける離散化支配方程式は、下式（２４）に示される。

　コントロールポイントｄにおける離散化支配方程式は、下式（２５）に示される。

　そして、式（２２）～（２５）を全て足し加えると、下式（２６）を得る。

　各コントロールポイントａ，ｂ，ｃ，ｄのコントロールボリュームの体積Ｖａ，Ｖｂ，Ｖｃ，Ｖｄが時間について一定であることを利用して時間微分項の中にくり入れると、式（２６）は、下式（２７）となる。

　今、各コントロールポイントａ，ｂ，ｃ，ｄのコントロールボリュームが占める全体積をＶａｂｃｄ、平均密度ρ（オーバーライン）ａｂｃｄとすると、全体積をＶａｂｃｄが下式（２８）で示され、平均密度ρ（オーバーライン）ａｂｃｄが下式（２９）で示される。

　したがって式（２７）は、下式（３０）のように示される。

　式（３０）は、コントロールポイントａ，ｂ，ｃ，ｄのコントロールボリュームが占める全領域に流入する質量流束と流出する質量流束の差が単位時間にコントロールポイントａ，ｂ，ｃ，ｄのコントロールボリュームが占める全領域の平均密度の時間変化（質量の時間変化）に等しいことを意味する。つまり、各コントロールポイントａ，ｂ，ｃ，ｄごとに離散化した質量保存の式は、全コントロールポイントのコントロールボリュームが占める領域に対しても成立する。すなわち、コントロールポイントが示すコントロールボリューム領域についての離散化支配方程式は、全てのコントロールポイントについて足し加えると、計算対象である解析領域の全領域に関する保存則を満足する方程式にならなくてはならない。

　続いて、コントロールポイントの全数をＮとし、式（４３）に示される質量保存の式を足し加えると、下式（３１）が得られる。

　式（３１）において、各コントロールポイント間の境界面の面積は、コントロールポイントａ側から見てもコントロールポイントｂ側から見ても等しいとすると、各コントロールポイント間の質量流束（ρ［ｎ］・［ｕ］）・Ｓは、コントロールポイントａ側とコントロールポイントｂ側で正負が逆で絶対値が等しくなるので差し引きゼロとなりキャンセルされる。つまり、式（３１）は、計算する全領域に対して、流入する質量と流出する質量との差が全領域での質量の単位時間変化に等しいことを示す。したがって、式（３１）は、解析領域全体における質量保存の式となる。よって、式（３１）が、計算する全領域についての質量保存則を満足するためには、２つのコントロールポイント間の境界面の面積と一致するという条件及び法線ベクトルが一方のコントロールポイント側から見た場合と他方のコントロールポイント側から見た場合とで絶対値が一致するという条件が必要である。

　また、質量保存則を満足するためには、下式（３２）に示される全コントロールポイントのコントロールボリュームの占める体積が解析領域の全体積と一致するという条件が必要である。

　これは、連続体の密度ρが、ρ１＝ρ２＝・・・＝ρと１つの変数で表されることを考えれば容易に理解される。
　このように質量保存則を満足するためには、全コントロールポイントのコントロールボリュームの体積の総和が解析領域の体積と一致するという条件が必要である。

　なお、ここでは、質量保存の式に対して説明を行ったが、保存則は、連続体の運動量やエネルギに対しても成立しなければならない。これらの物理量に対しても、後述の式（５０）、（５５）において全コントロールポイントに対して足し加えることによって、保存側が満足されるためには、全コントロールポイントのコントロールボリュームの占める体積が解析領域の全体積と一致するという条件と、２つのコントロールポイント間の境界面の面積が一致する条件及び法線ベクトルが一方のコントロールポイント側から見た場合と他方のコントロールポイント側から見た場合とで絶対値が一致する（正負逆符号）という条件とが必要であることが分かる。

　また、保存則を満たすためには、図７に示すように、コントロールポイントａの占めるコントロールボリュームを考えた場合に、コントロールポイントａを通り、任意の向きの単位法線ベクトル［ｎ］ｐを持つ無限に広い投影面Ｐを考えたときに下式（３３）が成り立つという条件が必要である。

　なお、図７及び式（３３）において、Ｓｉが境界面Ｅｉの面積、［ｎ］ｉが境界面Ｅｉの単位法線ベクトル、ｍがコントロールボリュームの面の総数を示す。
　式（３３）は、コントロールボリュームを構成する多面体が、閉包空間を構成することを示す。この式（３３）は、コントロールボリュームを構成する多面体の一部が凹んでいる場合であっても成立する。
　なお、図８に示すように、２次元における三角形についても式（３３）が成り立つ。また、多面体の１つの面を微小面ｄＳとし、ｍを∞とする極限を取ると、下式（３４）となり、図９に示すような閉包曲面体についても成り立つことが分かる。

　式（３３）が成り立つという条件は、ガウスの発散定理や、式（１４）に示す一般化されたグリーンの定理が成り立つために必要な条件である。
　そして、一般化されたグリーンの定理は、連続体の離散化のための基本となる定理である。したがって、グリーンの定理にしたがって体積分を面積分に変化し離散化させる場合において、保存則を満足させるためには式（３３）が成り立つという条件は必須となる。

　このように、前述の計算用データモデル及び物理量計算を用いて数値解析を行う際に、物理量の保存則が満足されるためには、以下の３つの条件が必要となる。
　（ａ）全コントロールポイントのコントロールボリュームの体積（全分割領域の体積）の総和が解析領域の体積と一致する。
　（ｂ）２つのコントロールポイント間の境界面の面積が一致する及び法線ベクトルが一方のコントロールポイント側（境界面を挟む一方の分割領域）から見た場合と他方のコントロールポイント側（境界面を挟む他方の分割領域）から見た場合とで絶対値が一致する。
　（ｃ）コントロールポイントを通り（分割領域を通り）、任意の向きの単位法線ベクトル［ｎ］ｐを持つ無限に広い投影面Ｐを考えたときに式（３３）が成り立つ。

　つまり、保存則を満足させる場合には、これらの条件を満足するように計算用データモデルを作成する必要がある。ただし、前述のように本数値解析手法においては、計算用データモデルの作成にあたり、セル形状を任意に変形することができることから、容易に上記３つの条件を満足するように計算用データモデルを作成することができる。

　次に、従来の粒子法であるＭＰＳ（Moving Particle Semi-Implicit）法が保存則を満足できない理由、計算負荷が増大する理由、及び本数値解析手法の粒子法に対する優位性について詳説する。

　ＭＰＳ法は、適切に設定された半径ｒｅの球内に存在する粒子を検出し、それらと結合関係を結んで計算する手法であり、例えば図１０に示すように、粒子ｉ周りに複数の粒子ｊが存在する場合に、粒子ｉにおけるラプラシアン（∇２ψ）ｉは、下式（３５）のように近似される。

　なお、図１０において、ψｉが粒子ｉにおける物理量を示し、ψｊが粒子ｊにおける物理量を示し、［ｒ］ｉｊが粒子ｉから粒子ｊへの距離ベクトルを示している。
　また、式（３５）におけるｄは次元数を示す定数であり、３次元である場合には３となる。また式（３５）におけるω（ｒ）は重み関数であり、下式（３６）のように示される。また、式（３５）におけるｍは、結合関係にある粒子の数を示す。

　一方、図１１に示すように、粒子ｉ，ｊをコントロールポイントとして捉え、粒子ｉのコントロールボリュームがＶｉ、粒子ｉと粒子ｊとの間の境界面の面積がＳｉｊ、粒子ｉと粒子ｊとの間の境界面の法線ベクトルが［ｎ］ｉｊ、粒子ｉから粒子ｊへの距離ベクトルが［ｒ］ｉｊである場合には、粒子ｉにおけるラプラシアン（∇２ψ）ｉは、下式（３７）のように近似される。

　そして、式（３６）と式（３７）とを比較すると、仮に［ｒ］ｉｊと［ｎ］ｉｊが同じ向きであるとした場合には、下式（３８），（３９）が得られる。

　式（３９）は、左辺と右辺の次元が（１／距離）で一致している。したがって、右辺のＭＰＳ法定式化は、２つの粒子ｉ，ｊ間の距離のみから、下式（４０）に示す（面積／体積）という量（すなわち式（９）で定義された比率）を算出する手法と解釈出来る。

　しかし、ｍ個の粒子の結合関係のみから、境界面の面積Ｓｉｊと体積Ｖｉを求めるには、関係式が不足しており、具体的な数値は決定できず、あくまでも、式（４０）の比率のみが決定される。このため、仮に、ＭＰＳ法の離散化支配方程式から、境界面の面積Ｓｉｊと体積Ｖｉを求めたとしても、前述した保存則が満足される条件（ａ）～（ｃ）を満足する保証は全くない。このことは、ＭＰＳ法が、保存則の満足性という点で大きな問題を有していることを示している。

　数値解析を工学上の問題、特に機械設計問題やプラント設計問題に適用した場合において、定量値（圧力、温度、熱量等）の評価は極めて重要であるが、数値解析が保存則を満足しない場合には、定量性が保証されなくなる。つまり、ＭＰＳ法では、質量保存、運動量保存、エネルギ保存という保存則が満足される保証がなく、定量性が保証されない。これに対して、本数値解析手法によれば、保存則を満足させることができ、定量性を保証することができる。

　なお、ＭＰＳ法は、前述のように、時刻の変化に伴って粒子が移動するため、例えば、その都度、前述の半径ｒｅの球内に存在する粒子を検出する近傍探索処理を行う必要があり、物理量計算における計算負荷が大きくなる。これに対して、本数値解析手法においてコントロールボリューム及びコントロールポイントは、時刻が変化した場合であっても移動しない。このため、予めコントロールボリュームやコントロールポイントの配置関係が分かっている場合は近傍探索処理を行うことなく物理量計算を行うことができる。このため、物理量計算における計算負荷をＭＰＳ法と比較して小さくすることができる。なお、コントロールボリューム及びコントロールポイントの配置関係が予め分かっていない場合であっても、最初に一度のみ、コントロールボリュームやコントロールポイントの配置関係を定める処理を行えば良い。

　なお、前述の説明においては、ナビエ・ストークスの式及び連続の式から重み付き残差積分法に基づいて導出した離散化支配方程式を用いる物理量の計算例について説明したが、本数値解析手法において用いられる離散化支配方程式はこれに限られるものではない。つまり、種々の方程式（質量保存の方程式、運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式、及び波動方程式等）から重み付き残差積分法に基づいて導出されると共に、VertexとConnectivityとを必要としない量のみを使用して物理量を算出可能な離散化支配方程式であれば本数値解析手法に用いることができる。

　そして、このような離散化支配方程式の特性によって、従来の有限要素法や有限体積法のようにいわゆるメッシュを必要としない、メッシュレスでの計算が可能となる。また、たとえ、プリ処理において、セルの幾何学的形状を規定するVertexとConnectivityとを利用するとしても、従来の有限要素法、有限体積法、ボクセル法のようなメッシュに対する制約がないため、計算用データモデルの作成に伴う作業負荷を低減できる。

　以下に、質量保存の方程式、運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式、及び波動方程式から、重み付き残差積分法に基づいて、VertexとConnectivityとを必要としない量のみを使用する離散化支配方程式が導出可能であること、すなわち本数値解析手法において他の支配方程式を用いることができることについて説明する。

（１）質量保存の方程式
　オイラー座標系における質量保存の方程式は、下式（４１）のように微分形式で示される。

　なお、式（４１）においてｔが時間、ｘｉ（ｉ＝１，２，３）がカーテシアン系における座標、ρが密度、ｕｉ（ｉ＝１，２，３）が変形速度の成分、添字ｉ（ｉ＝１，２，３）がカーテシアン座標系における各方向成分を示している。また、添字ｉに関しては総和規約に従うものとする。
　そして、式（４１）を、重み付き残差積分法に基づいて、コントロールポイントのコントロールボリュームの体積Ｖに対して積分すると下式（４２）と示される。

　さらに、図５に示すコントロールポイントａについて離散化し、代数方程式に変換すると、下式（４３）のように示される。

　ここで、添字ａｂが付く、ρａｂ、［ｎ］ａｂは、コントロールポイントａとコントロールポイントｂとの間の境界面Ｅ上における物理量であることを示す。また、ｍは、コントロールポイントａと結合関係（境界面を挟む関係）にある全てのコントロールポイントの数である。

　そして、式（４３）をコントロールポイントａのコントロールボリュームの体積であるＶａで割ると、下式（４４）が得られ、さらに下式（４５）とすると、質量保存の方程式が離散化された下式（４６）を得る。

　このような式（４６）は、重み付き残差積分法に基づいて導出されると共に、VertexとConnectivityとを必要としない量のみを使用する方程式であるため、本数値解析手法の離散化支配方程式として用いることができる。

　なお、本数値解析手法で用いられる離散化支配方程式は、従来の幾何学的形状を示す量を使用する方程式を重み付き残差積分法に基づいて導出する過程で敢えて途中にて留めることによって得られることは前述した。つまり式（４６）は、質量保存の方程式を、重み付き残差積分法に基づいて、Vertex等を使用する方程式を導出する過程で得られるものである。ここで、図１２は、２次元三角形状のセルを示す模式図である。図１２における三角形aの面積、辺の長さ、法線ベクトルを下表に示す。なお、下表において記号×は、外積を示している。

　図１２に示すように、セルが２次元三角形状である場合において、質量保存の方程式を、重み付き残差積分法に基づいて、Vertex等を使用する離散化支配方程式とすると、下式（４７）となる。

　ただし、式（４７）において、［ｒ］ｉが頂点（Vertex）ｉの位置ベクトル（座標）であり、記号×がベクトルの外積を示している。また、ρａｂ及びρを一定とし、［ｒ］ｉｊを［ｒ］ｊ－［ｒ］ｉとし、［ｒ］４を［ｒ］１とする。

（２）運動量保存の方程式
　オイラー座標系における運動量保存の方程式は、下式（４８）のように微分形式で示される。

　なお、式（４８）においてσｉｊ（ｉ，ｊ＝１，２，３）が連続体内部応力、ｆｉ（ｉ＝１，２，３）が連続体に作用する外力（例えば重力）を示し、他の量は式（４１）と同様である。また、添字ｊに関しては総和規約に従うものとする。
　この式（４７）は、構造体や材料、流体等の応力場の基礎方程式である。

　そして、式（４７）を、重み付き残差積分法に基づいて、コントロールポイントのコントロールボリュームの体積Ｖに対して積分すると下式（４９）と示される。

　さらに、図５に示すコントロールポイントａについて離散化し、代数方程式に変換すると、下式（５０）のように示される。

　式（５０）をコントロールポイントａのコントロールボリュームの体積であるＶａで割り、さらに式（４５）を導入すると、運動量保存の方程式が離散化された下式（５１）が得られる。

　なお、運動量保存の方程式において応力テンソルの対称性を考慮すると、角運動量保存の方程式も、運動量保存の方程式と同様に離散化できる。

（３）移流拡散方程式
　ある物質Ｃの連続体内への移流拡散現象は、下式（５２）の移流拡散方程式で示される。

　なお、式（５２）において、Ｃが物質Ｃの濃度、μＣが物質Ｃの拡散係数、ｑＣが物質Ｃのソース（シンク）項、ρが連続体の密度、ｕｉが連続体の変形速度を示す。
　そして、式（５２）を、重み付き残差積分法に基づいて積分し、さらに離散化して離散化支配方程式に変換すると、下式（５３）を得る。

　なお、Ｃａ等の添字ａが付く量がコントロールポイントａの物理量であり、Ｃａｂ等の添字ａｂが付く量がコントロールポイントａとコントロールポイントｂとの間の境界面における物理量である。

（４）エネルギ保存の方程式
　エネルギ保存則は、熱エネルギ保存と運動エネルギ保存に分けられるが、運動エネルギの保存は前述の運度量保存に含まれるため、ここでは熱エネルギ保存の方程式の一般形を下式（５４）に示す。

　なお、式（５４）において、Ｕが連続体の内部エネルギ、ｑｉが熱流束ベクトル、ｒが熱源，熱エネルギのソース（シンク）項、σｉｊが応力テンソル、Ｄｉｊが変形速度テンソルを示す。なお、σｉｊ・Ｄｉｊのテンソル２重積の項は、応力仕事率と呼ばれる。また、添字ｉ，ｊについては、総和規約を適用する。
　そして、式（５４）を、重み付き残差積分法に基づいて積分し、さらに離散化して離散化支配方程式に変換すると、下式（５５）を得る。

　そして、熱流束ベクトル［ｑ］に、熱流束に加えて電気的または化学的エネルギ等の全ての非力学的エネルギの流束を加えると、エネルギ保存の方程式は、非常に広い範囲のエネルギの保存を示す方程式となる。

（５）波動方程式
　前述の質量保存の方程式、運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式等の保存形で表される物理法則の方程式は、「放物型」と「楕円型」と呼ばれる偏微分方程式の性質を合わせ持った方程式である。これに対して、波の伝わりや振動の伝わりを表す波動方程式は「双曲型」と呼ばれ、一般形として下式（５６）と示される。

　なお、式（５６）において、ｕが振幅，変位、αが波の伝播速度を示している。
　そいして、式（５６）を、重み付き残差積分法に基づいて積分し、さらに離散化して離散化支配方程式に変換すると、下式（５７）が得られる。

　式（５６）から、空間方向にはコントロールボリュームの幾何学的形状を必要としない離散化手法がそのまま適用できることが分かる。ただし、時間方向には２階の導関数であるので、高精度な時間積分法が用いられる。

　以上のような式（４６），（５１），（５３），（５５），（５７）は、重み付き残差積分法に基づいて導出されると共に、VertexとConnectivityとを必要としない量のみを使用して物理量を算出可能な方程式であるため、本数値解析手法の離散化支配方程式として用いることができる。そして、本数値解析手法は、上記離散化支配方程式を用いることによって、定常及び非定常の、流体力学、熱伝導、移流拡散、構造力学、波動、及びこれらの物理現象が連成した現象における数値解析に適用することができる。

　次に、本数値解析手法を用いる場合には、異なる座標系で設計された部品を容易に組み上げた状態で数値解析を行うことができる。これは、本数値解析手法においては、有限要素法で用いられる要素と異なり、コントロールボリュームの具体的な幾何学的形状を必要とせず、また２つの隣接するコントロールポイント間の「距離」は、絶対座標系における距離である必要はなく、「計算上の距離」であれば良いためである。
　したがって、本数値解析手法によれば、図１３に示すように、異なる座標系で設計された部品Ａ，Ｂ，Ｃを組み上げた状態で数値解析したい場合において、各部品の座標系を一致させることなく数値解析を行うことができる。

　なお、一般的な有限要素法と呼ばれる数値解析手法は、図１４や図１５に示すような要素の交差を全く許さない手法である。このため、有限要素法において用いられるアプリケーションソフトにて要素の交差が検出されるとエラー出力がなされる。そして、このようなことから、有限要素において計算用データモデルの作成に非常に大きな作業負担がかかることは前述の通りである。

　一方、本数値解析手法によれば、図１４や図１５に示すように、コントロールポイント間の結合に交差を許すことが可能である。これは、各コントロールポイントが占めるコントロールボリュームの体積（分割領域の体積）、境界面の面積、及び境界面の法線ベクトルという情報量は、具体的な幾何学的形状を有している必要がない。したがって、コントロールポイント間の結合が交差する場合であっても物理量計算を実行することが可能である。

　したがって、計算用データモデルを作成する際の制約が減り、計算用データモデルの作成の自由度が大幅に拡大する。ただし、流体などにおける物理現象を考えた場合には、１つのコントロールポイントに近接するコントロールポイントからの情報によって、そのコントロールポイントの物理用がアップデートされるという性質を有しているため、遠方のコントロールポイントと結合させると計算精度が悪化する虞がある。このため、本数値解析手法においても適切な範囲においてコントロールポイント間結合を形成することが好ましい。

　本発明は、数値解析を行う対象物体を複数の多角形で表現したパッチデータが記憶されたパッチデータ記憶手段と、前記数値解析を行うのに必要なパラメータを記憶するパラメータ記憶手段と、前記対象物体を含む解析領域内を複数の直方体に分割したボクセルデータを記憶するボクセルデータ記憶手段と、前記パラメータ記憶手段に記憶されている前記パラメータに基づき、前記ボクセルデータを定義し、各ボクセルに対して、ボクセル属性を付与して、前記ボクセルデータ記憶手段に記憶するボクセルデータ生成手段と、前記解析領域内の領域分割を行うための初期点データを記憶する初期点データ記憶手段と、前記ボクセルの中心点を使用して、前記ボクセルの数より少ない前記初期点データを生成して、前記初期点データ記憶手段に記憶する初期点データ生成手段と、前記対象物体を複数の分割領域に分割した分割領域データを記憶する分割領域データ記憶手段と、前記ボクセルデータ記憶手段に記憶されている各ボクセルの属性と、前記初期点データ記憶手段に記憶されている前記初期点データとに基づき、前記対象物体内を、複数の前記ボクセルで構成する分割領域を定義し、定義した前記分割領域データを前記分割領域データ記憶手段に記憶する分割領域データ生成手段と、前記数値解析を行うための計算用データを記憶する計算用データ記憶手段と、前記分割領域データ記憶手段に記憶されている前記分割領域データに基づき、各分割領域の境界面データを生成し、該境界面データを前記計算用データとして前記計算用データ記憶手段に記憶する計算用データ生成手段とを備えることを特徴とする。

　本発明においては、前記分割領域データ生成手段は、前記初期点データ記憶手段に記憶されている前記初期点データの中から、各ボクセルの中心点に最も距離が近い初期点を選択して、各ボクセルが属する分割領域を特定することにより前記分割領域を生成することを特徴とする。

　本発明においては、前記分割領域データ生成手段は、所定の軸に垂直な前記ボクセルの境界面で構成するスライス面と、前記初期点のそれぞれを中心とする球体の交差する断面を定義し、該断面から前記初期点からの距離に応じて高さの異なるポテンシャル立体を定義し、該ポテンシャル立体を３次元の隠面処理を行って描画した画像データに基づいて、前記スライス面上の分割領域を定義する処理を、前記解析領域内の全ての前記スライス面に対して行うことにより前記分割領域を生成することを特徴とする。

　本発明においては、前記初期点データ生成手段は、前記パッチデータ記憶手段に記憶されている前記パッチデータと、前記ボクセルデータ記憶手段に記憶されている前記ボクセルデータとから、前記対象物体内部に位置する前記ボクセルの内側ボクセルを選択し、前記内側ボクセルのうち、隣接するボクセルが存在しない内側ボクセルに内接する球を定義し、前記球上に球上点を定義するともに、前記対象物体の壁上に壁上点を定義し、対となる前記球上点と前記壁上点との間に分割点を定義し、該分割点と前記ボクセルの中心点とを前記初期点とすることを特徴とする。

　本発明は、数値解析を行う対象物体を複数の多角形で表現したパッチデータが記憶されたパッチデータ記憶手段と、前記数値解析を行うのに必要なパラメータを記憶するパラメータ記憶手段と、前記対象物体を含む解析領域内を複数の直方体に分割したボクセルデータを記憶するボクセルデータ記憶手段と、前記解析領域内の領域分割を行うための初期点データを記憶する初期点データ記憶手段と、前記対象物体を複数の分割領域に分割した分割領域データを記憶する分割領域データ記憶手段と、前記数値解析を行うための計算用データを記憶する計算用データ記憶手段とを備える計算用データ生成装置における計算用データ生成方法であって、前記パラメータ記憶手段に記憶されている前記パラメータに基づき、前記ボクセルデータを定義し、各ボクセルに対して、ボクセル属性を付与して、前記ボクセルデータ記憶手段に記憶するボクセルデータ生成ステップと、前記ボクセルの中心点を使用して、前記ボクセルの数より少ない前記初期点データを生成して、前記初期点データ記憶手段に記憶する初期点データ生成ステップと、前記ボクセルデータ記憶手段に記憶されている各ボクセルの属性と、前記初期点データ記憶手段に記憶されている前記初期点データとに基づき、前記対象物体内を、複数の前記ボクセルで構成する分割領域を定義し、定義した前記分割領域データを前記分割領域データ記憶手段に記憶する分割領域データ生成ステップと、前記分割領域データ記憶手段に記憶されている前記分割領域データに基づき、各分割領域の境界面データを生成し、該境界面データを前記計算用データとして前記計算用データ記憶手段に記憶する計算用データ生成ステップとを有することを特徴とする。

　本発明においては、前記分割領域データ生成ステップは、前記初期点データ記憶手段に記憶されている前記初期点データの中から、各ボクセルの中心点に最も距離が近い初期点を選択して、各ボクセルが属する分割領域を特定することにより前記分割領域を生成することを特徴とする。

　本発明においては、前記分割領域データ生成ステップは、所定の軸に垂直な前記ボクセルの境界面で構成するスライス面と、前記初期点のそれぞれを中心とする球体の交差する断面を定義し、該断面から前記初期点からの距離に応じて高さの異なるポテンシャル立体を定義し、該ポテンシャル立体を３次元の隠面処理を行って描画した画像データに基づいて、前記スライス面上の分割領域を定義する処理を、前記解析領域内の全ての前記スライス面に対して行うことにより前記分割領域を生成することを特徴とする。

　本発明においては、前記初期点データ生成ステップは、前記パッチデータ記憶手段に記憶されている前記パッチデータと、前記ボクセルデータ記憶手段に記憶されている前記ボクセルデータとから、前記対象物体内部に位置する前記ボクセルの内側ボクセルを選択し、前記内側ボクセルのうち、隣接するボクセルが存在しない内側ボクセルに内接する球を定義し、前記球上に球上点を定義するともに、前記対象物体の壁上に壁上点を定義し、対となる前記球上点と前記壁上点との間に分割点を定義し、該分割点と前記ボクセルの中心点とを前記初期点とすることを特徴とする。

　本発明は、数値解析を行う対象物体を複数の多角形で表現したパッチデータが記憶されたパッチデータ記憶手段と、前記数値解析を行うのに必要なパラメータを記憶するパラメータ記憶手段と、前記対象物体を含む解析領域内を複数の直方体に分割したボクセルデータを記憶するボクセルデータ記憶手段と、前記解析領域内の領域分割を行うための初期点データを記憶する初期点データ記憶手段と、前記対象物体を複数の分割領域に分割した分割領域データを記憶する分割領域データ記憶手段と、前記数値解析を行うための計算用データを記憶する計算用データ記憶手段とを備える計算用データ生成装置上のコンピュータに計算用データ生成を行わせる計算用データ生成プログラムであって、前記パラメータ記憶手段に記憶されている前記パラメータに基づき、前記ボクセルデータを定義し、各ボクセルに対して、ボクセル属性を付与して、前記ボクセルデータ記憶手段に記憶するボクセルデータ生成ステップと、前記ボクセルの中心点を使用して、前記ボクセルの数より少ない前記初期点データを生成して、前記初期点データ記憶手段に記憶する初期点データ生成ステップと、前記ボクセルデータ記憶手段に記憶されている各ボクセルの属性と、前記初期点データ記憶手段に記憶されている前記初期点データとに基づき、前記対象物体内を、複数の前記ボクセルで構成する分割領域を定義し、定義した前記分割領域データを前記分割領域データ記憶手段に記憶する分割領域データ生成ステップと、前記分割領域データ記憶手段に記憶されている前記分割領域データに基づき、各分割領域の境界面データを生成し、該境界面データを前記計算用データとして前記計算用データ記憶手段に記憶する計算用データ生成ステップとを前記コンピュータに行わせることを特徴とする。

　本発明によれば、従来の数値解析手法である有限要素法、有限体積法、ボクセル法、ボクセル法の改良手法、及び粒子法の問題点を解決し、解析精度の悪化を伴わずにソルバ処理における計算負荷の低減を図ることができる数値解析装置に入力する計算用データ作成の作業負担を軽減することができるという効果が得られる。

従来の有限要素法による計算用データモデルの一例を示す概念図である。従来の有限体積法による計算用データモデルの一例を示す概念図である。本発明と従来の数値解析手法とを比較するための表である。本発明と従来の有限体積法とを詳細に比較するための図である。本数値解析手法の計算用データモデルの一例を示す概念図である。本数値解析手法において物理量の保存則が満足される条件について説明するための複数の分割領域を示す模式図である。コントロールポイントを通り、任意の向きの単位法線ベクトルを持つ無限に広い投影面を示す模式図である。２次元における三角形のコントロールボリュームを考えた場合において物理量の保存則が満足される条件について説明する模式図である。球のコントロールボリュームを考えた場合において物理量の保存則が満足される条件について説明する模式図である。粒子法による近傍探索について示す模式図である。本数値解析手法におけるコントロールボリュームを示す模式図である。２次元三角形状のセルを示す模式図である。異なる座標系で設計された部品を組み上げた状態を示す模式図である。要素が交差する様子を示す模式図である。コントロールポイント間の結合が交差する様子を示す模式図である。数値解析装置のハードウェア構成を概略的に示すブロック図である。数値解析方法を示すフローチャートである。数値解析方法にて行うプリ処理を示すフローチャートである。数値解析方法にて行うソルバ処理を示すフローチャートである。解析領域が移動境界を含む場合における数値解析方法を示すフローチャートである。第１の実施形態の構成を示すブロック図である。図２１に示すパラメータ記憶部１２のテーブル構造を示す図である。図２１に示す形状データ記憶部１５に記憶される形状データの一例を示す図である。解析領域と解析対象物体の位置関係を示す図である。図２１に示すパッチデータ記憶部１７に記憶されるパッチデータの一例を示す図である。図２１に示すパッチデータ記憶部１７のテーブル構造を示す図である。ボクセルの構成を示す図である。図２１に示すボクセルデータ記憶部１９のテーブル構造を示す図である。各ボクセルを定義するデータを求める方法を示す図である。図２１に示す初期点データ記憶部２１のテーブル構造を示す図である。図２１に示す分割領域データ記憶部２４のテーブル構造を示す図である。図２１に示す中間データ記憶部２６のテーブル構造を示す図である。図２１に示す計算用データ記憶部２７のテーブル構造を示す図である。図２１に示すボクセルデータ生成部１８がボクセルデータを生成する処理動作を示すフローチャートである。分割領域データ生成部２３が分割領域データを生成して、分割領域データ記憶部２４に記憶する処理動作を示すフローチャートである。計算用データ生成部２５が分割領域データから計算用データを生成して、計算用データ記憶部２７に記憶する処理動作を示すフローチャートである。第２の実施形態の構成を示すブロック図である。図３７に示すグラフィックボード２８の機能構成を示す図である。図３７に示すグラフィックボード２８の座標系を示す図である。ポテンシャル立体の構成を示す図である。グラフィックボード２８に入力する三角形データの構成を示す図である。グラフィックボード２８における三角形の描画処理を示す図である。図３７に示す分割領域データ生成部２３１の処理動作を示すフローチャートである。図３７に示すグラフィックボード２８内の三角形描画部の処理動作を示すフローチャートである。第３の実施形態の構成を示すブロック図である。図４５に示す初期点データ生成部２０の構成を示すブロック図である。図４６に示す内側ボクセルデータ生成部２０１の動作を示すフローチャートである。図４６に示す壁上点・球上点データ生成部２０３の動作を示すフローチャートである。図４６に示す分割点データ生成部２０６の動作を示すフローチャートである。図４６に示す初期点データ出力部２０８の動作を示すフローチャートである。図４６に示す構成の変形例を示すブロック図である。図５１に示す内側ボクセルデータ修正部２１０の動作を示すフローチャートである。内側ボクセルデータを示す説明図である。球上点の一例を示す説明図である。壁上点、球上点、ボクセルの中心点を示す説明図である。分割点を生成する動作を示す説明図である。境界層を形成するための点群の一例を示す説明図である。内側ボクセルデータを示す説明図である。分割領域データの一例を示す説明図である。

　始めに、従来の数値解析手法である有限要素法、有限体積法、ボクセル法、ボクセル法の改良手法、及び粒子法の問題点を解決し、解析精度の悪化を伴わずにソルバ処理における計算負荷の低減を図ることができる数値解析装置について、図面を参照して説明する。ここでは、本発明による計算用データ生成装置が生成した計算用データを入力して数値解析装置について説明する。また、以下の説明においては、車両の室内空間における空気の流速を数値解析によって求める場合について説明する。

　図１６は、数値解析装置Ａのハードウェア構成を概略的に示すブロック図である。この図に示すように、数値解析装置Ａは、パーソナルコンピュータやワークステーション等のコンピュータによって構成されるものであり、ＣＰＵ１、記憶装置２、ＤＶＤ（Digital Versatile Disc）ドライブ３、入力装置４、出力装置５、及び通信装置６を備えている。

　ＣＰＵ１は、記憶装置２、ＤＶＤドライブ３、入力装置４、出力装置５、及び通信装置６と電気的に接続されており、これらの各種装置から入力される信号を処理すると共に、処理結果を出力するものである。

　記憶装置２は、メモリ等の内部記憶装置及びハードディスクドライブ等の外部記憶装置によって構成されており、ＣＰＵ１から入力される情報を記憶すると共にＣＰＵ１から入力される指令に基づいて記憶した情報を出力するものである。
　そして、記憶装置２は、プログラム記憶部２ａとデータ記憶部２ｂとを備えている。

　プログラム記憶部２ａは、数値解析プログラムＰを記憶している。この数値解析プログラムＰは、所定のＯＳにおいて実行されるアプリケーションプログラムであり、コンピュータから構成される数値解析装置Ａを、数値解析を行うように機能させるものである。そして、数値解析プログラムＰは、数値解析装置Ａを、例えば計算用データモデル作成手段と物理量計算手段として機能させるものである。そして、図１６に示すように、数値解析プログラムＰは、プリ処理プログラムＰ１と、ソルバ処理プログラムＰ２と、ポスト処理プログラムＰ３とを有している。

　プリ処理プログラムＰ１は、ソルバ処理を実行するための前処理（プリ処理）を数値解析装置Ａに実行させるものであり、数値解析装置Ａを計算用データモデル作成手段として機能させることによって計算用データモデルを作成させる。また、プリ処理プログラムＰ１は、数値解析装置Ａに、ソルバ処理を実行するにあたり必要となる条件の設定を実行させ、さらには上記計算用データモデルや設定された条件を纏めたソルバ入力データファイルＦの作成を実行させる。

　そして、プリ処理プログラムＰ１は、数値解析装置Ａを計算用データモデル作成手段として機能させる場合に、まず数値解析装置Ａに対して、解析空間を含む３次元形状データを取得させ、この取得させた３次元形状データに含まれる解析空間を示す解析領域の作成を実行させる。

　その後、プリ処理プログラムＰ１は、数値解析装置Ａに対して、解析領域全体を有限個の領域に分割した分割領域を形成する。なお、後に詳説するが、ソルバ処理において、数値解析手法にて説明した離散化支配方程式（VertexとConnectivityとを必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化支配方程式）を用いる。このため、計算用データモデルの作成にあたり、保存則を満たす条件の下、分割領域の形状を任意に選択することができる。

　また、プリ処理プログラムＰ１は、数値解析装置Ａを計算用データモデル作成手段として機能させる場合に、数値解析装置Ａに対して、作成させた室内空間を示す解析領域に含まれる分割領域の各々の内部に対して１つのコントロールポイントを仮想的に配置する処理を実行させ、コントロールポイントの配置情報、及び各コントロールポイントが占めるコントロールボリュームの体積データを記憶させる。

　また、プリ処理プログラムＰ１は、数値解析装置Ａを計算用データモデル作成手段として機能させる場合に、数値解析装置Ａに対して、上記分割領域同士の境界面である境界面の面積及び法線ベクトルの算出を実行させ、これらの境界面の面積及び法線ベクトルを記憶させる。
　また、プリ処理プログラムＰ１は、数値解析装置Ａを計算用データモデル作成手段として機能させる場合に、コントロールボリューム（コントロールポイント）の結合情報（Link）を作成させ、このLinkを記憶させる。

　そして、プリ処理プログラムＰ１は、数値解析装置Ａに対して、上記各コントロールポイントが占めるコントロールボリュームの体積と、境界面の面積及び法線ベクトルと、コントロールポイント（すなわち分割領域）の配置情報（座標）と、Linkとを纏めさせて計算用データモデルを作成させる。

　また、プリ処理プログラムＰ１は、数値解析装置Ａに、前述のソルバ処理を実行するにあたり必要となる条件の設定を行わせる場合には、物性値の設定、境界条件の設定、初期条件の設定、計算条件の設定を行わせる。ここで、物性値とは、室内空間における空気の密度、粘性係数等である。境界条件とは、コントロールポイント間の物理量の交換の法則を規定するものであり、前述した式（１０）で示されるナビエ・ストークスの式に基づく離散化支配方程式、及び式（１１）で示される連続の式に基づく離散化支配方程式である。また、境界条件には、室内空間と外部空間との境界面に臨む分割領域を示す情報が含まれる。初期条件とは、ソルバ処理を実行する際の最初の物理量を示すものであり、各分割領域の流速の初期値である。計算条件とは、ソルバ処理における計算の条件であり、例えば反復回数や収束基準である。

　また、プリ処理プログラムＰ１は、数値解析装置Ａに、ＧＵＩ（Graphical User Interface）を形成させる。より詳細には、プリ処理プログラムＰ１は、出力装置５が備えるディスプレイ５ａに対してグラフィックを表示させると共に、入力装置４が備えるキーボード４ａやマウス４ｂによって操作が可能な状態とさせる。

　ソルバ処理プログラムＰ２（物理量計算プログラム）は、数値解析装置Ａにソルバ処理を実行させるものであり、数値解析装置Ａを物理量計算装置として機能させる。
　そして、ソルバ処理プログラムＰ２は、数値解析装置Ａを物理量計算手段として機能させる場合に、計算用データモデルが有するコントロールボリュームの体積と境界面の面積及び法線ベクトルとを含むソルバ入力データファイルＦを用いて、解析領域における物理量を物理量計算させる。

　そして、ソルバ処理プログラムＰ２は、数値解析装置Ａを物理量計算手段として機能させる場合に、数値解析装置Ａに対して、ソルバ入力データファイルＦに含まれるナビエ・ストークスの式及び連続の式の離散化係数行列の作成を実行させると共に、マトリックス形成用のデータテーブルの作成を実行させる。

　また、ソルバ処理プログラムＰ２は、数値解析装置Ａを物理量計算手段として機能させる場合に、数値解析装置Ａに対して、前述した式（１０）で示されるナビエ・ストークスの式に基づく離散化支配方程式、及び前述した式（１１）で示される連続の式に基づく離散化支配方程式から、下式（５８）で示すマトリックス計算用の大規模粗行列方程式の組み上げを実行させる。

　なお、式（５８）において［Ａ］が大規模粗行列を示し、［Ｂ］が境界条件ベクトルを示し、［Ｘ］が流速の解を示す。

　また、ソルバ処理プログラムＰ２は、上記離散化支配方程式に非圧縮性等の付帯条件が存在する場合には、数値解析装置Ａに対して、この付帯条件の行列方程式への組み上げを実行させる。そして、ソルバ処理プログラムＰ２は、数値解析装置Ａに対して、ＣＧ法（共役勾配法）等による行列方程式の解の計算、当該解の下式（５９）を用いた解のアップデート、収束条件の判定を実行させ、最終的な計算結果を取得させる。

　ポスト処理プログラムＰ３は、数値解析装置Ａにポスト処理を実行させるものであり、数値解析装置Ａに対して、ソルバ処理において取得された計算結果に基づく処理を実行させる。より詳細には、ポスト処理プログラムＰ３は、数値解析装置Ａに対して、計算結果の可視化処理、抽出処理を実行させる。ここで、可視化処理とは、例えば、断面コンタ表示、ベクトル表示、等値面表示、アニメーション表示を出力装置５に出力させる処理である。また、抽出処理とは、作業者が指定する領域の定量値を抽出して数値やグラフとして出力装置５に出力させる、あるいは作業者が指定する領域の定量値を抽出してファイル化したものの出力を実行させる処理である。また、ポスト処理プログラムＰ３は、数値解析装置Ａに対して、自動レポート作成、計算残差の表示及び分析を実行させる。

　データ記憶部２ｂは、図１６に示すように、計算用データモデルＭ、境界条件を示す境界条件データＤ１、計算条件を示す計算条件データＤ２、物性値を示す物性値データＤ３、及び初期条件を示す初期条件データＤ４を有するソルバ入力データファイルＦと、３次元形状データＤ５と、計算結果データＤ６等を記憶するものである。また、データ記憶部２ｂは、ＣＰＵ１の処理過程において生成される中間データを一時的に記憶する。

　ＤＶＤドライブ３は、ＤＶＤメディアＸを取り込み可能に構成されており、ＣＰＵ１から入力される指令に基づいて、ＤＶＤメディアＸに記憶されるデータを出力するものである。そして、ＤＶＤメディアＸに数値解析プログラムＰが記憶されており、ＤＶＤドライブ３は、ＣＰＵ１から入力される指令に基づいて、ＤＶＤメディアＸに記憶される数値解析プログラムＰを出力する。

　入力装置４は、数値解析装置Ａと作業者とのマンマシンインターフェイスであり、ポインティングデバイスであるキーボード４ａやマウス４ｂを備えている。出力装置５は、ＣＰＵ１から入力される信号を可視化して出力するものであり、ディスプレイ５ａ及びプリンタ５ｂを備えている。通信装置６は、数値解析装置ＡとＣＡＤ装置Ｃ等の外部装置との間においてデータの受け渡しを行うものであり、社内ＬＡＮ（Local Area Network）等のネットワークＢに対して電気的に接続されている。

　次に、このように構成された数値解析装置Ａを用いた数値解析方法（数値解析方法）について、図１７～図１９のフローチャートを参照して説明する。

　図１７のフローチャートに示すように、数値解析方法は、プリ処理（ステップＳ１）と、ソルバ処理（ステップＳ２）と、ポスト処理（ステップＳ３）とから構成されている。

　なお、数値解析方法を行うより前に、ＣＰＵ１は、ＤＶＤドライブ３に取り込まれたＤＶＤメディアＸに記憶された数値解析プログラムＰをＤＶＤメディアＸから取り出し、記憶装置２のプログラム記憶部２ａに記憶させる。そして、ＣＰＵ１は、入力装置４から数値解析の開始を指示する信号が入力されると、記憶装置２に記憶された数値解析プログラムＰに基づいて数値解析を実行する。より詳細には、ＣＰＵ１は、プログラム記憶部２ａに記憶されたプリ処理プログラムＰ１に基づいてプリ処理（ステップＳ１）を実行し、プログラム記憶部２ａに記憶されたソルバ処理プログラムＰ２に基づいてソルバ処理（ステップＳ２）を実行し、プログラム記憶部２ａに記憶されたポスト処理プログラムＰ３に基づいてポスト処理（ステップＳ３）を実行する。なお、このようにＣＰＵ１がプリ処理プログラムＰ１に基づくプリ処理（ステップＳ１）を実行することによって、数値解析装置Ａが計算用データモデル作成手段として機能される。また、ＣＰＵ１がソルバ処理プログラムＰ２に基づくソルバ処理（ステップＳ２）を実行することによって、数値解析装置Ａが物理量計算手段として機能される。

　図１８は、プリ処理（ステップＳ１）を示すフローチャートである。この図に示すように、プリ処理（ステップＳ１）が開始されると、ＣＰＵ１は、通信装置６に、ネットワークＢを介してＣＡＤ装置Ｃから車両の室内空間を含む３次元形状データＤ５を取得させる（ステップＳ１ａ）。ＣＰＵ１は、取得した３次元形状データＤ５を記憶装置２のデータ記憶部２ｂに記憶させる。

　続いて、ＣＰＵ１は、計算用データモデルの作成を実行する（ステップＳ１ｄ）。より詳細には、まずＣＰＵ１は、３次元形状データＤ５から、解析空間の全領域を含む分割領域にて分割された解析領域の作成を実行する。なお、解析領域を構成する分割領域は、前述の原理説明で詳説したように、VertexとConnectivityとを持たない計算用データモデルを作成するため、分割領域の幾何学的形状に制約を課すことなく計算用データモデル作成することができる。つまり、計算用データモデルを作成するにあたり、解析領域を構成する分割領域は、任意の形状を取ることができる。

　次にＣＰＵ１は、室内空間を示す解析領域に含まれる各分割領域内に１つのコントロールポイントを仮想的に配置する。ここでは、ＣＰＵ１は、分割領域の重心を算出し、各々の重心に対して１つのコントロールポイントを仮想的に配置する。そして、ＣＰＵ１は、コントロールポイントの配置情報、各コントロールポイントが占めるコントロールボリュームの体積（コントロールポイントが配置される分割領域の体積）を算出し、記憶装置２のデータ記憶部２ｂに一時的に記憶させる。また、ＣＰＵ１は、分割領域同士の境界面である境界面の面積及び法線ベクトルを算出し、これらの境界面の面積及び法線ベクトルを記憶装置２のデータ記憶部２ｂに一時的に記憶させる。また、ＣＰＵ１は、Linkを作成し、このLinkを記憶装置２のデータ記憶部２ｂに一時的に記憶させる。

　そして、ＣＰＵ１は、データ記憶部２ｂに記憶された、コントロールポイントの配置情報と、各コントロールポイントが占めるコントロールボリュームの体積と、境界面の面積及び法線ベクトルと、Linkとをデータベース化することによって計算用データモデルＭを作成し、作成した計算用データモデルを記憶装置２のデータ記憶部２ｂ内に記憶させる。

　なお、このようなステップＳ１ｄにおいては、室内空間を含む解析領域を同一形状の分割領域にて均等に分割し、さらに室内空間から食み出した分割領域を削除し、さらにその結果生じた解析領域と室内空間との隙間に新たな分割領域を充填することによって最終的な解析領域Ｋ２が作成される。このため、室内空間の全領域が重ならない分割領域によって充填された状態とされる。したがって、計算用データモデルは、前述した保存則を満足するための３つの条件（ａ）～（ｃ）を満たすものとされている。

　また、ステップＳ１ｄにおいて、先に分割領域を形成し、その後コントロールポイントを配置し、各コントロールポイントに対して、自らが配置された分割領域の体積を割り当てる構成を採用している。

　しかしながら、先にコントロールポイントを解析領域に配置し、各コントロールポイントに対して後から体積を割り当てることも可能である。具体的には、例えば、異なるコントロールポイントにぶつかるまでの半径や、結合関係にある（Linkで関連付けられた）コントロールポイントまでの距離に基づいて、各コントロールポイントに対して重み付けを行う。ここでコントロールポイントｉの重みをｗｉ、基準体積をＶ＋とし、コントロールポイントｉに割り当てられる体積Ｖｉを下式（６０）とする。

　さらに、各コントロールポイントの体積Ｖｉの総和は、解析領域の体積Ｖｔｏｔａｌと等しいため、下式（６１）が成り立つ。

　この結果、基準体積Ｖ＋は下式（６２）で求めることができる。

　したがって、各コントロールポイントに割り当てる体積は、式（６１），（６２）から求めることができる。このような方法を用いれば、プリ処理において、VertexとConnectivityとを用いることなく、計算用データモデルに持たせる分割領域の体積を求めることができる。

　また、当該計算用データモデルの作成（ステップＳ１ｄ）において、ＣＰＵ１は、ＧＵＩを形成し、ＧＵＩから指令（例えば分割領域の密度を示す指令や分割領域の形状を示す指令）が入力された場合には、当該指令を反映させた処理を実行する。したがって、作業者は、ＧＵＩを操作することによって、コントロールポイントの配置や分割領域の形状を任意に調節することが可能とされている。ただし、ＣＰＵ１は、数値解析プログラムに記憶された保存則を満足するための３つの条件に照らし合わせ、ＧＵＩから入力される指令が、当該条件から外れる場合には、その旨をディスプレイ５ａに表示させる。

　続いて、ＣＰＵ１は、物性値データの設定（ステップＳ１ｅ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に物性値の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される物性値を示す信号を物性値データＤ３としてデータ記憶部２ｂに一時的に記憶させることで物性値の設定を行う。なお、ここで言う物性値とは、室内空間における流体（すなわち空気）の特性値であり、空気の密度、粘性係数等である。

　続いて、ＣＰＵ１は、境界条件データの設定（ステップＳ１ｆ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に境界条件の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される境界条件を示す信号を境界条件データＤ１としてデータ記憶部２ｂに一時的に記憶させることで境界条件データの設定を行う。なお、ここで言う境界条件とは、室内空間の物理現象を支配する離散化支配方程式や、室内空間と外部空間との境界面に臨むコントロールポイントの特定情報、及び室内空間と外部空間と間における熱の伝熱条件等を示す。

　なお、数値解析方法においては、室内空間における流速を数値解析により求めることを目的とするものであるため、上記離散化支配方程式として、前述のナビエ・ストークスの式に基づく離散化支配方程式（１０）及び連続の式に基づく離散化支配方程式（１１）が用いられる。また、これらの離散化支配方程式は、例えば、数値解析プログラムＰに予め記憶された複数の離散化支配方程式をディスプレイ５ａ上に表示された複数の離散化支配方程式から作業者がキーボード４ａやマウス４ｂを用いることによって選択される。

　続いて、ＣＰＵ１は、初期条件データの設定（ステップＳ１ｇ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に初期条件の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される初期条件を示す信号を初期条件データＤ４としてデータ記憶部２ｂに一時的に記憶させることで初期条件データの設定を行う。なお、ここ言う初期条件とは、各コントロールポイント（各分割領域）における初期流速である。

　続いて、ＣＰＵ１は、計算条件データの設定（ステップＳ１ｈ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に計算条件の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される計算条件を示す信号を計算条件データＤ２としてデータ記憶部２ｂに一時的に記憶させることで計算条件データの設定を行う。なお、ここで言う計算条件とは、ソルバ処理（ステップＳ２）における計算の条件であり、例えば、反復回数や収束基準を示す。

　続いて、ＣＰＵ１は、ソルバ入力データファイルＦの作成（ステップＳ１ｉ）を行う。具体的には、ＣＰＵ１は、ステップＳ１ｄにて作成された計算用データモデルＭと、ステップＳ１ｅで設定された物性値データＤ３と、ステップＳ１ｆで設定された境界条件データＤ１と、ステップＳ１ｇで設定された初期条件データＤ４と、ステップＳ１ｈで設定された計算条件データＤ２とをソルバ入力データファイルＦに格納することによってソルバ入力データファイルＦを作成する。なお、このソルバ入力データファイルＦは、データ記憶部２ｂに記憶される。

　以上のようなプリ処理（ステップＳ１）が完了すると、ＣＰＵ１は、ソルバ処理プログラムＰ２に基づいて、図１９のフローチャートに示すソルバ処理（ステップＳ２）を実行する。図１９に示すように、ソルバ処理（ステップＳ２）が開始されると、ＣＰＵ１は、プリ処理（ステップＳ１）で作成されたソルバ入力データファイルＦを取得する（ステップＳ２ａ）。なお、前述した数値解析方法のように、単一の装置（数値解析装置Ａ）によってプリ処理及びソルバ処理を実行する場合には、既にデータ記憶部２ｂにソルバ入力データファイルＦが記憶されているため、ステップＳ２ａを省略することができる。ただし、プリ処理（ステップＳ１）とソルバ処理（ステップＳ２）とが異なる装置において実行される場合には、ネットワークやリムーバルディスクによって搬送されるソルバ入力データファイルＦを取得する必要があるため、本ステップＳ２ａを行う必要がある。

　続いて、ＣＰＵ１は、ソルバ入力データの整合性を判定する（ステップＳ２ｂ）。なお、ソルバ入力データとは、ソルバ入力データファイルＦに格納されたデータを示し、計算用データモデルＭ、境界条件データＤ１、計算条件データＤ２、物性値データＤ３及び初期条件データＤ４である。具体的には、ＣＰＵ１は、ソルバ処理において物理量計算を実行可能なソルバ入力データがソルバ入力データファイルＦに全て格納されているかを分析することによってソルバ入力データの整合性の判定を行う。
　そして、ＣＰＵ１は、ソルバ入力データが不整合であると判定した場合には、ディスプレイ５ａにエラーを表示させ（ステップＳ２ｂ）、さらには不整合である部分のデータを入力するための画面を表示させる。その後、ＣＰＵ１は、ＧＵＩから入力される信号に基づいてソルバ入力データの調整を行い（ステップＳ２ｃ）、再度ステップＳ２ａを実行する。

　一方、ＣＰＵ１は、ステップＳ２ｂにおいてソルバ入力データの整合性があると判定した場合には、初期計算処理（ステップＳ２ｅ）を実行する。具体的には、ＣＰＵ１は、境界条件データＤ１に記憶された離散化支配方程式（すなわちナビエ・ストークスの式に基づく離散化支配方程式（１０）及び連続の式に基づく離散化支配方程式（１１））から離散化係数行列を作成し、さらにマトリクス計算用のデータテーブルの作成を行うことによって初期計算処理を行う。

　続いて、ＣＰＵ１は、大規模粗行列方程式の組み上げ（ステップＳ２ｆ）を行う。具体的には、ＣＰＵ１は、ナビエ・ストークスの式に基づく離散化支配方程式（１０）及び連続の式に基づく離散化支配方程式（１１）から、前述の式（５８）で示すマトリックス計算用の大規模粗行列方程式の組み上げを行う。

　続いて、ＣＰＵ１は、離散化支配方程式に、非圧縮性や接触等の付帯条件が存在するかの判定を行う。この付帯条件は、例えば、境界条件データとしてソルバ入力データファイルＦに格納されている。そして、ＣＰＵ１は、離散化支配方程式に付帯条件が存在すると判定した場合には当該付帯条件の大規模行列方程式の組み込み（ステップＳ２ｈ）を実行した後に大規模行列方程式の計算（ステップＳ２ｉ）を実行する。一方、ＣＰＵ１は、離散化支配方程式に付帯条件が存在しないと判定した場合には付帯条件の大規模行列方程式の組み込み（ステップＳ２ｈ）を実行することなく大規模行列方程式の計算（ステップＳ２ｉ）を実行する。そしてＣＰＵ１は、大規模行列方程式を例えば、ＣＧ法（共役勾配法）によって解き、前述の式（５９）を用いて解のアップデート（ステップＳ２ｊ）を行う。

　続いて、ＣＰＵ１は、式（５９）の残差が収束条件に達したか否かの判定（ステップＳ２ｇ）を行う。具体的には、ＣＰＵ１は、式（５９）の残差を計算し、計算条件データＤ２に含まれる収束条件と比較し、これによって式（５９）の残差が収束条件に達したか否かの判定を行う。そして、残差が収束条件に達していないと判定した場合には、ＣＰＵ１は、物性値のアップデートを行った後、再度ステップＳ２ｇを実行する。つまり、ＣＰＵ１は、式（５９）の残差が収束条件に達するまで、物性値のアップデートを行いながらステップＳ２ｆ～Ｓ２ｇを繰り返し行う。

　一方、残差が収束条件に達したと判定した場合には、ＣＰＵ１は、計算結果の取得を行う（ステップＳ２ｌ）。具体的には、ＣＰＵ１は、直前のステップＳ２ｉにおいて計算された物理量の解を計算結果データとしてデータ記憶部２ｂに記憶させることによって計算結果の取得を行う。このようなソルバ処理（ステップＳ２）によって、室内空間における空気の流速が求められる。

　以上のようなソルバ処理（ステップＳ２）が完了すると、ＣＰＵ１は、ポスト処理プログラムＰ３に基づいてポスト処理（ステップＳ３）を実行する。具体的には、例えばＣＰＵ１は、ＧＵＩから入力される指令に基づいて、計算結果データから、例えば断面コンタデータ、ベクトルデータ、等値面データ、アニメーションデータを生成し、当該データを、出力装置５に可視化させる。

　また、ＣＰＵ１は、ＧＵＩから入力される指令に基づいて、室内空間の一部における定量値（計算結果）を抽出して数値やグラフとし、この数値やグラフを出力装置５に可視化させ、さらには数値やグラフをファイルとして纏めて出力する。また、ＣＰＵ１は、ＧＵＩから入力される指令に基づいて、例えば計算結果データから自動レポート作成、計算残差の表示及び分析を行ってその結果を出力する。

　以上のような数値解析装置Ａ、数値解析方法及び数値解析プログラムによれば、プリ処理にてコントロールボリュームの体積と境界面の面積及び法線ベクトルとを有する計算用データモデルＭが作成され、ソルバ処理にて計算用データモデルＭに含まれるコントロールボリュームの体積と境界面の面積及び法線ベクトルとを用いて各コントロールボリュームにおける物理量が計算される。

　つまり、前述の数値解析手法にて説明したように、数値解析装置Ａ、数値解析方法及び数値解析プログラムによれば、VertexとConnectivityを有する計算用データモデルを作成することなく数値解析を行うことが可能となる。このため、３次元形状データの修正あるいは変更作業に対する規制も大幅に緩和され、VertexとConnectivityを有する計算用データモデルと比較して計算用データモデルＭを遥かに容易に作成することが可能となる。したがって、計算用データモデルＭの作成における作業負担を軽減することが可能となる。

　また、数値解析装置Ａ、数値解析方法及び数値解析プログラムにおいては、従来の数値解析手法と異なり、ソルバ処理においてVertexとConnectivityを用いてコントロールボリュームの体積及び境界面の面積及び法線ベクトルを算出する必要がない。したがって、ソルバ処理における計算負荷を低減させることが可能となる。したがって、数値解析装置Ａ、数値解析方法及び数値解析プログラムによれば、計算用データモデルの作成における作業負担を軽減すると共に、ソルバ処理における計算負荷の低減を図ることが可能となる。

　また、数値解析装置Ａ、数値解析方法及び数値解析プログラムにおいては、解析領域に分割領域を重なることなく充填させている。このため、前述した保存側を満足するための３つの条件（ａ）～（ｃ）が満たされることとなり、保存則を満足して流速を計算することができる。

　なお、計算用データモデルＭは、従来の有限体積法や有限要素法に用いられるメッシュデータや、従来の粒子法に用いられる粒子データ、または単なる点群データから容易に計算用データモデルを作成することができる。この場合でも、ボクセル法のような外部空間との境界領域における分割領域の修正は必要ない。例えば、メッシュデータから計算用データモデルを作成する場合には、各要素を分割領域（コントロールポイントが占めるコントロールボリューム）として捉えることによって容易に行うことができる。また、粒子データから計算用データモデルを作成する場合には、解析領域に配置されたコントロールポイントを囲う閉空間を解析領域に充填して得られる閉空間を分割領域（コントロールポイントが占めるコントロールボリューム）として捉えることによって容易に行うことができる。

　また、数値解析装置Ａ、数値解析方法及び数値解析プログラムによれば、前述のように、プリ処理における計算用データモデルの作業負担が大幅に減少し、またソルバ処理における計算負荷を低減させることが可能となる。したがって、解析領域の形状が時系列的に変化する場合、すなわち解析領域が移動境界を含む場合であっても、図２０のフローチャートに示すように、解析領域が形状変化するたびにプリ処理とソルバ処理とを繰り返し行うことによって、現実的な時間内で物理量の計算が可能となる。

　また、物理量計算において、境界面上における流速を算出する場合に、当該境界面を挟むコントロールポイントにおける流速を当該コントロールポイントから境界面までの距離の比率に応じて重み付け平均した値を境界面上における流速としても良い。

　このような場合には、境界面における流速を、単純に当該境界面を挟むコントロールポイントの流速の平均とする場合と比較して、より正確に室内空間の流速を求めることが可能となる。ただし、このような場合には、計算用データモデルが、境界面がコントロールポイントとコントロールポイントとを結ぶ線分のどの内分点に存在するかを示す比率αをデータとして有している必要がある。このため、隣り合うコントロールボリュームの内部に配置された各コントロールポイントから当該コントロールボリュームに挟まれる境界面までの距離の比率をさらに含む計算用データモデルＭをプリ処理にて作成し、この比率に基づいて境界面における流速をソルバ処理にて計算する。

　なお、物理量計算において、コントロールポイント間を結ぶベクトルと当該コントロールポイントに挟まれた境界面が直交する場合には、法線ベクトルをコントロールポイント間を結ぶ距離ベクトルの単位ベクトルに置き換えて物理量計算することができる。ただし、このような場合には、計算用データモデルに上記距離ベクトルあるいは当該距離ベクトルを算出するためのコントロールポイントの位置座標を示すデータを備えさせる必要がある。このため、隣り合うコントロールボリュームの内部に配置された各コントロールポイントを結ぶ距離ベクトルあるいは該距離ベクトルを算出するためのコントロールポイントの座標をさらに含む計算用データモデルＭをプリ処理にて作成し、距離ベクトルが当該コントロールポイントに挟まれた境界面と直交する場合に、距離ベクトルを境界面の法線ベクトルに置き換えて流速を計算する。

　前述した説明においては、運動量保存の方程式の変形例であるナビエ・ストークスの式及び連続の式から導出した離散化支配方程式を用いて空気の流速を数値解析によって求める構成について説明した。しかしながら、これに限定されるものではなく、質量保存の方程式、運動量保存の方程式、角運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式及び波動方程式の少なくともいずれかから導出した離散化支配方程式を用いて物理量を数値解析によって求めることが可能である。

　また、前述した境界面特性量として、境界面の面積と境界面の法線ベクトルとを用いる構成について説明した。しかしながら、これに限定されるものではなく、境界面特性量として他の量（例えば境界面の周長）を用いることもできる。

　また、保存則を満足するために前述の３つの条件を満たすように計算用データモデルを作成する構成について説明した。しかしながら、これに限定されるものではなく、保存則を満足させる必要がない場合には、計算用データモデルを必ずしも前述の３つの条件を満たすように作成する必要はない。

　また、分割領域の体積を、当該分割領域の内部に配置されるコントロールポイントが占めるコントロールボリュームの体積として捉えた構成について説明した。しかしながら、これに限定されるものではなく、分割領域の内部に対してコントロールポイントを配置する必要は必ずしもない。このような場合には、コントロールポイントが占めるコントロールボリュームの体積を分割領域の体積に置き換えることによって数値解析を行うことができる。

　また、数値解析プログラムＰがＤＶＤメディアＸに記憶されて搬送可能な構成について説明した。しかしながら、これに限定されるものではなく、数値解析プログラムＰを他のリムーバブルメディアに記憶させて搬送可能とする構成を採用することもできる。また、プリ処理プログラムＰ１とソルバ処理プログラムＰ２とを別々のリムーバブルメディアに記憶させて搬送可能とすることもできる。また、数値解析プログラムＰは、ネットワークを介して伝達することも可能である。

＜第１の実施形態＞
　次に、本発明の第１の実施形態による計算用データ生成装置を説明する。図２１は同実施形態の構成を示すブロック図である。この図において、符号Ｄは、計算用データを生成する計算用データ生成装置であり、コンピュータ装置から構成する。計算用データ生成装置Ｄは、コンピュータネットワークＢと接続されており、このコンピュータネットワークに接続されているＣＡＤ装置Ｃとの間で情報通信が可能である。符号４は、キーボードやマウスから構成する入力部である。符号５は、液晶ディスプレイ装置等から構成する出力装置である。

　符号１１は、入力装置４から入力されるパラメータを読み込むパラメータ入力部である。符号１２は、パラメータ入力部１１が入力した数値計算に必要なパラメータを記憶するパラメータ記憶部である。符号１３は、コンピュータネットワークＢを介して、ＣＡＤ装置Ｃとの間で情報通信を行う通信部である。符号１４は、通信部１３を介して、ＣＡＤ装置Ｃから解析対象物体の形状データを入力する形状データ入力部である。符号１５は、形状データ入力部１４が入力した形状データを記憶する形状データ記憶部である。

　符号１６は、形状データ記憶部１５に記憶されている形状データを入力して、解析対象物体のパッチデータ（物体を多角形で表現したデータ）を生成するパッチデータ生成部である。符号１７は、パッチデータ生成部１６が生成したパッチデータを記憶するパッチデータ記憶部である。符号１８は、パラメータ記憶部１２に記憶されているパラメータと、パッチデータ記憶部１７に記憶されているパッチデータとを入力して、解析領域を定義するボクセルデータを生成するボクセルデータ生成部である。

　符号１９は、ボクセルデータ生成部１８が生成したボクセルデータを記憶するボクセルデータ記憶部である。符号２０は、ボクセルデータ記憶部１９に記憶されているボクセルデータに基づいて、初期点データを生成する初期点データ生成部である。符号２１は、初期点データ生成部２０において生成された初期点データを記憶する初期点データ記憶部である。符号２２は、初期点データ記憶部２１に記憶されている初期点を間引いて、初期点データ記憶部２１に記憶されている初期点データ更新する処理を行うデータ間引き部である。

　符号２３は、初期点データ記憶部２１に記憶されている初期点データと、ボクセルデータ記憶部１９に記憶されているボクセルデータとを入力して、解析対象物体の分割領域データを生成する分割領域データ生成部である。符号２４は、分割領域データ生成部２３が生成した分割領域データを記憶する分割領域データ記憶部である。符号２５は、分割領域データ記憶部２４に記憶されている分割領域データを入力して、計算用データを生成する計算用データ生成部である。符号２６は、計算用データ生成部２５が計算用データを生成する際の中間データを記憶する中間データ記憶部である。符号２７は、計算用データ生成部２５が生成した計算用データを記憶する計算用データ記憶部である。この計算用データ記憶部２７に記憶されたデータは、数値解析の入力データとなる。

　次に、図２２を参照して、図２１に示すパラメータ記憶部１２のテーブル構造を説明する。図２２は、図２１に示すパラメータ記憶部１２のテーブル構造を示す図である。パラメータ記憶部１２は、パラメータ入力部１１が入力装置４から入力されたパラメータを記憶するものであり、解析領域データ、解析領域基準点データ、ボクセル分割数データからなる。解析領域データとは、直方体からなる解析領域を定義するためのデータであり、直方体の２点の座標値から構成する。解析領域基準点データは、解析対象物体内に存在する点の座標値であり、物体の内側と外側を判定するためのものである。ボクセル分割数データとは、解析対象物体を完全に含む解析領域の直方体（点(Xmin,Ymin,Zmin)及び点(Xmax, Ymax,Zmax)の２点の解析領域データで定義される直方体）の直交座標系の各軸方向の分割数である。これらのデータは、解析作業者が、入力装置４から入力したものである。

　次に、図２３を参照して、図２１に示す形状データ記憶部１５に記憶される形状データについて説明する。図２３は、図２１に示す形状データ記憶部１５に記憶される形状データの一例を示す図である。形状データ記憶部１５には、ＣＡＤ装置Ｃによって設計された解析対象物体の形状を定義するデータが記憶される。図２３に示す例では、説明を簡単にするために、解析対象物体が円柱状の形状であるものとして説明する。形状データは、例えば、形状の稜線を定義するための座標値等から構成する。ここでは、物体の内側か外側かのデータは含まれておらず、解析領域基準点（Ｂｘ，Ｂｙ，Ｂｚ）を円柱の内側に定義することにより、解析領域基準点が存在する側を物体の内側とする。

　次に、図２４を参照して、解析領域と解析対象物体の位置関係について説明する。図２４は、解析領域と解析対象物体の位置関係を示す図である。解析対象物体は、解析領域データで定義される直方体の中に完全に含まれる位置関係を有している。

　次に、図２５を参照して、図２１に示すパッチデータ記憶部１７に記憶されるパッチデータについて説明する。図２５は、図２１に示すパッチデータ記憶部１７に記憶されるパッチデータの一例を示す図である。パッチデータとは、解析対象物体の外形を多角形（ポリゴン）で表現したデータである。ここでは、この多角形をパッチと称する。パッチデータ生成部１６は、形状データを読み込んで、公知の方法を用いて、パッチデータを生成して、パッチデータ記憶部１７に記憶する。パッチデータの生成方法は、公知の生成方法のいずれでも使用可能であるため、ここでは、パッチデータ生成処理の詳細な説明を省略する。図２５に示す例では、三角形のパッチデータを生成した例を示している。各多角形（ここでは三角形）は、各頂点（Ｐ１、Ｐ２、Ｐ３）の座標値によって定義されるものである。図２５においては、説明を簡単にするために、円柱を六角柱で表現しているが、実際には、所望の計算精度が得られる多角柱で表現されるものである。

　次に、図２６を参照して、図２１に示すパッチデータ記憶部１７のテーブル構造を説明する。図２６は、図２１に示すパッチデータ記憶部１７のテーブル構造を示す図である。パッチデータ記憶部１７には、解析対象物体上に生成されたパッチの数と、各パッチ（ここでは三角形）それぞれについて、頂点の３次元座標値が、パッチを一意に特定可能なパッチＩＤとが関係付けられて記憶される。

　次に、図２７を参照して、ボクセルデータについて説明する。図２７は、ボクセルの構成を示す図である。ボクセルは、解析領域データで定義される解析領域（直方体の領域）を、ボクセル分割数データ（Ｎｘ、Ｎｙ、Ｎｚ）に基づき分割して得られる小さな直方体である。図２７に示す例では、Ｎｘ＝６、Ｎｙ＝９、Ｎｚ＝６である場合の分割例を示している。

　次に、図２８を参照して、図２１に示すボクセルデータ記憶部１９のテーブル構造を説明する。図２８は、図２１に示すボクセルデータ記憶部１９のテーブル構造を示す図である。ボクセルデータ記憶部１９は、Ｎｘ×Ｎｙ×Ｎｚ個のボクセル属性データを記憶することができる３次元テーブルを備えており、ボクセルの分割と一致したテーブル構造になっている。ボクセル属性データは、３つの属性のいずれかを示す値が記憶される。１つ目の属性は、パッチと交差しているボクセルであること示す属性であり、交差しているパッチのパッチＩＤが記憶される。２つ目の属性は、解析対象物体の内部に存在するボクセルであることを示す属性であり、－２が記憶される。３つ目の属性は、解析対象物体の外部に存在するボクセルであることを示す属性であり、－１が記憶される。

　次に、図２９を参照して、各ボクセルを定義するデータについて説明する。図２９は、各ボクセルを定義するデータを求める方法を示す図である。各ボクセルは、ボクセルの中心点座標値と各軸方向の寸法によって定義される。各ボクセルの寸法（Ｓｘ、Ｓｙ、Ｓｚ：Ｘ方向の辺の長さ、Ｙ方向の辺の長さ、Ｚ方向の辺の長さ）は全て同一であり、分割領域の各軸方向の長さをボクセル分割数で除算して求めることができる。また、中心座標値（Ｘｃ，Ｙｃ，Ｚｃ）は、ボクセルの位置（ｉ，ｊ，ｋ：Ｘ方向の位置（何番目），Ｙ方向の位置，Ｚ方向の位置）と、分割領域のサイズと、ボクセル分割数とから簡単な計算式によって求めることができる。１つのボクセルを特定するためのデータ（中心座標値と寸法）は、簡単な計算式で求めることができため、必要な時に演算によって求めるようにして、ボクセルデータ記憶部１９には、属性データのみが記憶され、中心座標値と寸法については記憶されていない。

　次に、図３０を参照して、図２１に示す初期点データ記憶部２１のテーブル構造を説明する。図３０は、図２１に示す初期点データ記憶部２１のテーブル構造を示す図である。初期点データ記憶部２１には、初期点の個数と、初期点のそれぞれについて、各初期点を一意に識別可能な初期点ＩＤと、初期点の３次元座標値とが関係付けられて記憶される。この初期点は、各ボクセルの中心点から生成し、データ間引き部２２が初期点を間引き、不要な初期点を削除した後の初期点データが記憶されることになる。

　次に、図３１を参照して、図２１に示す分割領域データ記憶部２４のテーブル構造を説明する。図３１は、図２１に示す分割領域データ記憶部２４のテーブル構造を示す図である。分割領域データ記憶部２４は、Ｎｘ×Ｎｙ×Ｎｚ個（ボクセルと同数）の分割領域データを記憶することができる３次元テーブルを備えており、ボクセルの分割と一致したテーブル構造になっている。分割領域データは、２つの状態のいずれかを示す値が記憶される。一方の状態は、分割領域に含まれる状態を示すものであり、そのボクセルの中心座標から最も近い初期点の初期点ＩＤが記憶される。ここで記憶される初期点ＩＤは、分割領域ＩＤに相当する。他方の状態は、解析対象物体の外部に存在するボクセルである状態を示すものであり、－１が記憶される。分割領域は、解析対象物体内を複数の３次元領域に分割したものであり、解析対象物体の外部に存在するボクセルは、ボクセルデータ記憶１９に記憶されるボクセル属性がそのまま記憶されることになる。

　次に、図３２を参照して、図２１に示す中間データ記憶部２６のテーブル構造を説明する。図３２は、図２１に示す中間データ記憶部２６のテーブル構造を示す図である。中間データ記憶部２６には、分割領域ＩＤが異なる隣り合うボクセルそれぞれが属する分割領域ＩＤの組合せと、法線ベクトルと、法線ベクトルの長さからもとめた境界面の面積を記憶する領域を備えている。この記憶領域は、計算用データを生成するために、分割領域データ記憶部２４に記憶されている分割領域データを集計する際に得られる中間データを記憶するものである。

　次に、図３３を参照して、図２１に示す計算用データ記憶部２７のテーブル構造を説明する。図３３は、図２１に示す計算用データ記憶部２７のテーブル構造を示す図である。計算用データ記憶部２７には、分割領域の個数（初期点の個数と同数）と、分割領域を一意に識別可能な分割領域ＩＤ毎に、分割領域の重心座標と分割領域の体積が記憶される。また、分割領域の境界面の個数と、境界面を一意に識別することが可能な境界面ＩＤ毎に分割領域ＩＤの組合せと、法線ベクトルと、境界面の面積とが記憶される。

　次に、図２１に示す計算用データ生成装置Ｄの処理動作を説明する。まず、計算用データを作成しようとする作業者は、入力装置４から、形状データ入力部１４に対して形状データを入力する指示の操作を行う。形状データ入力部１４は、通信部１３を介して、ＣＡＤ装置Ｃから解析対象物体の形状データを入力して、入力した形状データを形状データ記憶部１５に記憶する。形状データが形状データ記憶部１５に記憶されると、パッチデータ生成部１６は、形状データ記憶部１５から形状データを読み込み、パッチデータを生成して、パッチデータ記憶部１７に記憶する。これにより、パッチデータ記憶部１７には、図２６に示すパッチデータが記憶されることになる。

　また、作業者は、入力装置４からパラメータの入力を行う。これを受けて、パラメータ入力部１１は、入力装置から入力されたパラメータを読み込み、パラメータ記憶部１２に記憶する。これにより、パラメータ記憶部１２には、図２２に示すパラメータが記憶されることになる。

　次に、作業者は、入力装置４から計算用データの生成を指示する操作を行う。これを受けて、ボクセルデータ生成部１８は、パラメータ記憶部１２に記憶されているパラメータと、パッチデータ記憶部１７に記憶されているパッチデータとを入力して、ボクセルデータの生成を行う。

　ここで、図３４を参照して、図２１に示すボクセルデータ生成部１８がボクセルデータを生成する処理動作を説明する。図３４は、図２１に示すボクセルデータ生成部１８がボクセルデータを生成する処理動作を示すフローチャートである。まず、ボクセルデータ生成部１８は、パラメータ記憶部１２に記憶されている解析領域データとボクセル分割数データを読み込み、これらのデータからボクセルを定義して、各ボクセルの中心座標と寸法を求める（ステップＳ１１）。そして、ボクセルデータ生成部１８は、定義したボクセルを１つ選択する（ステップＳ１２）。

　次に、ボクセルデータ生成部１８は、パッチデータ記憶部１７に記憶されているパッチデータを参照して、選択したボクセルとパッチが交差するかを幾何学的に計算して求める（ステップＳ１３）。そして、ボクセルデータ生成部１８は、選択したボクセルとパッチが交差するか否かを判定する（ステップＳ１４）。この判定の結果、交差しなければボクセルデータ生成部１８は、選択したボクセルと解析対象物体との位置関係を求め（ステップＳ１５）、選択したボクセルが解析対象物体の外部領域であるか否かを判定する（ステップＳ１６）。この判定の結果、ボクセルが解析対象物体の外部に位置するボクセルであれば、ボクセルデータ記憶部１９中の選択したボクセルの記憶領域に、－１（解析対象物体外部のボクセルを示す値）を記憶し（ステップＳ１７）、内部に位置するボクセルであれば、－２（解析対象物体内部のボクセルを示す値）を記憶する（ステップＳ１８）。

　一方、選択したボクセルとパッチが交差する場合、ボクセルデータ生成部１８は、ボクセルデータ記憶部１９中の選択したボクセルの記憶領域に、交差するパッチのパッチＩＤを記憶する（ステップＳ１９）。そして、ボクセルデータ生成部１８は、定義したボクセルの全てについて処理を行ったか否かを判定し（ステップＳ２０）、全てのボクセルについて属性データを記憶するまでステップＳ１２～Ｓ１９の処理動作を繰り返す。この処理動作によって、ボクセルデータ記憶部１９には、全てのボクセルについて属性データが関係付けられて記憶されたことになる。

　次に、ボクセルデータ記憶部１９に属性データが記憶されると、初期点データ生成部２０は、ボクセルデータ記憶部１９に記憶されているボクセルそれぞれについて、中心座標を求め、この中心座標値を初期点とし、各初期点に初期点ＩＤを付与して、初期点データ記憶部２１に記憶する。これにより、初期点データ記憶部２１（図３０）には、ボクセルと同数の初期点が記憶されたことになる。

　次に、データ間引き部２２は、初期点データ記憶部２１に記憶されている初期点の数を減らすために、初期点データの間引きを行う。ここで行う間引き処理は、公知の方法を用いて行う。例えば、所定の空間間隔で所定の数の初期点データを削除するなどの処理によって行えばよい。そして、データ間引き部２２は、削除するべき初期点データを初期点データ記憶部２１上から削除し、残った初期点データの数をカウントして、初期点個数として初期点データ記憶部２１に記憶する。

　次に、分割領域データ生成部２３は、初期点データ記憶部２１に記憶されている初期点データと、ボクセルデータ記憶部１９に記憶されているボクセルデータとを参照して、分割領域データを生成して、分割領域データ記憶部２４に記憶する。

　ここで、図３５を参照して、分割領域データ生成部２３が分割領域データを生成して、分割領域データ記憶部２４に記憶する処理動作を説明する。図３５は、分割領域データ生成部２３が分割領域データを生成して、分割領域データ記憶部２４に記憶する処理動作を示すフローチャートである。まず、分割領域データ生成部２３は、分割領域データ記憶部２４のテーブルに全て－１を記憶する（ステップＳ２１）。続いて、分割領域データ生成部２３は、ボクセルデータ記憶部１９に記憶されているボクセルの属性データが－２またはパッチＩＤであるボクセルの中からボクセルを１つ選択し（ステップＳ２２）、この選択したボクセルの中心座標を求める。そして、分割領域データ生成部２３は、初期点データ記憶部２１を参照して、求めたボクセルの中心座標に最も近い（直線距離が最も短い）初期点を選択し（ステップＳ２３）、選択した初期点の初期点ＩＤを分割領域データ記憶部２４中の選択したボクセルの記憶領域に記憶する（ステップＳ２４）。

　次に、分割領域データ生成部２３は、ボクセルの全てについて処理を行ったか否かを判定し（ステップＳ２５）、全てのボクセルについて分割領域データ（初期点ＩＤ）を記憶するまでステップＳ２２～Ｓ２４の処理動作を繰り返す。この処理動作によって、分割領域データ記憶部２４には、全てのボクセルについて、そのボクセルが属する分割領域のＩＤ（初期点ＩＤ）か、または－１のいずれかが記憶されたことになる。

　次に、計算用データ生成部２５は、分割領域データ記憶部２４に記憶されている分割領域データを読み込み、計算用データを生成して計算用データ記憶部２７に記憶する。

　ここで、図３６を参照して、計算用データ生成部２５が分割領域データから計算用データを生成して、計算用データ記憶部２７に記憶する処理動作を説明する。図３６は、計算用データ生成部２５が分割領域データから計算用データを生成して、計算用データ記憶部２７に記憶する処理動作を示すフローチャートである。まず、計算用データ生成部２５は、分割領域データ記憶部２４を参照して、分割領域データとして分割領域ＩＤが記憶されている（－１が記憶されていない）ボクセルの中から基準となるボクセルを１つ選択する（ステップＳ３１）。そして、計算用データ生成部２５は、分割領域データ記憶部２４を参照して、選択した基準ボクセルに接している６つの隣接ボクセル（上下左右前後のボクセル）を特定する（ステップＳ３２）。

　次に、計算用データ生成部２５は、特定した６つの隣接ボクセルから１つの隣接ボクセルを選択し（ステップＳ３３）、基準ボクセルと選択した隣接ボクセルのそれぞれが属する分割領域を特定する。この分割領域の特定は、分割領域データとして記憶されている分割領域ＩＤにより行う。そして、計算用データ生成部２５は、２つの分割領域ＩＤが同一か否かを判定する（ステップＳ３５）。この判定の結果、同じ分割領域ＩＤであれば、ステップＳ３３に戻り、次の隣接ボクセルを選択する。

　一方、２つの分割領域ＩＤが同じでなければ、計算用データ生成部２５は、中間データ記憶部２６を参照して、２つの分割領域ＩＤの組合せが既に登録されているかを探索する（ステップＳ３６）。そして、計算用データ生成部２５は、２つの分割領域ＩＤの組合せが既に登録されているか否かを判定し（ステップＳ３７）、登録されていなければ、この２つの分割領域ＩＤを中間データ記憶部２６に記憶する（ステップＳ３８）。判定の結果、２つの分割領域ＩＤの組合せが既に登録されている場合は、分割領域ＩＤの登録は行わない。

　次に、計算用データ生成部２５は、基準ボクセルと選択した隣接ボクセルとが接している面の法線ベクトルを中間データ記憶部２６に書き込む。このとき、既に法線ベクトルが書き込まれていた場合、計算用データ生成部２５は、書き込まれていた法線ベクトルに対して、新たに求めた法線ベクトルを加算して、新たに中間データ記憶部２６に書き込む（ステップＳ３９）。

　次に、計算用データ生成部２５は、６つの隣接ボクセル全ての処理を行ったか否かを判定し（ステップＳ４０）、６つの隣接ボクセルに対する処理が終わるまで、ステップＳ３３～Ｓ３９の処理を繰り返す。そして、６つの隣接ボクセルに対する処理が終了した時点で、全てのボクセルについて処理したか（全てのボクセルを基準ボクセルとして選択したか）を判定し（ステップＳ４１）、全てのボクセル（－１でないボクセル）に対する処理が終わるまで、ステップＳ３１～Ｓ４０の処理を繰り返す。この処理動作により、中間データ記憶部２６には、隣接する分割領域ＩＤの組合せについて法線ベクトルの加算値が得られたことになる。

　次に、計算用データ生成部２５は、２つの分割領域ＩＤの組合せそれぞれについて、法線ベクトルの加算値から法線ベクトルの長さを求め、この法線ベクトルの長さを境界面の面積として中間データ記憶部２６に記憶する（ステップＳ４２）。そして、計算用データ生成部２５は、中間データ記憶部２６に記憶されている２つの分割領域ＩＤの組合せから分割領域の個数と、境界面の個数とを求める（ステップＳ４３）。

　次に、計算用データ生成部２５は、分割領域データ記憶部２４に記憶されている分割領域データと、中間データ記憶部２６に記憶されている中間データとに基づき、計算用データを生成して、計算用データ記憶部２７に記憶する。このとき、計算用データ生成部２５は、分割領域データ記憶部２４に記憶されている分割領域データから分割領域それぞれの重心位置座標と分割領域の体積（ボクセルの個数にボクセルの体積を乗算したもの）を求めて計算用データ記憶部２７に記憶する。また、計算用データ生成部２５は、中間データ記憶部２６に記憶されている中間データを計算用データ記憶部２７に記憶するともに、先に求めた分割領域の個数と、境界面の個数とを計算用データ記憶部２７に記憶する。この処理動作によって、図３３に示す計算用データが生成されたことになる。

　このように、簡単な演算処理によって計算用データを生成することが可能になるため、解析作業者が行う作業負荷を軽減することができる。

＜第２の実施形態＞
　次に、本発明の第２の実施形態による計算用データ生成装置を説明する。図３７は同実施形態の構成を示すブロック図である。この図において、図２１に示す第１の実施形態による装置と同一の部分には同一の符号を付し、その説明を省略する。この図に示す装置が図２１に示す装置と異なる点は、カラーの３次元描画が行えるグラフィックボード２８が新たに設けられ、これに伴って分割領域生成部２３１の処理動作が異なる点である。

　次に、図３８を参照して、図３７に示すグラフィックボード２８の機能構成を説明する。図３８は、図３７に示すグラフィックボード２８の機能構成を示す図である。グラフィックボード２８は、縮尺調整データと三角形データを入力し、ピクセルデータを出力データとして出力する機能を有している。グラフィックボード２８内の三角形描画部は、縮尺調整データと三角形データを入力し、画像バッファに三角形を描画する。また、３次元形状である三角形の奥行き方向の情報を記憶するために、ピクセルのＺ座標（奥行き方向の座標）を記憶するＺバッファを備えている。そして、複数の三角形を描画することにより得られる画像バッファのピクセル色情報がそのまま分割領域のＩＤとして領域分割処理を実行するものである。

　縮尺調整データは、描画領域データ（解析領域データのｚ座標を除いたもの）と、ピクセル分割数データ（ボクセル分割数データのＮｚを除いたもの）と、奥行き調整データ（解析領域データのｚ座標値）とからなる。また、三角形データは、描画するべき三角形の個数分の三角形の３つの頂点の座標値からなる。

　次に、図３９を参照して、図３７に示すグラフィックボード２８の座標系について説明する。図３９は、図３７に示すグラフィックボード２８の座標系を示す図である。グラフィックボード２８内の画像バッファとＺバッファは、同じ座標系を有しており、分割領域データ生成部２３１から出力される描画領域データ（Ｘｍａｘ，Ｘｍｉｎ，Ｙｍａｘ，Ｙｍｉｎ）によって座標系が決定する。そして、各ピクセルデータは、描画領域データとピクセル分割数データとから決まり、ボクセルと同じ寸法となる。また、奥行き調整データによりＺ座標の許容最大値と許容最小値とが決定する。グラフィックボード２８において描画するべき図形は、図３９の符号Ｐで示す立体物であり、この立体物をここではポテンシャル立体と称する。

　次に、図４０を参照して、ポテンシャル立体について説明する。図４０は、ポテンシャル立体の構成を示す図である。ポテンシャル立体は、関数ｆ（ｄ）のグラフ形状を軸対象に回転した形状である。関数ｆ（ｄ）は、ｆ（ｄ）＝（Ｒ０－（Δｈ＊＊２＋ｄ＊＊２）＊＊０．５）／Ｒ０で定義される。ここで、Δｈは、初期点とＺスライス面との距離、Ｒ０は、球体の半径である。ここで球体の求め方について説明する。球体の中心座標は、初期点データ記憶部２１に記憶されている初期点の座標である。すなわち、球体は、初期点を中心として所定の半径を持つ球体である。

　半径は、次のように求める。まず、解析対象物体の内部の体積（ボクセル属性が－２であるボクセルの数にボクセルの体積を乗算した値）をＶ０とし、初期点の数をＮｐａｒｔとすると、初期点１つ分に割り当てられる平均体積ΔＶは、Ｖ０を初期点の数で除算することで求めることができる。この平均体積ΔＶが直方体とした場合、平均体積ΔＶの３乗根により直方体の１辺の長さが求められる。この１辺の長さを持つ直方体の頂点が内接する球体が求めるべき球体である。すなわち、直方体の中心から頂点までの距離が球体の半径である。

　この球体とＺスライス面における球体の断面がポテンシャル立体の底面となる。ポテンシャル立体は、山型をしており、初期点とＺスライス面との距離が短いほど山の高さが高くなうような立体である。Ｚスライス面とは、図２７に示すボクセルのＺ方向のボクセルの境界でスライスすることにより得られる面のことである。すなわち、ボクセル分割数Ｎｚと同数のＺスライス面が存在することになる。

　次に、図４１を参照して、グラフィックボード２８に入力する三角形データについて説明する。図４１は、グラフィックボード２８に入力する三角形データの構成を示す図である。三角形データは、ポテンシャル立体の表面を複数の三角形のパッチデータで表現したデータである。各三角形は、頂点１、２、３の３つの頂点の３次元座標値で表現され、これらの頂点座標値に加えて、この三角形を定義する際に用いた球体の中心に位置する初期点の初期点ＩＤが、三角形のＩＤとして付与される。すなわち、１つのポテンシャル立体を表現する三角形にはすべて同じ三角形ＩＤ（初期点ＩＤ）が付与される。

　グラフィックボード２８は、この三角形データに基づいて、描画処理を行う。ここで、図４２を参照して、グラフィックボード２８における三角形の描画処理について説明する。図４２は、グラフィックボード２８における三角形の描画処理を示す図である。まず、三角形描画部は、１個目の三角形を画像バッファに描画する。このとき、画像バッファには、三角形ＩＤの値を表示色として用いて描画する。また、三角形描画部は、三角形のＺ座標値をＺバッファに書き込む。そして、２個目の三角形を描画する場合、Ｚバッファに書き込まれている値と、新たに書き込む三角形の値と比較して、Ｚ座標値が大きい場合のみ、Ｚバッファの値を更新する。そして、Ｚバッファの更新を行ったピクセルのみを新たな三角形ＩＤの値を表示色として用いて画像バッファに書き込む。この処理を繰り返すと、隠面処理を行ったのと同等の結果が得られることになる。そして、画像バッファの内容を出力データとして出力する。このデータは、Ｚスライス面において、初期点毎の領域分割を行ったのと同等の結果を得ることができる。

　第２の実施形態におけるボクセルデータ記憶部１９にボクセルデータを記憶する処理動作と、初期点データ記憶部２１に初期点データを記憶する処理動作は、第１の実施形態における処理動作と同様であるので、詳細な処理動作の説明を省略する。以下の説明においては、第１の実施形態と異なる処理動作のみを説明する。

　次に、図４３を参照して、図３７に示す分割領域データ生成部２３１の処理動作を説明する。図４３は、図３７に示す分割領域データ生成部２３１の処理動作を示すフローチャートである。まず、分割領域データ生成部２３１は、Ｚスライス面を１つ選択する（ステップＳ５１）。そして、分割領域データ生成部２３１は、初期点を１つ選択し、この選択した初期点を中心とした球体を定義する（ステップＳ５３）。

　次に、分割領域データ生成部２３１は、この定義した球体が選択したＺスライス面と交差するか否かを判定する（ステップＳ５４）。この判定の結果、球体とＺスライス面が交差する場合、分割領域データ生成部２３１は、断面円の形状を求め（ステップＳ５５）、この断面円から前述したポテンシャル立体を定義する（ステップＳ５６）。一方、球体とＺスライス面が交差しない場合、分割領域データ生成部２３１は、ポテンシャル立体の定義（ステップＳ５５、Ｓ５６）を行わない。

　次に、分割領域データ生成部２３１は、全ての初期点について処理したかを判定し（ステップＳ５７）、全ての初期点についてポテンシャル立体の定義処理が終わるまで、ステップＳ５２～Ｓ５６の処理を繰り返す。

　次に、分割領域データ生成部２３１は、定義した全てのポテンシャル立体から定義した三角形データをグラフィックボード２８へ出力する。これを受けて、グラフィックボード２８は、ポテンシャル立体をＺスライス面に投影し、ピクセルにそれぞれの初期点ＩＤを記憶する。このとき、隠面処理を施すことにより、ポテンシャル立体同士が重なり合うときに、高い方の初期点ＩＤをピクセルに記憶する。そして、グラフィックボード２８は、全ての三角形データを描画した結果の画像バッファの内容を分割領域データ生成部２３１へ出力する。

　ここで、図４４を参照して、図３７に示すグラフィックボード２８内の三角形描画部の処理動作を説明する。図４４は、図３７に示すグラフィックボード２８内の三角形描画部の処理動作を示すフローチャートである。まず、三角形描画部は、全３角形に対する処理を繰り返し行う（ステップＳ６１）。そして、三角形描画部は、現在の３角形の３つの頂点を取り出し（ステップＳ６２）、全てのピクセルに対する処理を繰り返す（ステップＳ６３）。続いて、三角形描画部は、現在のピクセルの中心の点を通り、方向が（０，０，＋１）である直線と、三角形の面との交点座標を（Ｐｘ，Ｐｙ，Ｐｚ）とする（ステップＳ６４）。

　次に、三角形描画部は、そのピクセルは、３角形の内部または境界線上に位置しているか否かを判定する（ステップＳ６５）。この判定の結果、３角形の内部または境界線上に位置していれば、ＰｚはそのピクセルのＺバッファの値よりも大きいか否かを判定する（ステップＳ６６）。この判定の結果、ＰｚはそのピクセルのＺバッファの値よりも大きい場合、そのピクセルのＺバッファの値にＰｚを代入し、画像バッファにＩＤ番号を２４ＢＩＴの整数とした色を代入する（ステップＳ６７）。そして、全てのピクセルを処理するまで繰り返し（ステップＳ６８）、さらに、全ての３角形を処理するまで繰り返す（ステップＳ６９）。これにより、全ての三角形の描画を行うことが可能となる。この画像バッファの内容は、Ｚスライス面を初期点毎の分割を行ったのと同等の結果が得られるため、分割領域データ生成部２３１は、この結果のデータ（初期点ＩＤ）を分割領域データ記憶部２４に記憶する。

　次に、分割領域データ生成部２３１は、全てのスライス面について処理したか否かを判定し（ステップＳ５９）、全てのＺスライス面について処理が終わるまで、ステップＳ５１～Ｓ５８の処理を繰り返す。この処理によって、分割領域データ記憶部２４には、前述した第１の実施形態において得られた分割領域データと同じデータを得ることができる。計算用データ生成部２５における処理動作は、第１の実施形態における処理動作と同様であるので、詳細な説明を省略する。

　なお、第２の実施形態において、球体の半径が小さい場合、分割領域に含まれないボクセルが発生する場合があるが、この場合は、球体の半径を大きくして再度前述した処理動作を繰り返し行えばよい。また、第２の実施形態では、分割領域データ生成部２３１によって全ての分割領域データを生成してから計算用データ生成部２５によって計算用データを生成する処理動作を説明したが、１つのＺスライス面の分割領域データができた時点で、その分割領域データから計算用データを生成し、次のＺスライス面について分割領域データを生成し、この分割領域データから計算用データを生成するというように、Ｚスライス面ごとに領域分割と計算用データ生成を行うようにしてもよい。このようにすることにより、中間データ記憶部２６に記憶されるデータ量を少なくすることができ、処理を簡単にすることが可能となる。

＜第３の実施形態＞
　次に、本発明の第２の実施形態による計算用データ生成装置を説明する。図４５は同実施形態の構成を示すブロック図である。この図において、図２１に示す第１の実施形態による装置と同一の部分には同一の符号を付し、その説明を省略する。この図に示す装置が図２１に示す装置と異なる点は、初期点データ生成部２０がパッチデータ記憶部１７のデータを読み込んで初期点を生成する点である。特に、第３の実施形態における初期点データ生成部２０は、解析対象物体の壁近傍（境界層）に密な初期点を生成するのが特徴である。

　次に、図４６を参照して、図４５に示す初期点データ生成部２０の構成を説明する。図４６は、図４５に示す初期点データ生成部２０の構成を示す図である。図４６において、符号２０１は、解析対象物体の内側に存在するボクセルの内側ボクセルデータを生成する内側ボクセルデータ生成部である。符号２０２は、内側ボクセルデータ生成部２０１で生成された内側ボクセルデータを記憶する内側ボクセルデータ記憶部である。符号２０３は、壁上と球上の点データを生成する壁上点・球上点データ生成部である。符号２０４は、壁上点・球上点データ生成部２０３で生成された壁上点データを記憶する壁上点データ記憶部である。符号２０５は、壁上点・球上点データ生成部２０３で生成された球上点データを記憶する球上点データ記憶部である。符号２０６は、壁上点データと球上点データを参照して、分割点データを生成する分割点データ生成部である。符号２０７は、分割点データ生成部２０６が生成した分割点データを記憶する分割点データ記憶部である。符号２０８は、分割点データと内側ボクセルデータから初期点を生成して出力し、初期点データ記憶部２１に記憶する初期点データ出力部である。

　次に、図４７を参照して、図４６に示す内側ボクセルデータ生成部２０１の動作を説明する。まず、内側ボクセルデータ生成部２０１は、ボクセルデータ記憶部１９に記憶されているボクセルデータを参照して、ボクセルを１つ選択する（ステップＳ７１）。そして、内側ボクセルデータ生成部２０１は、パッチデータ記憶部１７を参照して、選択したボクセルとパッチが交差するかを求め（ステップＳ７２）、選択したボクセルとパッチが交差するか否かを判定する（ステップＳ７３）。この判定の結果、交差していれば内側ボクセルデータ生成部２０１は、内側ボクセルデータ記憶部２０２中の選択したボクセルの記憶領域に、－１（解析対象物体の内側のボクセルを示す値）を記憶する（ステップＳ７６）。

　一方、交差していなければ、内側ボクセルデータ生成部２０１は、ボクセルと解析対象物体の位置関係を求め（ステップＳ７４）、選択したボクセルが解析対象物体の外部領域であるか否かを判定する（ステップＳ７５）。この判定の結果、ボクセルが解析対象物体の外部に位置するボクセルであれば、内側ボクセルデータ記憶部２０２中の選択したボクセルの記憶領域に、－１（解析対象物体外部のボクセルを示す値）を記憶し（ステップＳ７６）、内部に位置するボクセルであれば、０（解析対象物体内部のボクセルを示す値）を記憶する（ステップＳ７７）。そして、内側ボクセルデータ生成部２０１は、ボクセルの全てについて処理を行ったか否かを判定し（ステップＳ７８）、全てのボクセルについてステップＳ７１～Ｓ７６の処理動作を繰り返す。この処理により、図５３に示すように、解析対象物体の内側に位置するボクセルのデータが内側ボクセルデータ記憶部２０２に記憶されたことになる。

　次に、図４８を参照して、図４６に示す壁上点・球上点データ生成部２０３の動作を説明する。まず、壁上点・球上点データ生成部２０３は、内側ボクセルデータ記憶部２０２を参照して、内側ボクセルを１つ選択し（ステップＳ８１）、隣り合うｘｙｚの６方向（正負を区別）に内側ボクセルが存在するかを調べ（ステップＳ８２）、全てに存在するかを判定する（ステップＳ８３）。この判定の結果、全てに存在していなければ、壁上点・球上点データ生成部２０３は、選択したボクセルに内接する球を定義する（ステップＳ８４）。そして、壁上点・球上点データ生成部２０３は、図５４に示すように球面上に２６点個の球上点を定義し（ステップＳ８５）、ボクセルの中心から球上点へ向けたベクトルの先に内側ボクセルが存在するかを調べ（ステップＳ８６）、存在するか否かを判定する（ステップＳ８７）。

　この判定の結果、存在しなければ壁上点・球上点データ生成部２０３は、各点から最も近い点を形状データから探索し、そこに壁上点を定義し（ステップＳ８８）、定義した壁上点、球上点をそれぞれ壁上点データ記憶部２０４、球上点データ記憶部２０５に記憶する（ステップＳ８９）。そして、壁上点・球上点データ生成部２０３は、ボクセルの全てについて処理を行ったか否かを判定し（ステップＳ９０）、全てのボクセルについてステップＳ８１～Ｓ８９の処理動作を繰り返す。この処理により、図５５に示すように、定義した壁上点、球上点をそれぞれが壁上点データ記憶部２０４、球上点データ記憶部２０５に記憶されることになる。

　次に、図４９を参照して、図４６に示す分割点データ生成部２０６の動作を説明する。まず、分割点データ生成部２０６は、壁上点データ記憶部２０４、球上点データ記憶部２０５を参照して、壁上点と球上点を通る直線を定義し（ステップＳ９１）、この直線上において、壁上点から距離ｄ１離れた位置に初期点１を生成する（ステップＳ９２）。

　次に、分割点データ生成部２０６は、定義した直線上において、初期点１から距離ｄ２離れた位置に初期点２を生成する（ステップＳ９３）。そして、分割点データ生成部２０６は、壁上点と生成された初期点とを結ぶ線分が球面と交差したか、または、ユーザが指定した数だけ直線上に初期点を生成したかを判定し（ステップＳ９４）、図５６に示すように、この条件のいずれかを満たすまで、ステップＳ９３の処理を繰り返し、初期点３、初期点４、・・・を生成し、生成した初期点群を分割点データ記憶部２０７に記憶する。なお、繰り返す条件は、予めいずれを採用するかを決めておく。この処理により、図５７に示すように、球上点と壁上点を結ぶ直線上に境界層を形成するための初期点群が生成されることになる。

　次に、図５０を参照して、図４６に示す初期点データ出力部２０８の動作を説明する。まず、初期点データ出力部２０８は、内側ボクセルデータ記憶部２０２に記憶されているボクセルデータを読み出し、このボクセルの中心点を初期点として、初期点データ記憶部２１に記憶する（ステップＳ１０１）。続いて、初期点データ出力部２０８は、分割点データ記憶部２０７に記憶されている分割点データを読み出し、この分割点データを初期点として初期点データ記憶部２１に記憶する（ステップＳ１０２）。この処理により、初期点データ記憶部２１には、解析対象物体の内側に位置するボクセルの中心点が初期点データとして記憶されるとともに、分割点が境界層を形成するための初期点データとして記憶されることになる。

　次に、図５１を参照して、図４６に示す構成の変形例を説明する。図４６に示す構成と同一の部分には同一の符号を付し、その説明を省略する。図５１に示す構成が、図４６に示す構成と異なる点は、内側ボクセル元データ記憶部２０９と内側ボクセルデータ修正部２１０が新たに設けられている点である。これは、境界層の厚さを厚くするための初期点を生成するためのものである。内側ボクセル元データ記憶部２０９は、図４６に示す内側ボクセルデータ記憶部２０２に相当し、解析対象物体の内側に位置するボクセルのボクセルデータを記憶するものである。内側ボクセルデータ修正部２１０は、境界層を厚くするために、内側ボクセル元データ記憶部２０９に記憶されているボクセルデータを境界層に近いボクセルを取り除く処理を行う。

　次に、図５２を参照して、図５１に示す内側ボクセルデータ修正部２１０の動作を説明する。まず、内側ボクセルデータ修正部２１０は、内側ボクセル元データ記憶部２０９に記憶されているボクセルデータを参照して、内側ボクセルを１つ選択し（ステップＳ１１１）、隣り合うｘｙｚの６方向（正負を区別）に内側ボクセルが存在するかを調べ（ステップＳ１１２）、全てに存在するかを判定する（ステップＳ１１３）。この判定の結果、全てに存在していなければ、内側ボクセルデータ記憶部２０２中の選択したボクセルの記憶領域に、－１を記憶する（ステップＳ１１４）。そして、内側ボクセルデータ修正部２１０は、ボクセルの全てについて処理を行ったか否かを判定し（ステップＳ１１５）、全てのボクセルについてステップＳ１１１～Ｓ１１４の処理動作を繰り返す。この処理により、解析対象物体の内側に位置するボクセルのデータが内側ボクセルデータ記憶部２０２に記憶されたことになる。ただし、境界層付近のボクセルが取り除かれている状態となる。これにより、図５８に示すように、解析対象物体の内側に位置するボクセルのデータが内側ボクセルデータ記憶部２０２に記憶されたことになる。ただし、図５３に示すボクセルデータに比べて体積（図面においては面積）の小さいボクセルデータの集合となり、内側ボクセルの集合の体積が小さくなり、その分境界層を定義する空間が大きくなる。これ以降の処理動作は、前述した動作と同様である。

　第３の実施形態により生成された初期点データを初期点データ記憶部２１に記憶し、この初期点データ記憶部２１に記憶された初期点データに基づき、分割領域データ生成部２３が前述した処理動作によって分割領域データを生成することになる。これにより、図５９に示すような分割領域データを定義することができる。このように、境界層の分割を密に行うことにより、境界層における解析精度の向上を図ることが可能となる。

　また、図２１、図３７、図４５における処理部の機能を実現するためのプログラムをコンピュータ読み取り可能な記録媒体に記録して、この記録媒体に記録されたプログラムをコンピュータシステムに読み込ませ、実行することにより計算用データ生成処理を行ってもよい。なお、ここでいう「コンピュータシステム」とは、ＯＳや周辺機器等のハードウェアを含むものとする。また、「コンピュータシステム」は、ホームページ提供環境（あるいは表示環境）を備えたＷＷＷシステムも含むものとする。また、「コンピュータ読み取り可能な記録媒体」とは、フレキシブルディスク、光磁気ディスク、ＲＯＭ、ＣＤ－ＲＯＭ等の可搬媒体、コンピュータシステムに内蔵されるハードディスク等の記憶装置のことをいう。さらに「コンピュータ読み取り可能な記録媒体」とは、インターネット等のネットワークや電話回線等の通信回線を介してプログラムが送信された場合のサーバやクライアントとなるコンピュータシステム内部の揮発性メモリ（ＲＡＭ）のように、一定時間プログラムを保持しているものも含むものとする。

　また、上記プログラムは、このプログラムを記憶装置等に格納したコンピュータシステムから、伝送媒体を介して、あるいは、伝送媒体中の伝送波により他のコンピュータシステムに伝送されてもよい。ここで、プログラムを伝送する「伝送媒体」は、インターネット等のネットワーク（通信網）や電話回線等の通信回線（通信線）のように情報を伝送する機能を有する媒体のことをいう。また、上記プログラムは、前述した機能の一部を実現するためのものであってもよい。さらに、前述した機能をコンピュータシステムにすでに記録されているプログラムとの組み合わせで実現できるもの、いわゆる差分ファイル（差分プログラム）であってもよい。

　本出願を詳細にまた特定の実施態様を参照して説明したが、本発明の精神と範囲を逸脱することなく様々な変更や修正を加えることができることは当業者にとって明らかである。
　本出願は、２０１０年８月２４日出願の日本特許出願（特願２０１０-１８７３０９）及び２０１０年１２月２１日出願の日本特許出願（特願２０１０-２８４８７８）に基づくものであり、その内容はここに参照として取り込まれる。

　従来の数値解析手法である有限要素法、有限体積法、ボクセル法、ボクセル法の改良手法、及び粒子法の問題点を解決し、解析精度の悪化を伴わずにソルバ処理における計算負荷の低減を図ることができる数値解析装置に入力する計算用データを生成することが不可欠な用途にも適用できる。

Ａ・・・数値解析装置（物理量計算装置）
Ｐ・・・数値解析プログラム
Ｐ１・・・プリ処理プログラム
Ｐ２・・・ソルバ処理プログラム（物理量計算プログラム）
Ｐ３・・・ポスト処理プログラム
Ｍ・・・計算用データモデル
Ｖ_ｉ……粒子iが仮想的に占める体積
Ｓ_ｉｊ……粒子iと粒子jの境界面積
ｒ_ｉｊ……ｉからjの距離ベクトル
ｎ_ｉｊ……Ｓ_ｉｊの法線ベクトル
１・・・ＣＰＵ
２・・・記憶装置
２ａ・・・プログラム記憶部
２ｂ・・・データ記憶部
３・・・ＤＶＤドライブ
４・・・入力装置
４ａ・・・キーボード
４ｂ・・・マウス
５・・・出力装置
５ａ・・・ディスプレイ
５ｂ・・・プリンタ
Ｄ・・・計算用データ生成装置
１１・・・パラメータ入力部
１２・・・パラメータ記憶部
１３・・・通信部
１４・・・形状データ入力部
１５・・・形状データ記憶部
１６・・・パッチデータ生成部
１７・・・パッチデータ記憶部
１８・・・ボクセルデータ生成部
１９・・・ボクセルデータ記憶部
２０・・・初期点データ生成部
２１・・・初期点データ記憶部
２２・・・データ間引き部
２３、２３１・・・分割領域データ生成部
２４・・・分割領域データ記憶部
２５・・・計算用データ生成部
２６・・・中間データ記憶部
２７・・・計算用データ記憶部
２８・・・グラフィックボード

Claims

　数値解析を行う対象物体を複数の多角形で表現したパッチデータが記憶されたパッチデータ記憶手段と、
　前記数値解析を行うのに必要なパラメータを記憶するパラメータ記憶手段と、
　前記対象物体を含む解析領域内を複数の直方体に分割したボクセルデータを記憶するボクセルデータ記憶手段と、
　前記パラメータ記憶手段に記憶されている前記パラメータに基づき、前記ボクセルデータを定義し、各ボクセルに対して、ボクセル属性を付与して、前記ボクセルデータ記憶手段に記憶するボクセルデータ生成手段と、
　前記解析領域内の領域分割を行うための初期点データを記憶する初期点データ記憶手段と、
　前記ボクセルの中心点を使用して、前記ボクセルの数より少ない前記初期点データを生成して、前記初期点データ記憶手段に記憶する初期点データ生成手段と、
　前記対象物体を複数の分割領域に分割した分割領域データを記憶する分割領域データ記憶手段と、
　前記ボクセルデータ記憶手段に記憶されている各ボクセルの属性と、前記初期点データ記憶手段に記憶されている前記初期点データとに基づき、前記対象物体内を、複数の前記ボクセルで構成する分割領域を定義し、定義した前記分割領域データを前記分割領域データ記憶手段に記憶する分割領域データ生成手段と、
　前記数値解析を行うための計算用データを記憶する計算用データ記憶手段と、
　前記分割領域データ記憶手段に記憶されている前記分割領域データに基づき、各分割領域の境界面データを生成し、該境界面データを前記計算用データとして前記計算用データ記憶手段に記憶する計算用データ生成手段と
　を備えることを特徴とする計算用データ生成装置。
　前記分割領域データ生成手段は、
　前記初期点データ記憶手段に記憶されている前記初期点データの中から、各ボクセルの中心点に最も距離が近い初期点を選択して、各ボクセルが属する分割領域を特定することにより前記分割領域を生成することを特徴とする請求項１に記載の計算用データ生成装置。
　前記分割領域データ生成手段は、
　所定の軸に垂直な前記ボクセルの境界面で構成するスライス面と、前記初期点のそれぞれを中心とする球体の交差する断面を定義し、該断面から前記初期点からの距離に応じて高さの異なるポテンシャル立体を定義し、該ポテンシャル立体を３次元の隠面処理を行って描画した画像データに基づいて、前記スライス面上の分割領域を定義する処理を、前記解析領域内の全ての前記スライス面に対して行うことにより前記分割領域を生成することを特徴とする請求項１に記載の計算用データ生成装置。
　前記初期点データ生成手段は、
　前記パッチデータ記憶手段に記憶されている前記パッチデータと、前記ボクセルデータ記憶手段に記憶されている前記ボクセルデータとから、前記対象物体内部に位置する前記ボクセルの内側ボクセルを選択し、前記内側ボクセルのうち、隣接するボクセルが存在しない内側ボクセルに内接する球を定義し、前記球上に球上点を定義するともに、前記対象物体の壁上に壁上点を定義し、対となる前記球上点と前記壁上点との間に分割点を定義し、該分割点と前記ボクセルの中心点とを前記初期点とすることを特徴とする請求項１から３のいずれかに記載の計算用データ生成装置。
　数値解析を行う対象物体を複数の多角形で表現したパッチデータが記憶されたパッチデータ記憶手段と、前記数値解析を行うのに必要なパラメータを記憶するパラメータ記憶手段と、前記対象物体を含む解析領域内を複数の直方体に分割したボクセルデータを記憶するボクセルデータ記憶手段と、前記解析領域内の領域分割を行うための初期点データを記憶する初期点データ記憶手段と、前記対象物体を複数の分割領域に分割した分割領域データを記憶する分割領域データ記憶手段と、前記数値解析を行うための計算用データを記憶する計算用データ記憶手段とを備える計算用データ生成装置における計算用データ生成方法であって、
　前記パラメータ記憶手段に記憶されている前記パラメータに基づき、前記ボクセルデータを定義し、各ボクセルに対して、ボクセル属性を付与して、前記ボクセルデータ記憶手段に記憶するボクセルデータ生成ステップと、
　前記ボクセルの中心点を使用して、前記ボクセルの数より少ない前記初期点データを生成して、前記初期点データ記憶手段に記憶する初期点データ生成ステップと、
　前記ボクセルデータ記憶手段に記憶されている各ボクセルの属性と、前記初期点データ記憶手段に記憶されている前記初期点データとに基づき、前記対象物体内を、複数の前記ボクセルで構成する分割領域を定義し、定義した前記分割領域データを前記分割領域データ記憶手段に記憶する分割領域データ生成ステップと、
　前記分割領域データ記憶手段に記憶されている前記分割領域データに基づき、各分割領域の境界面データを生成し、該境界面データを前記計算用データとして前記計算用データ記憶手段に記憶する計算用データ生成ステップと
　を有することを特徴とする計算用データ生成方法。
　前記分割領域データ生成ステップは、
　前記初期点データ記憶手段に記憶されている前記初期点データの中から、各ボクセルの中心点に最も距離が近い初期点を選択して、各ボクセルが属する分割領域を特定することにより前記分割領域を生成することを特徴とする請求項５に記載の計算用データ生成方法。
　前記分割領域データ生成ステップは、
　所定の軸に垂直な前記ボクセルの境界面で構成するスライス面と、前記初期点のそれぞれを中心とする球体の交差する断面を定義し、該断面から前記初期点からの距離に応じて高さの異なるポテンシャル立体を定義し、該ポテンシャル立体を３次元の隠面処理を行って描画した画像データに基づいて、前記スライス面上の分割領域を定義する処理を、前記解析領域内の全ての前記スライス面に対して行うことにより前記分割領域を生成することを特徴とする請求項５に記載の計算用データ生成方法。
　前記初期点データ生成ステップは、
　前記パッチデータ記憶手段に記憶されている前記パッチデータと、前記ボクセルデータ記憶手段に記憶されている前記ボクセルデータとから、前記対象物体内部に位置する前記ボクセルの内側ボクセルを選択し、前記内側ボクセルのうち、隣接するボクセルが存在しない内側ボクセルに内接する球を定義し、前記球上に球上点を定義するともに、前記対象物体の壁上に壁上点を定義し、対となる前記球上点と前記壁上点との間に分割点を定義し、該分割点と前記ボクセルの中心点とを前記初期点とすることを特徴とする請求項５から７のいずれかに記載の計算用データ生成方法。
　数値解析を行う対象物体を複数の多角形で表現したパッチデータが記憶されたパッチデータ記憶手段と、前記数値解析を行うのに必要なパラメータを記憶するパラメータ記憶手段と、前記対象物体を含む解析領域内を複数の直方体に分割したボクセルデータを記憶するボクセルデータ記憶手段と、前記解析領域内の領域分割を行うための初期点データを記憶する初期点データ記憶手段と、前記対象物体を複数の分割領域に分割した分割領域データを記憶する分割領域データ記憶手段と、前記数値解析を行うための計算用データを記憶する計算用データ記憶手段とを備える計算用データ生成装置上のコンピュータに計算用データ生成を行わせる計算用データ生成プログラムであって、
　前記パラメータ記憶手段に記憶されている前記パラメータに基づき、前記ボクセルデータを定義し、各ボクセルに対して、ボクセル属性を付与して、前記ボクセルデータ記憶手段に記憶するボクセルデータ生成ステップと、
　前記ボクセルの中心点を使用して、前記ボクセルの数より少ない前記初期点データを生成して、前記初期点データ記憶手段に記憶する初期点データ生成ステップと、
　前記ボクセルデータ記憶手段に記憶されている各ボクセルの属性と、前記初期点データ記憶手段に記憶されている前記初期点データとに基づき、前記対象物体内を、複数の前記ボクセルで構成する分割領域を定義し、定義した前記分割領域データを前記分割領域データ記憶手段に記憶する分割領域データ生成ステップと、
　前記分割領域データ記憶手段に記憶されている前記分割領域データに基づき、各分割領域の境界面データを生成し、該境界面データを前記計算用データとして前記計算用データ記憶手段に記憶する計算用データ生成ステップと
　を前記コンピュータに行わせることを特徴とする計算用データ生成プログラム。
　前記分割領域データ生成ステップは、
　前記初期点データ記憶手段に記憶されている前記初期点データの中から、各ボクセルの中心点に最も距離が近い初期点を選択して、各ボクセルが属する分割領域を特定することにより前記分割領域を生成することを特徴とする請求項９に記載の計算用データ生成プログラム。
　前記分割領域データ生成ステップは、
　所定の軸に垂直な前記ボクセルの境界面で構成するスライス面と、前記初期点のそれぞれを中心とする球体の交差する断面を定義し、該断面から前記初期点からの距離に応じて高さの異なるポテンシャル立体を定義し、該ポテンシャル立体を３次元の隠面処理を行って描画した画像データに基づいて、前記スライス面上の分割領域を定義する処理を、前記解析領域内の全ての前記スライス面に対して行うことにより前記分割領域を生成することを特徴とする請求項９に記載の計算用データ生成プログラム。
　前記初期点データ生成ステップは、
　前記パッチデータ記憶手段に記憶されている前記パッチデータと、前記ボクセルデータ記憶手段に記憶されている前記ボクセルデータとから、前記対象物体内部に位置する前記ボクセルの内側ボクセルを選択し、前記内側ボクセルのうち、隣接するボクセルが存在しない内側ボクセルに内接する球を定義し、前記球上に球上点を定義するともに、前記対象物体の壁上に壁上点を定義し、対となる前記球上点と前記壁上点との間に分割点を定義し、該分割点と前記ボクセルの中心点とを前記初期点とすることを特徴とする請求項９から１１のいずれかに記載の計算用データ生成プログラム。