JPWO2020054087A1 - シミュレーション方法、物理量計算プログラム及び物理量計算装置 - Google Patents

シミュレーション方法、物理量計算プログラム及び物理量計算装置 Download PDF

Info

Publication number
JPWO2020054087A1
JPWO2020054087A1 JP2019506743A JP2019506743A JPWO2020054087A1 JP WO2020054087 A1 JPWO2020054087 A1 JP WO2020054087A1 JP 2019506743 A JP2019506743 A JP 2019506743A JP 2019506743 A JP2019506743 A JP 2019506743A JP WO2020054087 A1 JPWO2020054087 A1 JP WO2020054087A1
Authority
JP
Japan
Prior art keywords
region
calculation
divided
area
vertices
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
JP2019506743A
Other languages
English (en)
Other versions
JP6504333B1 (ja
Inventor
恒洋 齊藤
恒洋 齊藤
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
AGC Inc
Original Assignee
Asahi Glass Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Asahi Glass Co Ltd filed Critical Asahi Glass Co Ltd
Priority claimed from PCT/JP2018/043836 external-priority patent/WO2020054087A1/ja
Application granted granted Critical
Publication of JP6504333B1 publication Critical patent/JP6504333B1/ja
Publication of JPWO2020054087A1 publication Critical patent/JPWO2020054087A1/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13—Differential equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Computing Systems (AREA)
  • Operations Research (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

【課題】物理現象を数値的に解析する数値解析において、ソルバ処理に要する時間の削減が可能な技術を提供すること。【解決手段】コンピュータが、解析領域を複数の分割領域に分割し、各分割領域の体積と隣り合う分割領域同士の特性を示す分割領域特性量とを分割領域の頂点の座標及び該頂点の連結情報を必要としない量として有する分割領域での計算用データモデルを生成し、分割領域を複数集合させることによって、要求される数の集合領域を生成し、各集合領域の体積と隣り合う集合領域同士の特性を示す集合領域特性量とを集合領域の頂点の座標及び該頂点の連結情報を必要としない量として有する集合領域での計算用データモデルを生成し、解析領域での物性値と、集合領域での計算用データモデルとに基づいて、集合領域での解析結果である物理量を計算する、ことを特徴とするシミュレーション方法。

Description

本発明は、シミュレーション方法、物理量計算プログラム及び物理量計算装置に関する。
従来から、流速分布、応力分布及び温度分布等を数値解析によって求めるための数値解析手法として、例えば、有限要素法、有限体積法、ボクセル法及び粒子法が知られている。
しかし、このような従来の数値解析手法は、大変よく知られているように、十分な解析精度を得ようとすると、計算用データモデルの作成とソルバ処理とにおいて、膨大な作業及び時間を必要とすることが問題であった。
このような問題を解決するため、特許文献１の数値解析の方法が提案されている。特許文献１の方法は、従来の数値解析手法にとって必須であったメッシュを必要としない。また、特許文献１の方法は、解析対象の物理現象における物理量の保存則を満たしつつ、物理現象を数値的に解析可能である。さらに、特許文献１の方法は、十分な解析精度を得つつ、計算用データモデルの生成に必要な作業及び時間を軽減することができる。
特許文献１の方法においては、十分な解析精度と計算用データモデルの生成に必要な作業の軽減とを維持して、ソルバ処理に要する時間のいっそうの削減が期待される。
国際公開第２０１０／１５０７５８号
本発明が解決しようとする課題は、物理現象を数値的に解析する数値解析において、ソルバ処理に要する時間の削減が可能な技術を提供することである。
本発明の一態様は、コンピュータによって物理現象での物理量を数値的に解析するシミュレーション方法であって、コンピュータが、解析領域を複数の分割領域に分割し、前記分割領域の頂点の座標（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とを持たず、集合領域の体積と、集合領域の境界面特性量と、その他補助データ（例えば、後述する、集合領域の係合情報や集合コントロールポイントの座標等）とを有するものとなる。
このように本実施形態を用いた場合には、前述のように、幾何学的形状を規定する量を必要としない量である集合領域の体積と上記境界面特性量に基づいて各領域における物理量が計算可能とされる。このため、第二計算用データモデルに、集合領域の幾何学的形状を規定する量を持たせることなく物理量を算出することが可能となる。したがって、本実施形態を用いることにより、プリ処理において、少なくとも集合領域の体積と境界面特性量（境界面の面積及び境界面の法線ベクトル）とを有する第二計算用データモデルを作成すれば良くなり、幾何学的形状を規定する量を有する計算用データモデルを作成することなく物理量の計算を行うことが可能となる。
幾何学的形状を規定する量を持たない第二計算用データモデルは、集合領域の幾何学的形状を規定する量を必要としないため、集合領域の幾何学的形状を規定する量に縛られることなく作成することができる。
このため、３次元形状データの修正作業に対する規制も大幅に緩和される。よって、幾何学的形状を規定する量を持たない第二計算用データモデルは、幾何学的形状を規定する量を有する計算用データモデルと比較して遥かに容易に作成できる。したがって、本実施形態によれば、計算用データモデルの作成における作業負担を軽減できる。
また、本実施形態を用いる場合であっても、プリ処理においては、幾何学的形状を規定する量を使用しても構わない。つまり、プリ処理において幾何学的形状を規定する量を用いて分割領域の体積や境界面特性値等を算出しても良い。このような場合であっても、ソルバ処理においては集合領域の体積や境界面特性値があれば物理量の計算が可能であるため、プリ処理において幾何学的形状を規定する量を利用するにしても、集合領域の幾何学的形状に対する制約、例えば分割領域の歪みや捩じれ等に起因する制約がなく、計算用データモデルの作成における作業負担の軽減が可能である。
また、本実施形態を用いることによって、プリ処理において、集合領域の幾何学的形状に対する制約がなくなるため、集合領域を任意の形状に変更させることができる。このため、集合領域の数を増やすことなく、解析領域を実際に解析したい領域に容易にフィッティングさせることが可能となり、計算負荷を増大させることなく解析精度を向上させることが可能となる。
さらに、本実施形態を用いることによって集合領域の分布密度も任意に変更することができるため、必要な範囲で計算負荷の増大を許容しながらさらに解析精度を向上させることも可能である。
本実施形態では、プリ処理において、まず、任意に配置された分割領域の体積、境界面特性量（境界面の面積、境界面の法線ベクトル）及び物理量のやり取りを行う分割領域同士を関連付ける情報を有する分割領域の第一計算用データモデルを作成する。例えば、この物理量のやり取りを行う分割領域同士を関連付ける情報は、隣り合う分割領域同士の結合情報（link）と、隣り合う分割領域同士の距離からなる。そして、このlinkによって関連付けられる分割領域は、必ずしも空間的に隣接されている必要はなく、空間的に離間していても構わない。このようなlinkは、幾何学的形状を規定する量と関連するものではなく、幾何学的形状を規定する量と比較すると、極めて短時間で作成できる。また、後に詳説するが、本実施形態では、必要に応じて、分割領域の内部に配置されるコントロールポイントの座標を、第一計算用データモデルに持たせる場合もある。次に、本実施形態では、プリ処理において、分割領域を複数集合させることによって、要求される数の集合領域を生成する。そして、集合領域の体積、境界面特性量（境界面の面積、境界面の法線ベクトル）及び物理量のやり取りを行う集合領域同士を関連付ける情報を有する第二計算用データモデルを作成する。この物理量のやり取りを行う集合領域同士を関連付ける情報、隣り合う集合領域同士の結合情報（link）も、分割領域に対する特徴と同様の特徴を有する。また、後に詳説するが、本実施形態では、必要に応じて、集合領域の内部に配置されるコントロールポイントの座標を、第二計算用データモデルに持たせる場合もある。
そして、本実施形態では、プリ処理から、分割領域の体積、境界面特性量及びlink並びにコントロールポイントの座標等を有する分割領域の第一計算用データモデルと、境界条件や初期条件等とをソルバ処理に受け渡す。ソルバ処理では、受け渡されたその第一計算用データモデルに含まれる分割領域の体積、境界面特性量等を使用して離散化支配方程式を解くことによって物理量の算出を行う。
また、本実施形態では、プリ処理から、集合領域の体積、境界面特性量及びlink並びにコントロールポイントの座標等を有する第二計算用データモデルと、境界条件や初期条件等とをソルバ処理に受け渡す。ソルバ処理では、受け渡されたその第二計算用データモデルに含まれる集合領域の体積、境界面特性量等を使用して離散化支配方程式を解くことによって物理量の算出を行う。
そして、本実施形態では、ソルバ処理において、幾何学的形状を規定する量を使用しないで物理量を計算している点が従来の有限体積法と大きく異なり、この点が本実施形態の大きな特徴である。このような特徴は、ソルバ処理にて、幾何学的形状を規定する量を必要としない量のみを使用する離散化支配方程式を用いることによって得られる。
この結果、本実施形態では、ソルバ処理に幾何学的形状を規定する量を受け渡す必要がなくなり、プリ処理において、幾何学的形状を規定する量を持たない計算用データモデルを作成すればよい。したがって、従来の有限体積法と比較して本実施形態では、遥かに容易に計算用データモデルを作成でき、計算用データモデルの作成における作業負担を軽減できる。
まず、本実施形態を用いた数値解析手法における数値解析の処理の流れを簡単に説明する。
前述したように、本実施形態を用いた数値解析手法は、解析領域を複数の分割領域に分割する処理（以下「複数の分割領域に分割する処理」という）を有する。さらに、本数値解析手法は、分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された分割領域での支配方程式に基づき、各分割領域の体積と隣り合う分割領域同士の特性を示す分割領域特性量とを分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する第一計算用データモデルを生成する処理（以下「分割領域での第一計算用データモデルを生成する処理」という）を有する。
さらに、本数値解析手法は、分割領域を複数集合させることによって、要求される数の集合領域を生成する処理（以下「集合領域を生成する処理」という）を有する。さらに、本数値解析手法は、集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された集合領域での支配方程式に基づき、各集合領域の体積と隣り合う集合領域同士の特性を示す集合領域特性量とを集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する第二計算用データモデルを生成する処理（以下「集合領域での第二計算用データモデルを生成する処理」という）を有する。さらに、本数値解析手法は、解析領域での物性値と、集合領域での計算データモデルとに基づいて、物理量を計算する処理を有する。
（複数の分割領域に分割する処理）
第一計算用データモデルを作成するための解析領域を複数の分割領域に分割する処理について説明する。解析領域を複数の分割領域に分割する処理では、解析領域を幾何学的形状を規定する量を使用しないセルで細かく分割する。
（分割領域での第一計算用データモデルを生成する処理）
分割領域での第一計算用データモデルを生成する処理について説明する。なお、この分割領域での第一計算用データモデルを生成する処理は、特許文献１で開示されている処理と同様である。
図１において、セルＲ１，Ｒ２，Ｒ３・・・は、解析領域を分割して得られる分割領域であり、各々が体積Ｖａ，Ｖｂ，Ｖｃ・・・を有している。また、境界面Ｅは、セルＲ１とセルＲ２との間において物理量の交換が行われる面であり、本実施形態における境界面に相当する。また、面積Ｓａｂは、境界面Ｅの面積を示し、本実施形態における境界面特性量の１つである。また、［ｎ］ａｂは、境界面Ｅの法線ベクトルを示し、本実施形態における境界面特性量の１つである。
また、コントロールポイントａ，ｂ，ｃ・・・は、各セルＲ１，Ｒ２，Ｒ３の内部に配置されており、図１においては各セルＲ１，Ｒ２，Ｒ３・・・の重心位置に配置されている。ただし、コントロールポイントａ，ｂ，ｃ・・・は、必ずしも各セルＲ１，Ｒ２，Ｒ３・・・の重心位置に配置される必要はない。また、αは、コントロールポイントａからコントロールポイントｂまでの距離を１とした場合におけるコントロールポイントａから境界面Ｅまでの距離を示し、境界面Ｅがコントロールポイントａとコントロールポイントｂとを結ぶ線分のどの内分点に存在するかを示す比率である。なお、コントロールポイントａからコントロールポイントｂまでの距離は、隣り合う分割領域同士の距離の一例である。
なお、境界面は、セルＲ１とセルＲ２との間のみに限らず、隣り合う全てのセル間に存在する。そして、境界面の法線ベクトル及び境界面の面積も、境界面ごとに与えられる。
そして、実際の分割領域での第一計算用データモデルは、各コントロールポイントａ，ｂ，ｃ・・・の配置データと、各コントロールポイントａ，ｂ，ｃ・・・が存在するセルＲ１，Ｒ２，Ｒ３・・・の体積Ｖａ，Ｖｂ，Ｖｃ・・・を示す体積データと、各境界面の面積を示す面積データと、各境界面の法線ベクトル（以下「法線ベクトル」という。）を示す法線ベクトルデータとを有するデータ群として構築されている。
このことは、本数値解析手法の分割領域での第一計算用データモデルは、セルＲ１，Ｒ２，Ｒ３・・・の体積Ｖａ，Ｖｂ，Ｖｃ・・・と、隣り合うセルＲ１，Ｒ２，Ｒ３・・・同士の境界面の特性を示す境界面特性量である境界面の面積と、隣り合うセルＲ１，Ｒ２，Ｒ３・・・同士の境界面の特性を示す境界面特性量である境界面の法線ベクトルとを有して定義されていることを意味する。
なお、各セルＲ１，Ｒ２，Ｒ３・・・は、コントロールポイントａ，ｂ，ｃ・・・を有している。このため、セルＲ１，Ｒ２，Ｒ３・・・の体積Ｖａ，Ｖｂ，Ｖｃ・・・は、コントロールポイントａ，ｂ，ｃ・・・が仮想的に占める空間（コントロールボリューム）の体積として捉えることができる。
また、本数値解析手法の第一計算用データモデルは、必要に応じて、境界面が挟まれたコントロールポイント同士を結ぶ線分のどの内分点に存在するかの比率αを示す比率データを有する。
以下では、前述の第一計算用データモデルを用いて解析領域の各セルにおける流速を求める物理量計算例について説明する。なお、ここでは、各コントロールポイントにおける流速を各セルにおける流速として求める。
まず、本物理量計算において本数値解析手法は、流体解析の場合、下式（１）で示すナビエ・ストークスの式と、下式（２）で示す連続の式とを用いる。
なお、式（１），（２）において、ｔが時間を示し、ｘｉ（ｉ＝１，２，３）がカーテシアン系における座標を示し、ρが流体密度を示し、ｕｉ（ｉ＝１，２，３）が流体の流速成分を示し、Ｐが圧力を示し、μが流体の粘性係数を示し、添字ｉ（ｉ＝１，２，３），ｊ（ｊ＝１，２，３）がカーテシアン座標系における各方向成分を示している。また、添字ｊに関しては総和規約に従うものとする。
そして、式（１），（２）を、重み付き残差積分法に基づいて、分割領域のコントロールボリュームの体積に対して積分して示すと、式（１）が下式（３）のように示され、式（２）が下式（４）のように示される。
なお、式（３），（４）において、Ｖがコントロールボリュームの体積を示し、∫ＶｄＶが体積Ｖに関する積分を示し、Ｓがコントロールボリュームの面積を示し、∫ＳｄＳが面積Ｓに関する積分を示し、［ｎ］がＳの法線ベクトルを示し、ｎｉ（ｉ＝１，２，３）が法線ベクトル［ｎ］の成分を示し、∂／∂ｎが法線方向微分を示している。
ここで、説明を簡単化するために、流体の密度ρと粘性係数μを定数とする。ただし、以下の定数化は、流体の物性値が、時間、空間、温度等によって変化する場合に対して拡張可能である。
そして、図１のコントロールポイントａについて、境界面Ｅの面積Ｓａｂについて離散化し、代数方程式による近似式に変換すると、式（３）が下式（５）、式（４）が下式（６）のように示される。
ここで、添字ａｂが付く、［ｎ］ａｂ，［ｕ］ａｂ，ｕｉａｂ，ｎｉａｂ，Ｐａｂ，（∂ｕｉ／∂ｎ）ａｂは、コントロールポイントａとコントロールポイントｂとの間の境界面Ｅ上における物理量であることを示す。また、ｎｉａｂは［ｎ］ａｂの成分である。また、ｍは、コントロールポイントａと結合関係（境界面を挟む関係）にある全てのコントロールポイントの数である。
そして、式（５），（６）をＶａ（コントロールポイントａのコントロールボリュームの体積）で割ると、式（５）が下式（７）のように示され、式（６）が下式（８）のように示される。
ここで、下式（９）とする。
すると、式（７）が下式（１０）のように示され、式（８）が下式（１１）のように示される。
式（１０），（１１）において、［ｕ］ａｂ、ｕｉａｂ、Ｐａｂ、（∂ｕｉ／∂ｎ）ａｂは、コントロールポイントａとコントロールポイントｂ上の物理量の重み付け平均（移流項については、風上性を考慮した重み付け平均）により近似的に求められ、コントロールポイントａ，ｂ間の距離及び向きと、その間に存在する境界面Ｅとの位置関係（上記比率α）と、境界面Ｅの法線ベクトルの向きに依存して決定される。ただし、［ｕ］ａｂ、ｕｉａｂ、Ｐａｂ、（∂ｕｉ／∂ｎ）ａｂは、境界面Ｅの幾何学的形状を規定する量には無関係な量である。 また、式（９）で定義されるφａｂも（面積／体積）という量であり、コントロールボリュームの幾何学的形状を規定する量には無関係な量である。
つまり、このような式（１０）、（１１）は、セル形状を規定する幾何学的形状を規定する量を必要としない量のみを使用して物理量が算出可能な、重み付き残差積分法に基づく演算式である。
このため、物理量計算（ソルバ処理）に先立って前述の第一計算用データモデルを作成し、物理量計算において当該第一計算用データモデルと、式（１０）、（１１）の離散化支配方程式とを用いることによって、物理量計算においてコントロールボリュームの幾何学的形状を規定する量を全く使用せずに、流速の計算を行うことができる。
このように、物理量計算において幾何学的形状を規定する量を全く使用せずに流速の計算が行えることから、第一計算用データモデルに幾何学的形状を規定する量を持たせる必要がなくなる。よって、第一計算用データモデルの作成にあたり、セルの幾何学的形状に縛られる必要がなくなるため、セルの形状を任意に設定することができる。このため、本数値解析手法によれば、前述のように３次元形状データの修正作業に対する規制を大幅に緩和することができる。
なお、実際に式（１０）、（１１）を解くにあたり、［ｕ］ａｂ、ｕｉａｂ、Ｐａｂ等の境界面Ｅ上の物理量は、通常、線形補間によって補間される。例えば、コントロールポイントａの物理量をψａ、コントロールポイントｂの物理用をψｂとすると、境界面Ｅ上の物理量ψａｂは、下式（１２）によって求めることができる。
また、物理量ψａｂは、境界面が挟まれたコントロールポイント同士を結ぶ線分のどの内分点に存在するかの比率αを用いることによって、下式（１３）によって求めることもできる。
したがって、第一計算用データモデルが比率αを示す比率データを有している場合には、式（１３）を用いて境界面Ｅ上の物理量をコントロールポイントａとコントロールポイントｂとからの離間する距離に応じた重み付け平均を用いて算出することができる。
また、連続体モデルの方程式（ナビエ・ストークスの式等）には、式（１）に示すように、１階の偏導関数（偏微分）が含まれる。
ここで、連続体モデルの方程式の微係数を部分積分、ガウスの発散定理、あるいは一般化されたグリーンの定理を利用して、体積分を面積分に変換し、微分の次数を下げる。これによって１次微分は０次微分（スカラー量またはベクトル量）とすることができる。
例えば一般化されたグリーンの定理では、物理量をψとすると、下式（１４）という関係が成り立つ。
なお、式（１４）において、ｎｉ（ｉ＝１、２、３）は、表面Ｓ上の単位法線ベクトル［ｎ］のｉ方向の成分である。
連続体モデルの方程式の１次微分項は、体積分から面積分の変換により、境界面上ではスカラー量またはベクトル量として取り扱われる。そして、これらの値は、前述の線形補間等によって、各コントロールポイント上の物理量から補間できる。
また、連続体モデルの方程式によっては、２階の偏導関数が含まれる場合もある。
式（１４）の被積分関数をさらに１階微分した式は下式（１５）となり、連続体モデルの方程式の２次微分項は、体積分から面積分の変換により境界面Ｅ上では下式（１６）となる。
なお、式（１５）において∂／∂ｎは法線方向微分を示し、式（１６）において∂／∂ｎａｂは、［ｎ］ａｂ方向微分を示す。
つまり、連続体モデルの方程式の２次微分項は、体積分から面積分の変換により、物理量ψの法線方向微分（Ｓａｂの法線［ｎ］ａｂ方向への微分）に、［ｎ］の成分ｎｉａｂ、ｎｊａｂを乗じた形となる。
ここで式（１６）中の∂ψ／∂ｎａｂは、下式（１７）と近似される。
なお、コントロールポイントａとコントロールポイントｂとのコントロールポイント間ベクトル［ｒ］ａｂは、コントロールポイントａの位置ベクトル［ｒ］ａとコントロールポイントｂの位置ベクトル［ｒ］ｂから下式（１８）のように定義される。
したがって、境界面Ｅの面積がＳａｂであるため、式（１６）は下式（１９）となり、これを利用して式（１６）を計算できる。
なお、式（１６）の導出にあたり、次のことがわかる。
すべての線形偏微分方程式は、定数と、１次、２次、その他の偏導関数に係数を乗じた項の線形和で表わされる。式（１５）から式（１９）において、物理量ψをψの１次偏導関数に置き換えると、より高次の偏導関数の体積分を、式（１４）のように低次の偏導関数の面積分により求めることができる。この手順を、低次の偏微分から順次繰り返すと、線形偏微分方程式を構成するすべての項の偏導関数は、コントロールポイントの物理量ψと、式（１２）又は式（１３）で計算される境界面上のψであるψａｂと、式（１８）で定義されるコントロールポイント間ベクトルから求められるコントロールポイント間距離と、境界面Ｅの面積Ｓａｂ、境界面Ｅの法線ベクトルの成分ｎｉａｂとｎｊａｂから、すべて求めることができる。
本数値解析手法において物理量計算にあたり幾何学的形状を規定する量を必要としないことは前述した。このため、第一計算用データモデルの作成にあたり、コントロールボリュームの体積と、境界面の面積及び法線ベクトルとを、幾何学的形状を規定する量を使用しないで求めれば、式（１０）と式（１１）との離散化支配方程式を用いて、コントロールボリュームの幾何学的形状であるセルの幾何学的形状を全く使用せずに、流速の計算を行うことができる。
ただし、本数値解析手法においては、必ずしも、コントロールボリュームの体積と、境界面の面積及び法線ベクトルとを、コントロールボリュームの具体的な幾何学的形状を使用しないで求める必要はない。このように、ソルバ処理において幾何学的形状を規定する量を利用しないので、コントロールボリュームの具体的な幾何学的形状、具体的にはVertexとConnectivityとを利用するとしても、従来の有限要素法、有限体積法のような分割領域に関わる制約、すなわち分割領域の歪みや捩じれに対する制約がないため、前述のように容易に計算用データモデルの作成ができる。
なお、前述の説明においては、ナビエ・ストークスの式及び連続の式から重み付き残差積分法に基づいて導出した離散化支配方程式を用いる物理量の計算例について説明したが、本数値解析手法において用いられる離散化支配方程式はこれに限られるものではない。
つまり、種々の方程式（質量保存の方程式、運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式、及び波動方程式等）から重み付き残差積分法に基づいて導出されると共に、幾何学的形状を規定する量を必要としない量のみを使用して物理量を算出可能な離散化支配方程式であれば本数値解析手法に用いることができる。
そして、このような離散化支配方程式の特性によって、従来の有限要素法や有限体積法のようにいわゆるメッシュを必要としない、メッシュレスでの計算が可能となる。また、たとえ、プリ処理において、セルの幾何学的形状を規定する量を利用するとしても、従来の有限要素法、有限体積法、ボクセル法のようなメッシュに対する制約がないため、第一計算用データモデルの作成に伴う作業負荷を軽減できる。
本実施形態では、質量保存の方程式、運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式、及び波動方程式から、重み付き残差積分法に基づいて、幾何学的形状を規定する量を必要としない量のみを使用する離散化支配方程式が導出可能である。そのため、本数値解析手法において他の支配方程式を用いることができる。これについては、特許文献１に記載した理由と同様であるため、説明を省略する。
（集合領域を生成する処理）
次に、本数値解析手法の第二計算用データモデルを作成するための分割領域から集合領域を生成する処理について説明する。集合領域を生成する処理では、セルの集合和により集合領域を作成する。以下、集合領域は、ドメインと呼ぶこともある。ドメインが作成されることによって、解析領域はドメイン分割された状態となる。
図２で、集合領域を生成する処理では、解析領域内に自動生成されたコントロールボリューム（セル）を集合させ、新たに設定されたコントロールボリュームを、ドメインと定義する。ドメインは、コントロールボリュームであり、セルの集合和である。
図３に示される複数のドメインのうち、任意のドメインをドメインＡとする。ドメインＡ内に存在するセルの総数をＮＶＡとした場合、ドメインＡの体積ＶＡは式（２０）で表され、ドメインＡのコントロールポイントの座標ベクトル［ｒ］Ａは式（２１）で表される。以下の式（２０）、式（２１）により、セルの集合和を計算し、ドメインを設定する。
式（２０）、（２１）において、Ａ，Ｂ，・・・は、ドメインを表す添え字である。
ここで、ドメインを設定した後に、設定したドメインを集合させることによって、新たにドメインを設定するようにしてもよい。
図４に示されるように、複数のドメインの集合和が設定される。図４に示される例では、細い線によってドメインが示されている。
図５に示されるように、ドメインを複数集合させることによって、新たにドメインが設定される。新たに設定されたドメインは、太線によって示されている。
新たに設定されたドメインも、コントロールボリュームであり、ドメインの集合和である。新たに設定されたドメインについても、式（２０）、式（２１）により、新たに設定されたドメインの体積と、新たに設定されたドメインのコントロールポイントの座標ベクトルが計算される。図６に示されるように、コントロールポイントが示される。そして、新たに設定されたドメインは、ドメインと同様に取り扱われる。
セルの集合和によって作成されたドメインをドメイン１とし、ドメイン１の集合和によって作成されたドメインをドメイン２とし、ドメイン２の集合和によって作成されたドメインをドメイン３とする。解析領域に自動生成されたコントロールボリューム（セル）に基づいて、ドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型にドメインを設定できる。
要求される計算精度に合わせて、コントロールボリューム（セル）から、ドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型にドメインを設定できる。ここで、ドメイン１から、ドメイン２やドメイン３を経ることなく、最終ドメインを設定してもよい。コントロールボリューム（セル）から、ドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型にドメインを設定した場合や、ドメイン１から、ドメイン２を経ることなく、最終ドメインを設定した場合に、最初の解析領域の境界形状のセル分割精度が失われない。
従来の有限体積法や有限要素法では、初期の分割メッシュを集合させると、集合させたドメインの境界形状は複雑な多面となるため、計算を実行できない。具体的には、現状では、有限体積法では二十面体が限度であり、有限要素法では６面体を超える多面体の要素内補間関数を定義できない。このように、従来の方法においては、初期の分割メッシュを集合するという思想自体が存在しない。このことから、本実施形態は、分割領域から集合領域を形成することの動機づけも示唆もないところから着想され、従来の方法では不可能で、顕著な効果を奏する。
自動生成されたコントロールボリューム（セル）では、セル間でのマスバランス（質量保存）や運動量保存、エネルギ保存などの物理量の保存則を満足させながら数値解析ができる。従って、コントロールボリューム（セル）から、ドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型に設定されたドメインでも、ドメイン１から、ドメイン２やドメイン３を経ることなく、最終ドメインと階層構造型に設定された場合でも、ドメイン間での物理量の保存則は満足される。
解析領域内に自動生成されたコントロールボリューム（セル）を集合させることによって、新たにドメインを設定するための、コントロールボリューム（セル）の集合方法について説明する。
図７に示されるように、解析領域を直交格子状の領域に粗く分割し、その直交格子内にコントロールポイントの座標が含まれるセルを集合させる。
図８に示されるように、解析領域内に、ドメインのコントロールポイントを設定する。ドメインのコントロールポイントから、予め指定された半径の球体内にコントロールポイントの座標が含まれるセルを集合させる。半径を徐々に拡大し、解析領域内のセルを全て、ドメインのいずれかに集合させる。
また、ボクセル法で、解析領域を包含する領域にボクセルを生成し、そのボクセルをドメインとしてもよい。この場合、ボクセル内にコントロールポイントの座標が含まれるセルを集合させる。
ここでは、セルの集合方法として、三例について説明したが、この例に限られない。例えば、ここに示した例以外の集合方法が使用されてもよい。
解析領域内に自動生成されたコントロールボリューム（セル）を、コントロールボリューム（セル）の集合方法にしたがって集合させることによって新たにドメインを設定した場合には、最初の解析領域の境界形状のセル分割精度が失われない。このため、解析領域内に自動生成されたコントロールボリューム（セル）を、コントロールボリューム（セル）の集合方法にしたがって集合させることによって新たに設定されたドメイン間で、物理量の保存則を満足させながら数値解析ができる。
（集合領域での第二計算用データモデルを生成する処理）
さらに、集合領域での第二計算用データモデルを生成する処理について説明する。
図９は、発明の数値解析手法の集合領域における境界面特性量の一例を示す概念図である。図９は解析領域を分割する複数の分割領域と、複数の集合領域とを表す。図９において、実線で囲まれたドメインＡ及びドメインＢは、集合領域である。図９において、破線で囲まれた図形は分割領域である。例えば、セルＲ２０１〜セルＲ２０８は、分割領域である。ドメインＡは、セルＲ２０１〜セルＲ２０８を含む複数の分割領域を集合して得られる集合領域である。境界面ＥＡＢは、ドメインＡとドメインＢとの間において物理量の交換が行われる面であり、第二計算用データモデルにおける境界面に相当する。また、面積ＳＡＢは、境界面ＥＡＢの面積を示し、本実施形態における集合領域の境界面特性量の１つである。
［ｎ］ａ１〜［ｎ］ａ８の各々は、境界面ＥＡＢで互いに接するセル同士の境界面の特性を示す量である法線ベクトルである。
図１０のドメインＡ及びドメインＢは、図９のドメインＡ及びドメインＢである。コントロールポイントＡｃ及びコントロールポイントＢｃは、それぞれドメインＡ及びドメインＢの内部に配置されている。
ドメインＢに属し境界面ＥＡＢに接するセルｂの境界面の総数をＮＳＡＢとすると、ドメインＡとドメインＢとの間の境界面の面積ＳＡＢと、ドメインＡとドメインＢとの間の境界面の法線ベクトル［ｎ］ＡＢは、式（２２）、式（２３）で計算される。ここで、Ｓａｂは、ドメインＢに属し境界面ＥＡＢに接するセルｂの境界面の面積である。
図１１に示される例では、ドメインＡが解析領域のうちの外部空間との境界に接している場合を示す。この場合、図１２に示されるように、ドメインＢが解析領域の外部に存在すると仮定して、ドメインＢに接するセルａの境界面の集合和を計算することによって、式（２２）、式（２３）と同様に、ドメインＡとドメインＢとの間の境界面の面積ＳＡＢと、ドメインＡとドメインＢとの間の境界面の法線ベクトル［ｎ］ＡＢとを導出できる。
前述した数値解析手法では、解析領域内に自動生成されたコントロールボリューム（セル）を集合させることによって、ドメインを生成した。
さらに、数値解析手法では、セルａとセルｂとの間の境界面Ｓａｂにおける境界面特性量（境界面の面積Ｓａｂ、境界面の法線ベクトル［ｎ］ａｂ）を用いて、ドメインＡとドメインＢとの間の境界面ＳＡＢに関して、集合和を計算することによって、ドメインＡとドメインＢとの間の境界面ＳＡＢにおける境界面特性量（境界面の面積ＳＡＢ、境界面の法線ベクトル［ｎ］ＡＢ）を求めた。つまり、幾何学的形状を規定する量を使用しないセルによる連続体の数値解析手法を、ドメインに対しても全く同じ計算手法として適用した。
さらに、解析領域内に自動生成されたコントロールボリューム（セル）に基づいて、ドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型にドメインを設定することができる。したがって、連続体の数値解析手法を、ドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型に設定されたドメイン全てに対して適用できる。
コントロールボリューム（セル）から一度に最終ドメインを設定したときも、コントロールボリューム（セル）からドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型にドメインを設定したときも、最初の解析領域の境界形状のセル分割精度は失われない。これは、熱流体解析のように伝熱面積を重要視する数値解析では、非常に大きなメリットである。
解析領域内に自動的に生成されたコントロールボリューム（セル）の分割数が、数１０個〜数１０００個のオーダー以下の粗さの場合、その粗い分割数のセルを用いた数値解析は、解析結果に誤差を多く含むおそれがある。
一方、数１０００万個〜数億個以上のオーダーのセル分割による数値解析は、非常に高精度であるが、スーパーコンピュータを用いて行う計算レベルであり、大型のメモリ容量と、ＨＤＤ容量との計算資源が必要であり、長大な計算時間、解析結果処理のための多大な計算コストの増大を招く。
しかし、数１０００万個〜数億セル以上のオーダーのセル分割に関しては、セルの自動生成だけなら、比較的低容量のメモリを備えるコンピュータで、短時間で実行できる。
そこで、数１０００万個〜数億個以上のオーダーのセル分割を、解析領域内で行い、そのセル分割に基づいて、ドメイン１、ドメイン２、ドメイン３、・・・、最終ドメインと階層構造型にドメイン分割を行うこともできる。ドメイン分割数が、数千〜数万個のオーダーであれば、比較的低容量のメモリを備えるコンピュータで短時間に数値解析を実行できる。
この場合、セルを集合させたドメインの境界形状は非常に複雑な多面となるが、幾何学的形状を規定する量を使用しないセルによる連続体の数値解析手法を用いることによって、解析領域の境界形状のセル分割精度を維持した状態で、連続体の数値解析を、比較的低メモリのコンピュータで、解析結果に含まれる誤差を問題とならない程度に抑制しながら、短時間に実行できる。
次に、前述の第二計算用データモデルを用いて解析領域の集合領域における流速を求める物理量計算例について説明する。なお、ここでは、各集合コントロールポイントにおける流速を各集合領域における流速として求める。
第二計算用データモデルを用いる場合であっても、第一計算用データモデルを用いる場合と同様に、まずは、式（１）及び式（２）によって表されるナビエ・ストークスの式及び連続の式を重み付き残差積分法に基づいて集合領域の体積に対して積分する。
その結果、式（３）及び（４）が導出される。ただし、式（３）及び（４）において、第一計算用データモデルを用いる場合と異なり、Ｖが集合領域の体積を示し、∫ＶｄＶが集合領域の体積Ｖに関する積分を示し、Ｓが集合領域の面積を示し、∫ＳｄＳが集合領域の境界面Ｓに関する積分を示し、［ｎ］が集合領域の境界面Ｓの法線ベクトルを示し、ｎｉ（ｉ＝１，２，３）が集合領域の境界面Ｓの法線ベクトル［ｎ］の成分を示し、∂／∂ｎが集合領域の境界面Ｓの法線方向微分を示す。
第一計算用データモデルの場合と同様に、以下、説明を簡単化するために、流体の密度ρと粘性係数μを定数とする。ただし、第一計算用データモデルの場合と同様に、以下の定数化は、流体の物性値が、時間、空間、温度等によって変化する場合に対して拡張可能である。
集合コントロールポイントについて、図９と図１０に示すように、境界面ＥＡＢの面積ＳＡＢについて離散化し、代数方程式による近似式に変換すると、式（３）が下式（２４）、式（４）が下式（２５）のように示される。
ここで、添字ＡＢが付く、［ｎ］ＡＢ，［ｕ］ＡＢ，ｕｉＡＢ，ｎｉＡＢ，ＰＡＢ，（∂ｕｉ／∂ｎ）ＡＢは、集合領域Ａと集合領域Ｂとの間の境界面ＥＡＢ上における物理量であることを示す。また、ｎｉＡＢは［ｎ］ＡＢの成分である。また、ｍＡは、集合領域Ａと結合関係（境界面を挟む関係）にある全ての集合コントロールポイントの数である。
そして、式（２４），（２５）をＶＡ（集合領域Ａの体積）で割ると、式（２４）が下式（２６）のように示され、式（２５）が下式（２７）のように示される。
ここで、第一計算用データモデルの場合と同様に、ここで、下式（２８）の置換をおこなう。
すると、式（２６）が下式（２９）のように示され、式（２７）が下式（３０）のように示される。
式（２９），（３０）において、［ｕ］ＡＢ、ｕｉＡＢ、ＰＡＢ、（∂ｕｉ／∂ｎ）ＡＢは、集合コントロールポイントＡと集合コントロールポイントＢ上の物理量の重み付け平均（移流項については、風上性を考慮した重み付け平均）により近似的に求められ、集合コントロールポイントＡ，Ｂ間の距離及び向きと、その間に存在する集合領域の境界面ＥＡＢとの位置関係（上記比率α）と、集合領域の境界面ＥＡＢの法線ベクトルの向きに依存して決定される。ただし、［ｕ］ＡＢ、ｕｉＡＢ、ＰＡＢ、（∂ｕｉ／∂ｎ）ＡＢは、集合領域の境界面ＥＡＢの幾何学的な形状を規定する量には無関係な量である。
また、式（２８）で定義されるφＡＢも（面積／体積）という量であり、集合領域のコントロールボリュームの幾何学的形状を規定する量には無関係な量である。
つまり、このような式（２９）、（３０）は、集合領域の形状を規定するVertexとConnectivityを必要としない量のみを使用して物理量が算出可能な、重み付き残差積分法に基づく演算式である。
このため、物理量計算（ソルバ処理）に先立って前述の第二計算用データモデルを作成し、物理量計算において当該計算用データモデルと、式（２９），（３０）の離散化支配方程式とを用いることによって、物理量計算において集合領域の幾何学的形状を規定する量を全く使用せずに、流速の計算を行うことが可能である。
このように、物理量計算において集合領域の形状を規定するVertexとConnectivityとを全く使用せずに流速の計算が行えることから、第一計算用データモデルと同様に、第二計算用データモデルにおいても幾何学的形状を規定する量を持たせる必要がなくなる。
なお、実際に式（２９），（３０）を解くにあたり、［ｕ］ＡＢやＰＡＢ等の境界面ＥＡＢ上の物理量は、第一計算用データモデルでの計算式（１２）〜（１９）において、物理量Ψを分割領域のコントロールポイントの物理量から集合領域のコントロールポイントの物理量に置き換え、分割領域の境界面の面積と法線ベクトルを集合領域の境界面の面積（２２）式と法線ベクトル（２３）式に置き換え、分割領域のコントロールポイント間の距離（１８）式において分割領域のコントロールポイントの位置（座標）ベクトルを集合領域のコントロールポイントの位置（座標）ベクトル（２１）式に置き換えることで、第一計算用データモデルと同様に計算することができる。
次に、物理量計算にあたり、物理量の保存則が満足される条件について説明する。以下、まずは本数値解析における第一計算用データモデルを用いた計算において物理量の保存則が満足される条件を説明する。
第一計算用データモデルを用いた計算において、分割領域の質量保存則を表す連続の式（６）式におけるコントロールポイント間の分割領域の境界面の面積は、コントロールポイントａ側から見てもコントロールポイントｂ側から見ても等しいとすると、コントロールポイント間の質量流束（ρ［ｎ］・［ｕ］）・Ｓは、コントロールポイントａ側とコントロールポイントｂ側で正負が逆で絶対値が等しくなる。従って、（６）式を解析領域内の全分割領域に対して総和を取ると、分割領域間の質量流束は差し引きゼロとなりキャンセルされ、計算する解析領域全体に対して、流入する質量と流出する質量が等しいことを示す。
よって、計算する解析領域全体についての質量保存則を満足するためには、分割領域の２つのコントロールポイント間の境界面の面積が一致するという条件及び分割領域の境界面の法線ベクトルが一方のコントロールポイント側から見た場合と他方のコントロールポイント側から見た場合とで正負が逆で絶対値が一致するという条件が必要である。
また、質量保存則を満足するためには、下式（３１）に示される分割領域のコントロールボリュームの占める体積の総和が解析領域の全体積Ｖtotalと一致するという条件が必要である。なお、解析領域内の分割領域の総数はＮである。
なお、ここでは、質量保存の式に対して説明を行ったが、保存則は、連続体の運動量やエネルギに対しても成立しなければならない。これらの物理量に対しても、解析領域内の全分割領域に対して総和を取ることによって、保存側が満足されるためには、解析領域内の全分割領域のコントロールボリュームの体積の総和が解析領域の全体積と一致するという条件と、２つのコントロールポイント間の境界面の面積が一致する条件及び分割領域の境界面の法線ベクトルが一方のコントロールポイント側から見た場合と他方のコントロールポイント側から見た場合とで絶対値が一致する（正負逆符号）という条件とが必要であることが分かる。
また、保存則を満たすためには、図１３に示すように、分割領域のコントロールポイントａの占めるコントロールボリュームを考えた場合に、コントロールポイントａを通り、任意の向きの単位第一法線ベクトル［ｎ］ｐを持つ無限に広い投影面Ｐを考えたときに下式（３２）が成り立つという条件が必要である。単位法線ベクトルは、単位長さの法線ベクトルである。
なお、図１３及び式（３２）において、Ｓｉが境界面Ｅｉの面積、［ｎ］ｉが境界面Ｅｉの単位第一法線ベクトル、ｍがコントロールボリュームの面の総数を示す。
式（３２）は、コントロールボリュームを構成する多面体が、閉包空間を構成することを示す。この式（３２）は、コントロールボリュームを構成する多面体の一部が凹んでいる場合であっても成立する。
また、多面体の１つの面を微小面ｄＳとし、ｍを∞とする極限を取ると、下式（３３）となり、図１４に示すような閉包曲面体についても成り立つことが分かる。
式（３２）が成り立つという条件は、ガウスの発散定理や、式（１４）に示す一般化されたグリーンの定理が成り立つために必要な条件である。
そして、一般化されたグリーンの定理は、連続体の離散化のための基本となる定理である。したがって、グリーンの定理にしたがって体積分を面積分に変形し離散化させる場合において、保存則を満足させるためには式（３２）が成り立つという条件は必須となる。
このように、前述の第一計算用データモデル及び物理量計算を用いて数値解析を行う際に、物理量の保存則が満足されるためには、以下の３つの条件（以下「第一計算用条件」という。）が必要となる。
（ａ１）全コントロールポイントのコントロールボリュームの体積（全分割領域の体積）の総和が解析領域の体積と一致する。
（ｂ１）２つのコントロールポイント間の境界面の面積が一致する及び第一法線ベクトルが一方のコントロールポイント側（境界面を挟む一方の分割領域）から見た場合と他方のコントロールポイント側（境界面を挟む他方の分割領域）から見た場合とで絶対値が一致する。
（ｃ１）コントロールポイントを通り（分割領域を通り）、任意の向きの単位第一法線ベクトル［ｎ］ｐを持つ無限に広い投影面Ｐを考えたときに式（３２）が成り立つ。
第二計算用データモデルを用いた計算において、集合領域の質量保存則を表す連続の式（２５）式における集合コントロールポイント間の集合領域の境界面の面積は、集合コントロールポイントＡ側から見てもコントロールポイントＢ側から見ても等しいとすると、集合コントロールポイント間の質量流束は、集合コントロールポイントＡ側と集合コントロールポイントＢ側で正負が逆で絶対値が等しくなる。従って、（２５）式を解析領域内の全集合領域に対して総和を取ると、集合領域間の質量流束は差し引きゼロとなりキャンセルされ、計算する解析領域全体に対して、流入する質量と流出する質量が等しいことを示す。
よって、計算する解析領域全体についての質量保存則を満足するためには、集合領域の２つのコントロールポイント間の境界面の面積が一致するという条件及び集合領域の境界面の法線ベクトルが一方のコントロールポイント側から見た場合と他方のコントロールポイント側から見た場合とで正負が逆で絶対値が一致するという条件が必要である。
また、質量保存則を満足するためには、下式（３４）に示される集合領域のコントロールボリュームの占める体積の総和が解析領域の全体積Ｖtotalと一致するという条件が必要である。解析領域内の集合領域の総数をＮＤとする。
なお、ここでは、質量保存の式に対して説明を行ったが、保存則は、連続体の運動量やエネルギに対しても成立しなければならない。これらの物理量に対しても、全集合コントロールポイントに対して足し加えることによって、保存側が満足されるためには、全集合コントロールポイントの集合コントロールボリュームの占める体積が解析領域の全体積と一致するという条件と、２つの集合コントロールポイント間の境界面の面積が一致する条件及び第二法線ベクトルが一方の集合コントロールポイント側から見た場合と他方の集合コントロールポイント側から見た場合とで絶対値が一致する（正負逆符号）という条件とが必要であることが分かる。
また、保存則を満たすためには、図１５に示すように、集合コントロールポイントＡの占める集合コントロールボリュームを考えた場合に、集合コントロールポイントＡを通り、任意の向きの単位第二法線ベクトル［Ｎ］Ｐを持つ無限に広い投影面Ｐｄを考えたときに下式（３５）が成り立つという条件が必要である。単位第二法線ベクトルは、単位長さの第二法線ベクトルである。
なお、式（３５）において、Ｑｉがドメインのひとつの境界面Ｅｄｉの面積、［Ｎ］ｉが境界面Ｅｄｉの単位第二法線ベクトル、Ｍが集合コントロールボリュームの面の総数を示す。なお、添え字のｉは１〜Ｍの整数である。
式（３５）は、集合コントロールボリュームを構成する多面体が、閉包空間を構成することを示す。この式（３５）は、集合コントロールボリュームを構成する多面体の一部が凹んでいる場合であっても成立する。
また、多面体の１つの面を微小面ｄＱとし、Ｍを∞とする極限を取ると、下式（３６）となり、図１６に示すような閉包曲面体についても成り立つことが分かる。
式（３５）が成り立つという条件は、ガウスの発散定理や、式（１４）に示す一般化されたグリーンの定理が成り立つために必要な条件である。
そして、一般化されたグリーンの定理は、連続体の離散化のための基本となる定理である。したがって、グリーンの定理にしたがって体積分を面積分に変形し離散化させる場合において、保存則を満足させるためには式（３５）が成り立つという条件は必須となる。
このように、前述の第二計算用データモデル及び物理量計算を用いて数値解析を行う際に、物理量の保存則が満足されるためには、以下の３つの条件（以下「第二計算用条件」という。）が必要となる。
（ａ２）全集合領域の体積の総和が解析領域の体積と一致する。
（ｂ２）２つの集合コントロールポイント間の境界面の面積が一致する及び第二法線ベクトルが一方の集合コントロールポイント側（境界面を挟む一方の集合領域）から見た場合と他方の集合コントロールポイント側（境界面を挟む他方の集合領域）から見た場合とで絶対値が一致する。
（ｃ２）集合コントロールポイントを通り（集合領域を通り）、任意の向きの単位集合第二法線ベクトル［Ｎ］Ｐを持つ無限に広い投影面Ｐｄを考えたときに上式（３５）が成り立つ。
つまり、保存則を満足させる場合には、これらの条件を満足するように第一及び第二計算用データモデルを作成する必要がある。ただし、前述のように本数値解析手法においては、計算用データモデルの作成にあたり、分割領域のセル形状を任意に変形することができ、その分割領域セルを集合和させることにより集合領域を作成することから、容易に上記３つの条件を満足するように第一及び第二計算用データモデルを作成することができる。
このように、本数値解析手法においては、第一計算用データモデルだけでなく、第一計算用データモデルよりも解析領域の分割数が少ない第二計算用データモデルも保存則を満足する。そのため、本数値解析手法は、解析領域の分割数が少なくても解析精度が悪くならない、という特徴を有する。
なお、前述の説明においては、ナビエ・ストークスの式及び連続の式から重み付き残差積分法に基づいて導出した離散化支配方程式を用いる物理量の計算例について説明したが、本数値解析手法において用いられる離散化支配方程式はこれに限られるものではない。
つまり、種々の方程式（質量保存の方程式、運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式、及び波動方程式等）から重み付き残差積分法に基づいて導出されると共に、幾何学的形状を規定する量を必要としない量のみを使用して物理量を算出可能な離散化支配方程式であれば本数値解析手法に用いることができる。
そして、このような離散化支配方程式の特性によって、従来の有限要素法や有限体積法のようにいわゆるメッシュを必要としない、メッシュレスでの計算が可能となる。また、たとえ、プリ処理において、セルの幾何学的形状を規定する量を利用するとしても、従来の有限要素法、有限体積法、ボクセル法のようなメッシュに対する制約がないため、計算用データモデルの作成に伴う作業負荷を軽減できる。具体的には、従来の有限要素法や有限体積法等に基づくソフトによる幾何学的形状を規定する量に基づいて、各分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量を求めてもよい。
本実施形態では、質量保存の方程式、運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式、及び波動方程式から、重み付き残差積分法に基づいて、幾何学的形状を規定する量を必要としない量のみを使用する離散化支配方程式が導出可能である。そのため、本数値解析手法において他の支配方程式を用いることができる。これについては、特許文献１に記載した理由と同様であるため、説明を省略する。
なお、特許文献１では、解析領域が分割領域によって分割される場合について離散化支配方程式が導出可能であることを説明するが、解析領域が集合領域によって分割される場合についても、同様に離散化支配方程式が導出可能である。
（適用例）
以下の説明においては、本実施形態に係る物理量計算方法を含む数値解析方法と、本実施形態に係る物理量計算プログラムを含む数値解析プログラムと、本実施形態に係る物理量計算装置を含む数値解析装置との具体的な適用例について説明する。
また、以下の具体的な適用例においては、自動車のキャビン空間における空気の流速を数値解析によって求める場合について説明する。
図１７は、本実施形態の数値解析装置Ａのハードウェア構成を概略的に示すブロック図である。
この図に示すように、本実施形態の数値解析装置Ａは、パーソナルコンピュータやワークステーション等のコンピュータによって構成され、ＣＰＵ１、記憶装置２、ＤＶＤドライブ３、入力装置４、出力装置５、及び通信装置６を備えている。
ＣＰＵ１は、記憶装置２、ＤＶＤドライブ３、入力装置４、出力装置５、及び通信装置６と電気的に接続されており、これらの各種装置から入力される信号を処理すると共に、処理結果を出力する。なお、ＣＰＵ１は、演算部１ｏｐの具体例である。
記憶装置２は、メモリ等の内部記憶装置及びハードディスクドライブ等の外部記憶装置によって構成されており、ＣＰＵ１から入力される情報を記憶すると共にＣＰＵ１から入力される指令に基づいて記憶した情報を出力する。
そして、本実施形態において記憶装置２は、プログラム記憶部２ａとデータ記憶部２ｂとを備えている。
プログラム記憶部２ａは、数値解析プログラムＰを記憶している。この数値解析プログラムＰは、所定のＯＳ(Operating System)において実行されるアプリケーションプログラムであり、コンピュータから構成される本実施形態の数値解析装置Ａを、数値解析を行うように機能させる。そして、演算部１ｏｐが数値解析プログラムＰを実行することによって、本実施形態の数値解析装置Ａの各機能が実現される。
そして、図１７に示すように、数値解析プログラムＰは、プリ処理プログラムＰ１と、ソルバ処理プログラムＰ２と、ポスト処理プログラムＰ３とを有している。
プリ処理プログラムＰ１は、ソルバ処理を実行するための前処理（プリ処理）を本実施形態の数値解析装置Ａに実行させ、本実施形態の数値解析装置Ａを第一計算用データモデル作成部として機能させることによって第一計算用データモデルを作成させる。また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第二計算用データモデル作成部として機能させることによって第二計算用データモデルを作成させる。また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに、ソルバ処理を実行するにあたり必要となる条件の設定を実行させ、さらには上記計算用データモデルや設定された条件を纏めたソルバ入力データファイルＦの作成を実行させる。なお、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第一計算用データモデル作成部として機能させ、第一計算用データモデルが作成された後、本実施形態の数値解析装置Ａを第二計算用データモデル作成部として機能させる。
プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第一計算用データモデル作成部及び第二計算用データモデル作成部として機能させる場合に、まず本実施形態の数値解析装置Ａに対して、自動車のキャビン空間を含む３次元形状データを取得させ、この取得させた３次元形状データに含まれる自動車のキャビン空間を示す解析領域の作成を実行させる。
なお、後に詳説するが、本実施形態においては、ソルバ処理において、前述の本実施形態を用いた数値解析手法にて説明した離散化支配方程式を用いる。前述の本実施形態を用いた数値解析手法にて説明した離散化支配方程式とは、具体的には、幾何学的形状を規定する量を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化支配方程式である。このため、第一計算用データモデル及び第二計算用データモデルの作成にあたり、保存則を満たす条件の下、分割領域の形状、分割領域に基づく集合領域の形状及び解析領域の形状を任意に変更できる。よって、３次元形状データに含まれる自動車のキャビン空間の修正あるいは変更作業は簡易的なもので充分となる。そこで、本実施形態においてプリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに対して、取得させた３次元形状データに含まれる自動車のキャビン空間に存在する穴や隙間に、微小な閉曲面を覆いかぶせるラッピング処理によって修繕する処理を実行させる。
その後、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに対して、複数の分割領域に分割する処理で説明したように分割領域を形成し、ラッピング処理等により修繕されたキャビン空間の全領域を含む解析領域の作成を実行させる。続いて、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに対して作成された分割領域のうちキャビン空間から食み出した領域をカットすることによって、キャビン空間を示す解析領域の作成を実行させる。ここでも、ソルバ処理において前述の離散化支配方程式を用いることから、解析領域のうちキャビン空間から食み出した領域を容易にカットすることができる。
これにより、ボクセル法のように、外部空間との境界が階段状になることがなく、また、ボクセル法のカットセル法のような外部空間の境界付近の解析領域の形成に対して、経験や試行錯誤を要する非常に膨大な手作業を伴う特別な修正または処理を必要としない。そのため、本実施形態では、ボクセル法で問題となる外部空間との境界の処理に関わる問題がない。
なお、本実施形態においては、後述のようにキャビン空間とカットした領域との隙間に新たな任意形状の分割領域を充填することによって、直交格子形状のみによらない分割領域で解析領域が構成されるようにし、さらには解析領域に分割領域を重なることなく充填させている。
次に、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第一計算用データモデル作成部として機能させる場合に、本実施形態の数値解析装置Ａに対して、作成させたキャビン空間を示す解析領域に含まれる分割領域の各々の内部に対して１つのコントロールポイントを仮想的に配置する処理を実行させ、コントロールポイントの配置情報、及び各コントロールポイントが占めるコントロールボリュームの体積データを記憶させる。
また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第一計算用データモデル作成部として機能させる場合に、本実施形態の数値解析装置Ａに対して、上記分割領域同士の境界面である境界面の面積及び第一法線ベクトルの算出を実行させ、これらの境界面の面積及び第一法線ベクトルを記憶させる。
また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第一計算用データモデル作成部として機能させる場合に、コントロールボリューム又はコントロールポイントの結合情報（link）を作成させ、このlinkを記憶させる。
そして、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに対して、上記各コントロールポイントが占めるコントロールボリュームの体積と、境界面の面積及び第一法線ベクトルと、分割領域の配置情報を表すコントロールポイントの配置情報と、linkとを纏めさせて第一計算用データモデルを作成させる。配置情報が示す配置は、例えば座標を用いて示されてもよい。
プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第二計算用データモデル作成部として機能させる場合に、本実施形態の数値解析装置Ａに対して、作成された第一計算用データモデルを用いて、図７あるいは図８に示す方法でコントロールボリューム（セル）を集合させ、キャビン空間を示す解析領域内に、集合領域を生成する。
次に、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第二計算用データモデル作成部として機能させる場合に、本実施形態の数値解析装置Ａに対して、作成させたキャビン空間を示す解析領域に含まれる集合領域の各々の内部に対して１つの集合コントロールポイントを仮想的に配置する処理を実行させ、集合コントロールポイントの配置情報、及び各集合コントロールポイントが配置されたドメインの体積データを記憶させる。
プリ処理プログラムＰ１は、上記集合コントロールポイントを仮想的に配置する処理において、本実施形態の数値解析装置Ａにたいして、集合コントロールポイント算出処理を実行させることで、各集合コントロールポイントを配置する箇所を示す情報を取得する。集合コントロールポイント算出処理は、本実施形態の数値解析装置Ａが第一計算用データモデルが有するコントロールボリュームの体積とコントロールポイントの配置情報とを取得し、式（２０）及び式（２１）の計算を行うことで、集合コントロールポイントの配置箇所を算出する処理である。
また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第二計算用データモデル作成部として機能させる場合に、本実施形態の数値解析装置Ａに対して、上記集合領域同士の境界面の面積及び第二法線ベクトルの算出（以下「集合領域境界面特性量算出処理」という。）を実行させ、これらの境界面の面積及び第二法線ベクトルを記憶させる。
集合領域境界面特性量算出処理は、本実施形態の数値解析装置Ａが第一計算用データモデルが有する境界面の面積と第一法線ベクトルとの情報を取得し、式（２２）及び式（２３）の計算を実行することで、集合領域同士の境界面の面積及び第二法線ベクトルを算出する処理である。
また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａを第二計算用データモデル作成部として機能させる場合に、ドメイン又は集合コントロールポイントの結合情報（link）を作成させ、このlinkを記憶させる。
そして、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに対して、上記各集合コントロールポイントが配置されたドメインの体積と、境界面の面積及び第二法線ベクトルと、集合領域の配置情報を表す集合コントロールポイントの配置情報と、linkとを纏めさせて第二計算用データモデルを作成させる。配置情報が示す配置は、例えば座標を用いて示されてもよい。
また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに、前述のソルバ処理を実行するにあたり必要となる条件の設定を行わせる場合には、物性値の設定、境界条件の設定、初期条件の設定、計算条件の設定を行わせる。
ここで、物性値とは、キャビン空間における空気の密度、粘性係数等である。
境界条件とは、コントロールポイント間の物理量の交換の法則を規定する条件であり、、本実施形態においては前述した式（１０）で示されるナビエ・ストークスの式に基づく離散化支配方程式、及び式（１１）で示される連続の式に基づく離散化支配方程式である。
また、境界条件には、キャビン空間と外部空間との境界面に臨む集合領域を示す情報が含まれる。
初期条件とは、ソルバ処理を実行する際の最初の物理量を示すものであり、各分割領域及び集合領域の流速の初期値である。
計算条件とは、ソルバ処理における計算の条件であり、例えば反復回数や収束基準である。
また、プリ処理プログラムＰ１は、本実施形態の数値解析装置Ａに、ＧＵＩ（Graphical User Interface）を形成させる。より詳細には、プリ処理プログラムＰ１は、出力装置５が備えるディスプレイ５ａに対してグラフィックを表示させると共に、入力装置４が備えるキーボード４ａやマウス４ｂによって操作が可能な状態とさせる。
ソルバ処理プログラムＰ２（物理量計算プログラム）は、本実施形態の数値解析装置Ａにソルバ処理を実行させるものであり、本実施形態の数値解析装置Ａを物理量計算装置として機能させる。
そして、ソルバ処理プログラムＰ２は、本実施形態の数値解析装置Ａを物理量計算部として機能させる場合に、ソルバ入力データファイルＦを用いて、解析領域における物理量を数値計算させる。
そして、ソルバ処理プログラムＰ２は、本実施形態の数値解析装置Ａを物理量計算部として機能させる場合に、本実施形態の数値解析装置Ａに対して、ソルバ入力データファイルＦに含まれるナビエ・ストークスの式及び連続の式の離散化係数行列の作成を実行させると共に、マトリックス形成用のデータテーブルの作成を実行させる。
また、ソルバ処理プログラムＰ２は、本実施形態の数値解析装置Ａを物理量計算部として機能させる場合に、本実施形態の数値解析装置Ａに対して、前述した式（２９）で示されるナビエ・ストークスの式に基づく離散化支配方程式、及び前述した式（３０）で示される連続の式に基づく離散化支配方程式から、下式（３７）で示すマトリックス計算用の大規模疎行列方程式の組み上げを実行させる。
なお、ソルバ処理プログラムＰ２は、集合領域から作成される第二計算用データモデルを入力とする数値計算だけでなく、分割領域から作成される第一計算用データモデルを入力とする数値計算を実行することもできる。この場合は、本実施形態の数値解析装置Ａに対して、前述した式（１０）で示されるナビエ・ストークスの式に基づく離散化支配方程式、及び前述した式（１１）で示される連続の式に基づく離散化支配方程式から、下式（３７）で示すマトリックス計算用の大規模疎行列方程式の組み上げを実行させる。
なお、式（３７）において［Ａ］が大規模疎行列を示し、［Ｂ］が境界条件ベクトルを示し、［Ｘ］が流速の解を示す。
また、ソルバ処理プログラムＰ２は、上記離散化支配方程式に非圧縮性等の付帯条件が存在する場合には、本実施形態の数値解析装置Ａに対して、この付帯条件の行列方程式への組み上げを実行させる。
そして、ソルバ処理プログラムＰ２は、本実施形態の数値解析装置Ａに対して、ＣＧ法（共役勾配法）等による行列方程式の解の計算、当該解の下式（３８）を用いた解のアップデート、収束条件の判定を実行させ、最終的な計算結果を取得させる。
ポスト処理プログラムＰ３は、本実施形態の数値解析装置Ａにポスト処理を実行させるものであり、本実施形態の数値解析装置Ａに対して、ソルバ処理において取得された計算結果に基づく処理を実行させる。
より詳細には、ポスト処理プログラムＰ３は、本実施形態の数値解析装置Ａに対して、計算結果の可視化処理、抽出処理を実行させる。
ここで、可視化処理とは、例えば、断面コンター表示、ベクトル表示、等値面表示、アニメーション表示を出力装置５に出力させる処理である。また、抽出処理とは、作業者が指定する領域の定量値を抽出して数値やグラフとして出力装置５に出力させる、あるいは作業者が指定する領域の定量値を抽出してファイル化したものの出力を実行させる処理である。
また、ポスト処理プログラムＰ３は、本実施形態の数値解析装置Ａに対して、自動レポート作成、計算残差の表示及び分析を実行させる。
データ記憶部２ｂは、図１７に示すように、第一計算用データモデルＭ１、第二計算用データモデルＭ２、境界条件を示す境界条件データＤ１、計算条件を示す計算条件データＤ２、物性値を示す物性値データＤ３、及び初期条件を示す初期条件データＤ４を有するソルバ入力データファイルＦと、３次元形状データＤ５と、計算結果データＤ６等を記憶するものである。また、データ記憶部２ｂは、ＣＰＵ１の処理過程において生成される中間データを一時的に記憶する。
ＤＶＤドライブ３は、ＤＶＤメディアＸを取り込み可能に構成されており、ＣＰＵ１から入力される指令に基づいて、ＤＶＤメディアＸに記憶されるデータを出力するものである。そして、本実施形態においては、ＤＶＤメディアＸに数値解析プログラムＰが記憶されており、ＤＶＤドライブ３は、ＣＰＵ１から入力される指令に基づいて、ＤＶＤメディアＸに記憶される数値解析プログラムＰを出力する。
入力装置４は、本実施形態の数値解析装置Ａと作業者とのマンマシンインターフェイスであり、ポインティングデバイスであるキーボード４ａやマウス４ｂを備えている。
出力装置５は、ＣＰＵ１から入力される信号を可視化して出力するものであり、ディスプレイ５ａ及びプリンタ５ｂを備えている。
通信装置６は、本実施形態の数値解析装置ＡとＣＡＤ装置Ｃ等の外部装置との間においてデータの受け渡しを行い、社内ＬＡＮ（Local Area Network）等のネットワークＢに対して電気的に接続されている。
次に、このように構成された本実施形態の数値解析装置Ａを用いた数値解析方法（本実施形態の数値解析方法）について、図１８〜図２０のフローチャートを参照して説明する。
図１８のフローチャートに示すように、本実施形態の数値解析方法は、プリ処理（ステップＳ１）と、ソルバ処理（ステップＳ２）と、ポスト処理（ステップＳ３）とから構成されている。
なお、本実施形態の数値解析方法を行うより前に、ＣＰＵ１は、ＤＶＤドライブ３に取り込まれたＤＶＤメディアＸに記憶された数値解析プログラムＰをＤＶＤメディアＸから取り出し、記憶装置２のプログラム記憶部２ａに記憶させる。
そして、ＣＰＵ１は、入力装置４から数値解析の開始を指示する信号が入力されると、記憶装置２に記憶された数値解析プログラムＰに基づいて数値解析を実行する。より詳細には、ＣＰＵ１は、プログラム記憶部２ａに記憶されたプリ処理プログラムＰ１に基づいてプリ処理（ステップＳ１）を実行し、プログラム記憶部２ａに記憶されたソルバ処理プログラムＰ２に基づいてソルバ処理（ステップＳ２）を実行し、プログラム記憶部２ａに記憶されたポスト処理プログラムＰ３に基づいてポスト処理（ステップＳ３）を実行する。なお、このようにＣＰＵ１がプリ処理プログラムＰ１に基づくプリ処理（ステップＳ１）を実行することによって、本実施形態の数値解析装置Ａが計算用データモデル作成部として機能される。また、ＣＰＵ１がソルバ処理プログラムＰ２に基づくソルバ処理（ステップＳ２）を実行することによって、本実施形態の数値解析装置Ａが物理量計算部として機能される。
図１９は、プリ処理（ステップＳ１）を示すフローチャートである。この図に示すように、プリ処理（ステップＳ１）が開始されると、ＣＰＵ１は、通信装置６に、ネットワークＢを介してＣＡＤ装置Ｃから自動車のキャビン空間を含む３次元形状データＤ５を取得させる（ステップＳ１ａ）。ＣＰＵ１は、取得した３次元形状データＤ５を記憶装置２のデータ記憶部２ｂに記憶させる。
続いて、ＣＰＵ１は、ステップＳ１ａで取得した３次元形状データＤ５に基づいて、第一計算用データモデルの作成を実行する（ステップＳ１ｂ）。
具体的には、ＣＰＵ１は、キャビン空間を示す解析領域に含まれる各分割領域内に１つのコントロールポイントを仮想的に配置する。ここでは、ＣＰＵ１は、分割領域の重心を算出し、各々の重心に対して１つのコントロールポイントを仮想的に配置する。そして、ＣＰＵ１は、コントロールポイントの配置情報、各コントロールポイントが占めるコントロールボリュームの体積（コントロールポイントが配置される分割領域の体積）を算出し、記憶装置２のデータ記憶部２ｂに一時的に記憶させる。
また、ＣＰＵ１は、分割領域同士の境界面である境界面の面積及び第一法線ベクトルを算出し、これらの境界面の面積及び第一法線ベクトルを記憶装置２のデータ記憶部２ｂに一時的に記憶させる。
また、ＣＰＵ１は、分割領域間のlinkを作成し、この分割領域間のlinkを記憶装置２のデータ記憶部２ｂに一時的に記憶させる。
そして、ＣＰＵ１は、データ記憶部２ｂに記憶された、コントロールポイントの配置情報と、各コントロールポイントが占めるコントロールボリュームの体積と、境界面の面積及び法線ベクトルと、linkとをデータベース化することによって第一計算用データモデルを作成し、作成した第一計算用データモデルを記憶装置２のデータ記憶部２ｂ内に記憶させる。
また、本実施形態では、ステップＳ１ｂにおいて、先に分割領域を形成し、その後コントロールポイントを配置し、各コントロールポイントに対して、自らが配置された分割領域の体積を割り当てる構成を採用している。
しかしながら、本実施形態においては、先にコントロールポイントを解析領域に配置し、各コントロールポイントに対して後から体積を割り当てることも可能である。
具体的には、例えば、異なるコントロールポイントにぶつかるまでの半径や、結合関係にある（linkで関連付けられた）コントロールポイントまでの距離に基づいて、各コントロールポイントに対して重み付けを行う。
ここでコントロールポイントｉの重みをｗｉ、基準体積をＶ＋とし、コントロールポイントｉに割り当てられる体積Ｖｉを下式（３９）とする。
さらに、各コントロールポイントの体積Ｖｉの総和は、解析領域の体積Ｖtotalと等しいため、下式（４０）が成り立つ。
この結果、基準体積Ｖ＋は下式（４１）で求めることができる。
したがって、各コントロールポイントに割り当てる体積は、式（４０），（４１）から求めることができる。
このような方法を用いれば、プリ処理において、幾何学的形状を規定する量を用いることなく、第一計算用データモデルに持たせる分割領域の体積を求めることができる。
また、当該第一計算用データモデルの作成（ステップＳ１ｂ）において、ＣＰＵ１は、ＧＵＩを形成し、ＧＵＩから指令（例えば分割領域の密度を示す指令や分割領域の形状を示す指令）が入力された場合には、当該指令を反映させた処理を実行する。したがって、作業者は、ＧＵＩを操作することによって、コントロールポイントの配置や分割領域の形状を任意に調節することが可能とされている。
ただし、ＣＰＵ１は、数値解析プログラムに記憶された保存則を満足するための３つの条件に照らし合わせ、ＧＵＩから入力される指令が、当該条件から外れる場合には、その旨をディスプレイ５ａに表示させる。
続いて、ＣＰＵ１は、ステップＳ１ｂで作成した第一計算用データモデルに基づいて、第二計算用データモデルの作成を実行する（ステップＳ１ｃ）。
具体的には、ＣＰＵ１は、作成された第一計算用データモデルを用いて、図７あるいは図８に示す方法でコントロールボリューム（セル）を集合させ、キャビン空間を示す解析領域内に、集合領域を生成する。次に、ＣＰＵ１は、キャビン空間を示す解析領域に含まれる各集合領域内に１つの集合コントロールポイントを仮想的に配置する。ここでは、ＣＰＵ１は、第一計算用データモデルに基づいて、式（２０）及び（２１）によって、集合コントロールポイントを算出し、各々の集合領域に対して１つの集合コントロールポイントを仮想的に配置する。そして、ＣＰＵ１は、集合コントロールポイントの配置情報、各集合コントロールポイントが占める集合コントロールボリュームの体積（集合コントロールポイントが配置される集合領域の体積）を算出し、記憶装置２のデータ記憶部２ｂに一時的に記憶させる。
また、ＣＰＵ１は、集合領域同士の境界面である境界面の面積及び第二法線ベクトルを算出し、これらの境界面の面積及び第二法線ベクトルを記憶装置２のデータ記憶部２ｂに一時的に記憶させる。
また、ＣＰＵ１は、集合領域のlinkを作成し、この集合領域のlinkを記憶装置２のデータ記憶部２ｂに一時的に記憶させる。
そして、ＣＰＵ１は、データ記憶部２ｂに記憶された、集合コントロールポイントの配置情報と、各集合コントロールポイントを有する集合領域の体積と、境界面の面積及び第二法線ベクトルと、linkとをデータベース化することによって第二計算用データモデルを作成し、作成した第二計算用データモデルを記憶装置２のデータ記憶部２ｂ内に記憶させる。
また、当該第二計算用データモデルの作成（ステップＳ１ｃ）において、ＣＰＵ１は、ＧＵＩを形成し、ＧＵＩから指令（例えば分割領域の密度を示す指令や分割領域の形状を示す指令）が入力された場合には、当該指令を反映させた処理を実行する。したがって、作業者は、ＧＵＩを操作することによって、集合コントロールポイントの配置や集合領域の形状を任意に調節することが可能とされている。なお、集合領域は、あくまでもコントロールボリューム(セル)の集合なので、集合領域に対してコントロールボリューム(セル)を無視して形状変更はできない。
ただし、ＣＰＵ１は、数値解析プログラムに記憶された保存則を満足するための３つの条件に照らし合わせ、ＧＵＩから入力される指令が、当該条件から外れる場合には、その旨をディスプレイ５ａに表示させる。
続いて、ＣＰＵ１は、物性値データの設定（ステップＳ１ｄ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に物性値の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される物性値を示す信号を物性値データＤ３としてデータ記憶部２ｂに一時的に記憶させることで物性値の設定を行う。なお、ここで言う物性値とは、キャビン空間における流体である空気の特性値であり、空気の密度、粘性係数等である。
続いて、ＣＰＵ１は、境界条件データの設定（ステップＳ１ｅ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に境界条件の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される境界条件を示す信号を境界条件データＤ１としてデータ記憶部２ｂに一時的に記憶させることで境界条件データの設定を行う。なお、ここで言う境界条件とは、キャビン空間の物理現象を支配する離散化支配方程式や、キャビン空間と外部空間との境界面に臨む集合コントロールポイントの特定情報、及びキャビン空間と外部空間との間における熱の伝熱条件等を示す。
なお、本実施形態の数値解析方法は、キャビン空間における流速を数値解析により求めることを目的とするため、上記離散化支配方程式として、前述のナビエ・ストークスの式に基づく離散化支配方程式（２９）及び連続の式に基づく離散化支配方程式（３０）が用いられる。
なお、これらの離散化支配方程式は、例えば、数値解析プログラムＰに予め記憶された複数の離散化支配方程式をディスプレイ５ａ上に表示された複数の離散化支配方程式から作業者がキーボード４ａやマウス４ｂを用いることによって選択される。
続いて、ＣＰＵ１は、初期条件データの設定（ステップＳ１ｆ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に初期条件の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される初期条件を示す信号を初期条件データＤ４としてデータ記憶部２ｂに一時的に記憶させることで初期条件データの設定を行う。なお、ここ言う初期条件とは、各コントロールポイント（各分割領域）における初期流速、及び各集合コントロールポイント（各集合領域）における初期流速である。
続いて、ＣＰＵ１は、計算条件データの設定（ステップＳ１ｇ）を行う。具体的には、ＣＰＵ１は、ＧＵＩを用いて、ディスプレイ５ａ上に計算条件の入力画面を表示し、キーボード４ａあるいはマウス４ｂから入力される計算条件を示す信号を計算条件データＤ２としてデータ記憶部２ｂに一時的に記憶させることで計算条件データの設定を行う。なお、ここで言う計算条件とは、ソルバ処理（ステップＳ２）における計算の条件であり、例えば、反復回数や収束基準を示す。
続いて、ＣＰＵ１は、ソルバ入力データファイルＦの作成（ステップＳ１ｈ）を行う。
具体的には、ＣＰＵ１は、ステップＳ１ｂにて作成された第一計算用データモデルＭ１と、ステップＳ１ｃにて作成された第二計算用データモデルＭ２と、ステップＳ１ｄで設定された物性値データＤ３と、ステップＳ１ｅで設定された境界条件データＤ１と、ステップＳ１ｆで設定された初期条件データＤ４と、ステップＳ１ｇで設定された計算条件データＤ２とをソルバ入力データファイルＦに格納することによってソルバ入力データファイルＦを作成する。なお、このソルバ入力データファイルＦは、データ記憶部２ｂに記憶される。
以上のようなプリ処理（ステップＳ１）が完了すると、ＣＰＵ１は、ソルバ処理プログラムＰ２に基づいて、図１８のフローチャートに示すソルバ処理（ステップＳ２）を実行する。
図２０に示すように、ソルバ処理（ステップＳ２）が開始されると、ＣＰＵ１は、プリ処理（ステップＳ１）で作成されたソルバ入力データファイルＦを取得する（ステップＳ２ａ）。なお、本実施形態に示す数値解析方法のように、単一の装置（本実施形態の数値解析装置Ａ）によってプリ処理及びソルバ処理を実行する場合には、既にデータ記憶部２ｂにソルバ入力データファイルＦが記憶されているため、ステップＳ２ａを省略することができる。ただし、プリ処理（ステップＳ１）とソルバ処理（ステップＳ２）とが異なる装置において実行される場合には、ネットワークやリムーバルディスクによって搬送されるソルバ入力データファイルＦを取得する必要があるため、本ステップＳ２ａを行う必要がある。
続いて、ＣＰＵ１は、ソルバ入力データの整合性を判定する（ステップＳ２ｂ）。なお、ソルバ入力データとは、ソルバ入力データファイルＦに格納されたデータを示し、第一計算用データモデルＭ１、第二計算用データモデルＭ２、境界条件データＤ１、計算条件データＤ２、物性値データＤ３及び初期条件データＤ４である。
具体的には、ＣＰＵ１は、ソルバ処理において物理量計算を実行可能なソルバ入力データがソルバ入力データファイルＦに全て格納されているかを分析することによってソルバ入力データの整合性の判定を行う。
そして、ＣＰＵ１は、ソルバ入力データが不整合であると判定した場合には、ディスプレイ５ａにエラーを表示させ（ステップＳ２ｂ＋）、さらには不整合である部分のデータを入力するための画面を表示させる。その後、ＣＰＵ１は、ＧＵＩから入力される信号に基づいてソルバ入力データの調整を行い（ステップＳ２ｃ）、再度ステップＳ２ａを実行する。
一方、ＣＰＵ１は、ステップＳ２ｂにおいてソルバ入力データの整合性があると判定した場合には、初期計算処理（ステップＳ２ｅ）を実行する。
具体的には、ＣＰＵ１は、境界条件データＤ１に記憶された離散化支配方程式から離散化係数行列を作成し、さらにマトリクス計算用のデータテーブルの作成を行うことによって初期計算処理を行う。なお、境界条件データＤ１に記憶された離散化支配方程式とは、例えば、ナビエ・ストークスの式に基づく離散化支配方程式（２９）及び連続の式に基づく離散化支配方程式（３０）である。
なお、分割領域から作成される第一計算用データモデルをソルバ入力データとする数値計算を実行する場合は、ＣＰＵ１は、境界条件データＤ１に記憶された離散化支配方程式から離散化係数行列を作成し、さらにマトリクス計算用のデータテーブルの作成を行うことによって初期計算処理を行う。なお、境界条件データＤ１に記憶された離散化支配方程式とは、例えば、ナビエ・ストークスの式に基づく離散化支配方程式（１０）及び連続の式に基づく離散化支配方程式（１１）である。
続いて、ＣＰＵ１は、大規模疎行列方程式の組み上げ（ステップＳ２ｆ）を行う。具体的には、ＣＰＵ１は、ナビエ・ストークスの式に基づく離散化支配方程式（２９）及び連続の式に基づく離散化支配方程式（３０）から、前述の式（３７）で示すマトリックス計算用の大規模疎行列方程式の組み上げを行う。
なお、分割領域から作成される第一計算用データモデルをソルバ入力データとする数値計算を実行する場合は、ＣＰＵ１は、ナビエ・ストークスの式に基づく離散化支配方程式（１０）及び連続の式に基づく離散化支配方程式（１１）から、前述の式（３７）で示すマトリックス計算用の大規模疎行列方程式の組み上げを行う。
続いて、ＣＰＵ１は、離散化支配方程式に、非圧縮性や接触等の付帯条件が存在するかの判定を行う。この付帯条件は、例えば、境界条件データとしてソルバ入力データファイルＦに格納されている。
そして、ＣＰＵ１は、離散化支配方程式に付帯条件が存在すると判定した場合には当該付帯条件の大規模行列方程式の組み込み（ステップＳ２ｈ）を実行した後に大規模行列方程式の計算（ステップＳ２ｉ）を実行する。一方、ＣＰＵ１は、離散化支配方程式に付帯条件が存在しないと判定した場合には付帯条件の大規模行列方程式の組み込み（ステップＳ２ｈ）を実行することなく大規模行列方程式の計算（ステップＳ２ｉ）を実行する。
そしてＣＰＵ１は、大規模行列方程式を例えば、ＣＧ法（共役勾配法）によって解き、前述の式（３８）を用いて解のアップデート（ステップＳ２ｊ）を行う。
続いて、ＣＰＵ１は、式（３８）の残差が収束条件に達したか否かの判定（ステップＳ２ｋ）を行う。具体的には、ＣＰＵ１は、式（３８）の残差を計算し、計算条件データＤ２に含まれる収束条件と比較し、これによって式（３８）の残差が収束条件に達したか否かの判定を行う。
そして、残差が収束条件に達していないと判定した場合には、ＣＰＵ１は、物性値のアップデートを行った後、再度ステップＳ２ｇを実行する。つまり、ＣＰＵ１は、式（３８）の残差が収束条件に達するまで、物性値のアップデートを行いながらステップＳ２ｆ〜Ｓ２ｋを繰り返し行う。
一方、残差が収束条件に達したと判定した場合には、ＣＰＵ１は、計算結果の取得を行う（ステップＳ２ｌ）。具体的には、ＣＰＵ１は、直前のステップＳ２ｉにおいて計算された物理量の解を計算結果データとしてデータ記憶部２ｂに記憶させることによって計算結果の取得を行う。
このようなソルバ処理（ステップＳ２）によって、キャビン空間における空気の流速が求められる。なお、このようなソルバ処理（ステップＳ２）は、本実施形態の物理量計算方法に相当する。
以上のようなソルバ処理（ステップＳ２）が完了すると、ＣＰＵ１は、ポスト処理プログラムＰ３に基づいてポスト処理（ステップＳ３）を実行する。
具体的には、例えばＣＰＵ１は、ＧＵＩから入力される指令に基づいて、計算結果データから、例えば断面コンタデータ、ベクトルデータ、等値面データ、アニメーションデータを生成し、当該データを、出力装置５に可視化させる。
また、ＣＰＵ１は、ＧＵＩから入力される指令に基づいて、キャビン空間の一部における定量値（計算結果）を抽出して数値やグラフとし、この数値やグラフを出力装置５に可視化させ、さらには数値やグラフをファイルとして纏めて出力する。また、ＣＰＵ１は、ＧＵＩから入力される指令に基づいて、例えば計算結果データから自動レポート作成、計算残差の表示及び分析を行ってその結果を出力する。
さらに、本実施形態の数値解析装置Ａ、数値解析方法及び数値解析プログラムの利用者が、ポスト処理の結果に応じて３次元形状データを変更して再度、図１８のフローチャートの処理を繰り返す場合にも、現実的な時間内で物理量の計算が可能となる。
即ち、利用者は、ポスト処理により得られた計算結果データを評価することにより、解析対象である３次元形状データにより所望の結果が得られたと判断する場合には、シミュレーションを終了してよい。また利用者は、ポスト処理により得られた計算結果データを評価することにより、解析対象である３次元形状データにより所望の結果が得られていないと判断する場合には、３次元形状データを修正してから再度シミュレーションを実行してよい。
上記の動作において、シミュレーションが所望の結果を示す場合、その解析対象であった３次元形状データが表現する物理的実体(閉鎖された空間を構成する自動車キャビン、コックピット、住宅、若しくは電気機器や産業機器の内部等、又はガラスや鉄鋼等の製造装置等)の設計が満足できるものと判断し、当該物理的実体を製造・生産してよい。またシミュレーションが所望の結果を示さない場合、その解析対象で、あった３次元形状データが表現する物理的実体の設計が満足できないものと判断し、当該物理的実体の設計を変更し、この設計変更後の３次元形状データに基づいて再度シミュレーションを実行することになる。
以上のような本実施形態の数値解析装置Ａ、数値解析方法及び数値解析プログラムによれば、プリ処理にて集合コントロールボリュームの体積と境界面の面積及び第二法線ベクトルとを有する第二計算用データモデルＭ２が作成され、ソルバ処理にて第二計算用データモデルＭ２に含まれる集合コントロールボリュームの体積と境界面の面積及び第二法線ベクトルとを用いて各集合コントロールボリュームにおける物理量が計算される。
このように、本実施形態の数値解析方法を用いることで、物理量を計算することができる。そのため、本実施形態の数値解析方法は、物理現象を数値的に解析する物理現象解析方法である。
なお、本実施形態の数値解析装置Ａ、数値解析方法及び数値解析プログラムにおいては、解析領域に分割領域を重なることなく充填させている。このため、前述した保存側を満足するための６つの条件（ａ１）〜（ｃ１）及び（ａ２）〜（ｃ２）が満たされることとなり、保存則を満足して流速を計算することができる。
このように構成された本実施形態の数値解析装置Ａは、幾何学的形状を規定する量を有さない計算用データモデルであって、保存則を満たす計算用データモデルを作成するため、解析領域の分割数を少なくしてソルバ処理における計算負荷を軽くすることと、解析領域の分割数を少なくしても解析精度を悪くしないこととを両立することが可能である。
また、本実施形態の数値解析装置Ａ、数値解析方法及び数値解析プログラムによれば、前述のように、プリ処理における第一及び第二計算用データモデルの作業負担が大幅に減少し、またソルバ処理における計算負荷を軽減させることが可能となる。
したがって、解析領域が移動境界を含み解析領域の形状が時系列的に変化する場合であっても、本実施形態によれば、図２１のフローチャートに示すように、解析領域が形状変化するたびにプリ処理とソルバ処理とを繰り返し行うことによって、現実的な時間内で物理量の計算が可能となる。 さらに、本実施形態の数値解析装置Ａ、数値解析方法及び数値解析プログラムの利用者が、ポスト処理の結果に応じて３次元形状データを変更して再度、図２１のフローチャートの処理を繰り返す場合にも、現実的な時間内で物理量の計算が可能となる。
即ち、利用者は、ポスト処理により得られた計算結果データを評価することにより、解析対象である３次元形状データにより所望の結果が得られたと判断する場合には、シミュレーションを終了してよい。また利用者は、ポスト処理により得られた計算結果データを評価することにより、解析対象である３次元形状データにより所望の結果が得られていないと判断する場合には、３次元形状データを修正してから再度シミュレーションを実行してよい。
上記の動作において、シミュレーションが所望の結果を示す場合、その解析対象であった３次元形状データが表現する物理的実体(閉鎖された空間を構成する自動車キャビン、コックピット、住宅、若しくは電気機器や産業機器の内部等、又はガラスや鉄鋼等の製造装置等)の設計が満足できるものと判断し、当該物理的実体を製造・生産してよい。またシミュレーションが所望の結果を示さない場合、その解析対象で、あった３次元形状データが表現する物理的実体の設計が満足できないものと判断し、当該物理的実体の設計を変更し、この設計変更後の３次元形状データに基づいて再度シミュレーションを実行することになる。
移動境界は、解析領域の中で、対象とする物体が移動することによって変化する物体の境界をいう。解析領域が移動境界を含み解析領域の形状が時系列的に変化する場合としては、例えば、キャビンにおいて、人がいない状態から、人がキャビン内に入るまでの現象を再現する場合が挙げられる。その他、例えば、加熱炉の中を加熱対象物が移動する現象を再現する場合が挙げられる。
以上、添付図面を参照しながら本実施形態の好適な実施形態について説明したが、本実施形態は、上記実施形態に限定されないことは言うまでもない。前述した実施形態において示した各構成部材の諸形状や組み合わせ等は一例であって、本実施形態の主旨から逸脱しない範囲において設計要求等に基づき種々変更可能である。
上記実施形態においては、運動量保存の方程式の変形例であるナビエ・ストークスの式及び連続の式から導出した離散化支配方程式を用いて空気の流速を数値解析によって求める構成について説明した。
しかしながら、本実施形態はこれに限定されるものではなく、質量保存の方程式、運動量保存の方程式、角運動量保存の方程式、エネルギ保存の方程式、移流拡散方程式及び波動方程式の少なくともいずれかから導出した離散化支配方程式を用いて物理量を数値解析によって求めることが可能である。
また、上記実施形態においては、本実施形態の境界面特性量として、境界面の面積と境界面の法線ベクトルとを用いる構成について説明した。
しかしながら、本実施形態はこれに限定されるものではなく、境界面特性量として他の量（例えば境界面の周長）を用いることもできる。
また、上記実施形態においては、保存則を満足するために前述の６つの条件を満たすように計算用データモデルを作成する構成について説明した。
しかしながら、本実施形態はこれに限定されるものではなく、保存則を満足させる必要がない場合には、計算用データモデルを必ずしも前述の６つの条件を満たすように作成する必要はない。
また、上記実施形態においては、分割領域の体積を、当該分割領域の内部に配置されるコントロールポイントが占めるコントロールボリュームの体積として捉えた構成について説明した。
しかしながら、本実施形態はこれに限定されるものではなく、分割領域の内部に対してコントロールポイントを配置する必要は必ずしもない。このような場合には、コントロールポイントが占めるコントロールボリュームの体積を分割領域の体積に置き換えることによって数値解析を行うことができる。
また、上記実施形態においては、解析領域をまずは複数のコントロールボリューム（セル）に分割して第一計算用データモデルを作成し、次に、コントロールボリューム（セル）を集合させた集合領域（ドメイン）によって解析領域を分割して第二計算用データモデルを作成した。そして、そのソルバ処理においては、第二計算用データモデルを用いて計算を行った。しかし、第二計算用データモデルは、セルを集合させたドメインに基づく計算用データモデルであればよく、セルを集合させたドメインをさらに集合させた新たなドメインでの計算用データモデルであってもよい。
図２２〜図３８は、本実施形態における第一計算用データモデル（分割領域）及び第二計算用データモデル（集合領域）を使用して実施した熱流体シミュレーション結果の一例を示す図である。
本実施形態では、流体解析の基礎方程式として、式（１）で示すナビエ・ストークスの式と、式（２）で示す連続の式を用いた数値計算方法を説明したが、前述の通り、本実施形態では、ナビエ・ストークスの式だけでなく、移流拡散方程式に対しても重み付き残差積分法に基づいて、幾何学的形状を規定する量を必要としない量のみを使用する離散化支配方程式が導出可能であることを説明した。図２２〜図３８に示す熱流体シミュレーションの実行では、熱流体解析の基礎方程式として、式（１）で示すナビエ・ストークスの式と、式（２）で示す連続の式に加え、熱の移流拡散方程式を基礎方程式として使用した。熱の移流拡散方程式に対する本実施形態での重み付き残差積分法に基づく幾何学的形状を規定する量を必要としない量のみを使用する離散化支配方程式の導出については、特許文献１に詳細が記載されているためここでは説明を省略する。
熱流体シミュレーションでは、解析領域の一例は自動車のキャビンであり、境界条件の一例は夏季で、冷房空調条件である。具体的には、車外参照温度は３５度、車外熱伝達率は４０Ｗ／ｍ２Ｋ、乗車人員数は４名、空調吹出風速は５ｍ／ｓ、空調吹出温度は８℃である。なお、必要に応じて、境界条件として、エンジンルームの温度、トランクルームの温度、床下フロアーの温度、ダッシュボードの内側の温度、天井の温度、その他の少なくとも一つが含まれる場合等も解析可能である。
数値解析装置Ａは、解析領域（自動車のキャビン）を、前述したように、頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としないセルに分割する。
解析領域（自動車のキャビン）を、約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の一例を、図２２〜図２４に示す。図２２は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、図２３は、その垂直断面上の７点のサンプリング点での温度値（以降、温度値の単位はＫ（ケルビン）で表示する。）を示す。また、図２４は、解析領域であるキャビン空間内の運転席中央での垂直断面上の流速コンター図であり、流速の向きを矢印の向き、流速の大きさを矢印の大きさで示す。図２２〜図２４に示す結果は、定常状態の計算結果を得るまでに、インテル社製ＣＰＵ Ｘｅｏｎ（２．６ＧＨｚ）を搭載したＰＣを用いて、約３０時間の計算時間を要した。
本実施形態では、数値解析装置Ａは、解析領域（自動車のキャビン）に対して約４５０万個のセルを自動生成する。数値解析装置Ａは、前述したセルから集合領域を生成する処理によって、図２５〜図３２に示す約４５０万個のセルから２７個の集合領域（ドメイン）を生成する場合について説明する。２７個の集合領域の各々を生成した結果の一例を、図２５〜図３０に示す。図２５〜図３０の各々において、ＤＣＰは集合領域（ドメイン）のコントロールポイントであり、ＣＣＰはセルのコントロールポイントの一つである。図２５〜図３０は、順番にドメイン１からドメイン６までを示す。また、図２５〜図３０には、図示されていないが、数値解析装置Ａは、各集合領域について、外表面を取得する。数値解析装置Ａは、取得した外表面と、その外表面に接している部材との間の境界条件を取得する。境界条件は、外表面が接する部材の材料によって異なる場合がある。
図３１と図３２は、前述の２７個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の一例を示す。
図３１は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その垂直断面上の７点のサンプリング点での温度値を示す。図３１の７点のサンプリング点での温度値は、それぞれ上段の温度値と下段の括弧内の温度値に分けて示しているが、上段は前述の２７個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値であり、下段の括弧内の温度値は図２３に示す約４５０万個の分割領域（セル）に分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値を示す。上段と下段の温度値を比較すると、数度の温度値の差はあるものの平均的には良く一致している。
図３２は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その同じ垂直断面上に流速ベクトルを重ねて示した図である。キャビン空間内に空調吹出流れによる循環流が形成されていることが分かるが、図２４に示す約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果と比較すると、図３２ではキャビン空間内の流れの詳細な分布までは計算できていないことが分かる。
図３１と図３２に示す２７個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションは、定常状態の計算結果を得るまでに、インテル社製ＣＰＵ Ｘｅｏｎ（２．６ＧＨｚ）を搭載したＰＣを用いて、1秒以下の計算時間を要した。約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションに約３０時間の計算時間を要したことと比較すると数値計算は非常に高速である。
図３３と図３４は、７９２個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の一例を示す。
図３３は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その垂直断面上の７点のサンプリング点での温度値を示す。図３１と同様に、７点のサンプリング点での温度値は、それぞれ上段の温度値と下段の括弧内の温度値に分けて示しているが、上段は７９２個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値であり、下段の括弧内の温度値は図３５に示す約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値を示す。上段と下段の温度値を比較すると、数度の温度値の差はあるものの平均的には良く一致している。
図３４は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その同じ垂直断面上に流速ベクトルを重ねて示した図である。キャビン空間内に空調吹出流れによる循環流が形成されていることが分かるが、図３２と比較すると、図３４ではキャビン空間内の流れの詳細度が向上していることが分かる。
図３３と図３４に示す７９２個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションは、定常状態の計算結果を得るまでに、インテル社製ＣＰＵ Ｘｅｏｎ（２．６ＧＨｚ）を搭載したＰＣを用いて、約２０秒の計算時間を要した。約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションに約３０時間の計算時間を要したことと比較すると数値計算は非常に高速である。
図３５と図３６は、１６０５５個の集合領域を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の一例を示す。
図３５は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その垂直断面上の７点のサンプリング点での温度値を示す。図３１と同様に、７点のサンプリング点での温度値は、それぞれ上段の温度値と下段の括弧内の温度値に分けて示しているが、上段は１６０５５個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値であり、下段の括弧内の温度値は図３５に示す約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値を示す。上段と下段の温度値を比較すると、温度値の差は１度程度であり非常に良く一致している。
図３６は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その同じ垂直断面上に流速ベクトルを重ねて示した図である。キャビン空間内に空調吹出流れによる循環流が形成されていることが分かるが、図３２、図３４と比較すると、図３６ではキャビン空間内の流れの詳細度が大きく向上しており、キャビン空間内の局所的な流れや渦も計算できていることが分かる。
図３５と図３６に示す１６０５５個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションは、定常状態の計算結果を得るまでに、インテル社製ＣＰＵ Ｘｅｏｎ（２．６ＧＨｚ）を搭載したＰＣを用いて、約３分の計算時間を要した。約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションに約３０時間の計算時間を要したことと比較すると数値計算は非常に高速である。
図３７と図３８は、６６２５７個の集合領域を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の一例を示す。
図３７は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その垂直断面上の７点のサンプリング点での温度値を示す。図３１と同様に、７点のサンプリング点での温度値は、それぞれ上段の温度値と下段の括弧内の温度値に分けて示しているが、上段は６６２５７個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値であり、下段の括弧内の温度値は図２３に示す約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の温度値を示す。上段と下段の温度値を比較すると、温度値の差はゼロであり非常に良く一致している。
図３８は、解析領域であるキャビン空間内の運転席中央での垂直断面上の温度コンター図であり、その同じ垂直断面上に流速ベクトルを重ねて示した図である。キャビン空間内に空調吹出流れによる循環流が形成されていることが分かるが、図３２、図３４と比較すると、図３８ではキャビン空間内の流れの詳細度が大きく向上しており、キャビン空間内の局所的な流れや渦も計算できており、図２４に示す約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果と比較しても精度上差異のない結果である。
図３７と図３８に示す６６２５７個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションは、定常状態の計算結果を得るまでに、インテル社製ＣＰＵ Ｘｅｏｎ（２．６ＧＨｚ）を搭載したＰＣを用いて、約１２分の計算時間を要した。約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションに約３０時間の計算時間を要したことと比較すると数値計算は非常に高速である。
前述の通り、集合領域（ドメイン）の個数が、２７個、７９２個、１６０５５個、６６２５７個の集合領域（ドメイン）を含む第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の一例を示した。また、約４５０万個のセルに分割した第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果と比較して、第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーション結果の精度を示した。さらに、第二計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションは、第一計算用データモデルを使用して実施した３Ｄ熱流体シミュレーションの計算速度と比較して、非常に高速であることを示した。解析する目的や必要とされる精度に応じて、集合領域（ドメイン）の数を適切に選択することによって、特許文献１に相当するセルのみを使用した計算に比べ非常に高速に解析結果を得ることができることを示した。
また、上記実施形態においては、数値解析プログラムＰがＤＶＤメディアＸに記憶されて搬送可能な構成について説明した。
しかしながら、本実施形態はこれに限定されるものではなく、数値解析プログラムＰを他のリムーバブルメディアに記憶させて搬送可能とする構成を採用することもできる。
また、プリ処理プログラムＰ１とソルバ処理プログラムＰ２とを別々のリムーバブルメディアに記憶させて搬送可能とすることもできる。また、数値解析プログラムＰは、ネットワークを介して伝達することも可能である。
なお、本実施形態は、例えば、自動車ボディの形状、エアコンなどのHeating Ventilation and Air Conditioning（HVAC）での消費エネルギ、ガラス、人の存在、外部日射エネルギ、湿度、車速等をシミュレーションモデルに反映した自動車のキャビンの温熱解析に用いられてもよい。
また、本実施形態は、自動車のキャビンの温熱解析以外にも、自動車のエンジンの燃焼解析や、自動車の燃焼ガスの排気効率の解析や、自動車のエンジンルームの温熱解析等に用いられてもよい。
また、本実施形態は、自動車以外の分野における温熱解析に用いられてもよい。例えば、航空機、船舶、宇宙船、宇宙ステーションのキャビン、コックピット等の内部空間の温熱解析に用いられてもよい。また、例えば、住宅、ビル、アトリウム等の内部空間の温熱解析に用いられてもよい。また、例えば、電気機器や産業機器の内部の温熱解析に用いられてもよい。また、例えば、ガラス、鉄鋼等の製造設備の装置自体の温熱解析に用いられてもよいし、それら製造設備の装置周辺の温熱解析に用いられてもよい。
なお、第一計算用データモデルは、分割領域の計算用データモデルの一例である。また、第二計算用データモデルは、集合領域の計算用データモデルの一例である。また、分割領域の境界面特性量は、分割領域特性量の一例である。また、集合領域の境界面特性量は、集合領域特性量の一例である。また、第一法線ベクトルは、分割領域の境界面の法線ベクトルの一例である。また、第二法線ベクトルは、集合領域の境界面の法線ベクトルの一例である。なお、重み付き残差積分法に基づいて導出された離散化支配方程式と集合領域の体積と境界面特性量とによって算出された物理量と、空気の流速とは、解析結果である物理量の一例である。なお、数値解析方法はシミュレーション方法の一例である。
本出願を詳細にまた特定の実施態様を参照して説明したが、本発明の精神と範囲を逸脱することなく様々な変更や修正を加えることができることは当業者にとって明らかである。
上述した実施の形態に関し、さらに以下の付記を開示する。
（付記１）
コンピュータによって物理現象での物理量を数値的に解析するシミュレーション方法であって、コンピュータが、外部装置から解析領域の３次元形状データを取得し、
前記解析領域を複数の分割領域に分割し、
前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成し、
前記分割領域を複数集合させることによって、要求される数の集合領域を生成し、
前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成し、
前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算し、
前記解析結果である物理量の可視化データを生成して出力装置に表示し、
当該シミュレーション方法の利用者が、前記出力装置に表示された内容に応じて前記解析領域の前記形状を変更して再度、前記解析領域の前記形状変更後の３次元形状データから前記出力装置への表示までを繰り返すシミュレーション方法。
（付記２）
コンピュータによって物理現象での物理量を数値的に解析するシミュレーション方法であって、
コンピュータが、外部装置から解析領域の３次元形状データを取得して前記解析領域を複数の分割領域に分割し、
前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成し、
前記分割領域を複数集合させることによって、要求される数の集合領域を生成し、
前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成し、
前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算し、
前記解析領域が移動境界を含む場合には、前記解析領域の形状が時系列的に変化したかどうかを判断し、変化した場合には、前記前記分割領域での計算用データモデルの生成、前記集合領域での計算用データモデルの生成、前記集合領域での解析結果である物理量の計算を繰り返し、
変化しない場合には、前記解析結果である物理量の可視化データを生成して出力装置に表示する、シミュレーション方法。
（付記３）
コンピュータによって物理現象での物理量を数値的に解析するシミュレーション方法であって、
コンピュータが、外部装置から解析領域の３次元形状データを取得して前記解析領域を複数の分割領域に分割し、
前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成し、
前記分割領域を複数集合させることによって、要求される数の集合領域を生成し、
前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成し、
前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算し、前記解析結果である物理量の可視化データを生成して出力装置に表示する、シミュレーション方法。
（付記４）
前記分割領域での計算用データモデルの生成において、
全分割領域の体積の総和が解析領域の体積と一致するという条件と、
前記分割領域同士の境界面の面積が一致するという条件及び法線ベクトルが当該境界面を挟む一方の分割領域から見た場合と他方の分割領域から見た場合とで絶対値が一致するという条件と、
前記分割領域を通る無限に広い投影面Ｐの任意の向きの単位法線ベクトルが［ｎ］ｐ、境界面の面積がＳｉ、境界面の単位法線ベクトルが［ｎ］ｉ、分割領域の面の総数がｍ、前記［］で囲まれた文字がベクトルを示す太字であるときに下式（１）が成り立つという条件と、
が満足されるように前記分割領域を形成する、付記１から３のいずれか一項に記載の方法。
（付記５）
前記集合領域での計算用データモデルの生成において、
全集合領域の体積の総和が解析領域の体積と一致するという条件と、
隣り合う前記集合領域同士の境界面の面積が一致するという条件及び法線ベクトルが当該境界面を挟む一方の集合領域から見た場合と他方の集合領域から見た場合とで絶対値が一致するという条件と、
前記集合領域を通る無限に広い投影面Ｐの任意の向きの単位法線ベクトルが［N］P、境界面の面積がＱｉ、境界面の単位法線ベクトルが［N］ｉ、集合領域の面の総数がM、前記［］で囲まれた文字がベクトルを示す太字であるときに下式（２）が成り立つという条件と、
が満足されるように前記集合領域を形成する付記１から４のいずれか一項に記載の方法。
（付記６）
前記分割領域特性量は、隣り合う前記分割領域同士の境界面の特性を示す境界面特性量と、隣り合う前記分割領域同士の結合情報と、隣り合う前記分割領域同士の距離とからなり、
前記集合領域特性量は、隣り合う前記集合領域同士の境界面の特性を示す境界面特性量と、隣り合う前記集合領域同士の結合情報と、隣り合う前記集合領域同士の距離とからなる、
付記１から５のいずれか一項に記載の方法。
（付記７）
前記隣り合う前記分割領域同士の境界面の特性を示す境界面特性量は、前記隣り合う前記分割領域同士の境界面の面積と前記境界面との法線ベクトルであり、
前記隣り合う前記集合領域同士の境界面の特性を示す前記境界面特性量は、前記隣り合う前記集合領域同士の境界面の面積と前記境界面の法線ベクトルとである、
付記６に記載の方法。
（付記８）
前記分割領域での計算用データモデルの生成において、前記分割領域の体積と隣り合う前記分割領域同士の特性を示す前記分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）から得る、付記１から７のいずれか一項に記載の方法。
（付記９）
コンピュータに、外部装置から解析領域の３次元形状データを取得させて前記解析領域を複数の分割領域に分割させ、
前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成させ、
前記分割領域を複数集合させることによって、要求される数の集合領域を生成させ、
前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成させ、
前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算させる、ことを特徴とする物理量計算プログラム。
（付記１０）
物理現象での物理量を数値的に解析する物理量計算装置であって、
データを表示する出力装置と、
外部装置との間においてデータの受け渡しを行う通信装置と、
前記通信装置を介して前記外部装置から解析領域の３次元形状データを取得し、前記解析領域を複数の分割領域に分割し、前記分割領域を複数集合させることによって、要求される数の集合領域を生成する演算部と、
前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法によって導出された離散化された分割領域での支配方程式と、前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法によって導出された離散化された集合領域での支配方程式とを記憶する記憶部とを備え、
前記演算部は、前記記憶部に記憶された前記分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成し、前記記憶部に記憶された前記集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成し、前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算し、
前記解析結果である物理量の可視化データを生成して前記出力装置に表示することを特徴とする物理量計算装置。
Ａ……数値解析装置（物理量計算装置）
Ｐ……数値解析プログラム
Ｐ１……プリ処理プログラム
Ｐ２……ソルバ処理プログラム（物理量計算プログラム）
Ｐ３……ポスト処理プログラム
Ｍ１……第一計算用データモデル
Ｍ２……第二計算用データモデル
１……ＣＰＵ
２……記憶装置
２ａ……プログラム記憶部
２ｂ……データ記憶部
３……ＤＶＤドライブ
４……入力装置
４ａ……キーボード
４ｂ……マウス
５……出力装置
５ａ……ディスプレイ
５ｂ……プリンタ

Claims (8)

  1. コンピュータによって物理現象での物理量を数値的に解析するシミュレーション方法であって、
    コンピュータが、解析領域を複数の分割領域に分割し、
    前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成し、
    前記分割領域を複数集合させることによって、要求される数の集合領域を生成し、
    前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成し、
    前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算する、ことを特徴とするシミュレーション方法。
  2. 前記分割領域での計算用データモデルの生成において、
    全分割領域の体積の総和が解析領域の体積と一致するという条件と、
    前記分割領域同士の境界面の面積が一致するという条件及び法線ベクトルが当該境界面を挟む一方の分割領域から見た場合と他方の分割領域から見た場合とで絶対値が一致するという条件と、
    前記分割領域を通る無限に広い投影面Ｐの任意の向きの単位法線ベクトルが［ｎ］ｐ、境界面の面積がＳｉ、境界面の単位法線ベクトルが［ｎ］ｉ、分割領域の面の総数がｍ、前記［］で囲まれた文字がベクトルを示す太字であるときに下式（１）が成り立つという条件と、
    が満足されるように前記分割領域を形成する、請求項１に記載の方法。
  3. 前記集合領域での計算用データモデルの生成において、
    全集合領域の体積の総和が解析領域の体積と一致するという条件と、
    隣り合う前記集合領域同士の境界面の面積が一致するという条件及び法線ベクトルが当該境界面を挟む一方の集合領域から見た場合と他方の集合領域から見た場合とで絶対値が一致するという条件と、
    前記集合領域を通る無限に広い投影面Ｐの任意の向きの単位法線ベクトルが［N］P、境界面の面積がＱｉ、境界面の単位法線ベクトルが［N］ｉ、集合領域の面の総数がM、前記［］で囲まれた文字がベクトルを示す太字であるときに下式（２）が成り立つという条件と、
    が満足されるように前記集合領域を形成する請求項１又は２に記載の方法。
  4. 前記分割領域特性量は、隣り合う前記分割領域同士の境界面の特性を示す境界面特性量と、隣り合う前記分割領域同士の結合情報と、隣り合う前記分割領域同士の距離とからなり、
    前記集合領域特性量は、隣り合う前記集合領域同士の境界面の特性を示す境界面特性量と、隣り合う前記集合領域同士の結合情報と、隣り合う前記集合領域同士の距離とからなる、
    請求項１から３のいずれか一項に記載の方法。
  5. 前記隣り合う前記分割領域同士の境界面の特性を示す境界面特性量は、前記隣り合う前記分割領域同士の境界面の面積と前記境界面との法線ベクトルであり、
    前記隣り合う前記集合領域同士の境界面の特性を示す前記境界面特性量は、前記隣り合う前記集合領域同士の境界面の面積と前記境界面の法線ベクトルとである、
    請求項４に記載の方法。
  6. 前記分割領域での計算用データモデルの生成において、前記分割領域の体積と隣り合う前記分割領域同士の特性を示す前記分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）から得る、請求項１から５のいずれか一項に記載の方法。
  7. コンピュータに、解析領域を複数の分割領域に分割させ、
    前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成させ、
    前記分割領域を複数集合させることによって、要求される数の集合領域を生成させ、
    前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法に基づいて導出された離散化された集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成させ、
    前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算させる、ことを特徴とする物理量計算プログラム。
  8. 物理現象での物理量を数値的に解析する物理量計算装置であって、
    解析領域を複数の分割領域に分割し、前記分割領域を複数集合させることによって、要求される数の集合領域を生成する演算部と、
    前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法によって導出された離散化された分割領域での支配方程式と、前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量のみを使用すると共に重み付き残差積分法によって導出された離散化された集合領域での支配方程式とを記憶する記憶部とを備え、
    前記演算部は、前記記憶部に記憶された前記分割領域での支配方程式に基づき、各前記分割領域の体積と隣り合う前記分割領域同士の特性を示す分割領域特性量とを前記分割領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記分割領域での計算用データモデルを生成し、前記記憶部に記憶された前記集合領域での支配方程式に基づき、各前記集合領域の体積と隣り合う前記集合領域同士の特性を示す集合領域特性量とを前記集合領域の頂点の座標（Vertex）及び該頂点の連結情報（Connectivity）を必要としない量として有する前記集合領域での計算用データモデルを生成し、前記解析領域での物性値と、前記集合領域での計算用データモデルとに基づいて、前記集合領域での解析結果である物理量を計算する、ことを特徴とする物理量計算装置。
JP2019506743A 2018-09-12 2018-11-28 シミュレーション方法、物理量計算プログラム及び物理量計算装置 Active JP6504333B1 (ja)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
JP2018170766 2018-09-12
JP2018170766 2018-09-12
PCT/JP2018/043836 WO2020054087A1 (ja) 2018-09-12 2018-11-28 シミュレーション方法、物理量計算プログラム及び物理量計算装置

Publications (2)

Publication Number Publication Date
JP6504333B1 JP6504333B1 (ja) 2019-04-24
JPWO2020054087A1 true JPWO2020054087A1 (ja) 2020-10-22

Family

ID=66324129

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2019506743A Active JP6504333B1 (ja) 2018-09-12 2018-11-28 シミュレーション方法、物理量計算プログラム及び物理量計算装置

Country Status (2)

Country Link
US (1) US20200082035A1 (ja)
JP (1) JP6504333B1 (ja)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE102017127276A1 (de) * 2017-08-30 2019-02-28 Taiwan Semiconductor Manufacturing Co., Ltd. Standardzellen und abwandlungen davon innerhalb einer standardzellenbibliothek

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102804187B (zh) * 2009-06-25 2016-06-29 旭硝子株式会社 物理量计算方法、数值解析方法、物理量计算装置及数值解析装置
EP2610772A4 (en) * 2010-08-24 2015-04-22 Asahi Glass Co Ltd CALCULATION DATA GENERATION APPARATUS, CALCULATION DATA GENERATION METHOD, AND CALCULATION DATA GENERATION PROGRAM

Also Published As

Publication number Publication date
JP6504333B1 (ja) 2019-04-24
US20200082035A1 (en) 2020-03-12

Similar Documents

Publication Publication Date Title
JP6516081B1 (ja) シミュレーション方法、ｍｂｄプログラムによるシミュレーション方法、数値解析装置、ｍｂｄ用数値解析システム、数値解析プログラムおよびｍｂｄプログラム
JP5454557B2 (ja) 物理量計算方法、数値解析方法、物理量計算プログラム、数値解析プログラム、物理量計算装置及び数値解析装置
Gano et al. Hybrid variable fidelity optimization by using a kriging-based scaling function
Kroll et al. DLR project Digital-X: towards virtual aircraft design and flight testing based on high-fidelity methods
Bueno-Orovio et al. Continuous adjoint approach for the Spalart-Allmaras model in aerodynamic optimization
WO2020054086A1 (ja) シミュレーション方法、ｍｂｄプログラムによるシミュレーション方法、数値解析装置、ｍｂｄ用数値解析システム、数値解析プログラムおよびｍｂｄプログラム
Ramlan et al. Comparison between solidworks and ansys flow simulation on aerodynamic studies
JP6504333B1 (ja) シミュレーション方法、物理量計算プログラム及び物理量計算装置
CN114154294A (zh) 一种航空发动机高空模拟试验的虚拟试验架构
CN111950210A (zh) 一种基于cfd的商用车发动机舱热管理模拟方法
Ciampa et al. CST parametrization for unconventional aircraft design optimization
WO2020054087A1 (ja) シミュレーション方法、物理量計算プログラム及び物理量計算装置
Hajdik et al. Aerodynamic shape optimization with CAD-based geometric parameterization
Poole et al. Optimal domain element shapes for free-form aerodynamic shape control
Mosavi Multiobjective Optimization of Spline Curves Using Mode Frontier
Han et al. Robust functional response-based metamodel optimization considering both location and dispersion effects for aeronautical airfoil designs
Rausa et al. Multi-fidelity assessment of the aerodynamic performances of unconventional aircraft tail configurations
Tarsia Morisco et al. Extension of the Vertex-Centered Mixed-Element-Volume MUSCL scheme to mixed-element meshes
Mengaldo et al. Immersed boundary lattice Green function methods for external aerodynamics
Liatsikouras et al. A coupled CAD‐free parameterization‐morphing method for adjoint‐based shape optimization
EP3394773A1 (en) Composite design direction
Barutcu et al. A Parametric Modeling Approach for Prediction of Load Distribution due to Fluid Structure Interaction on Aircraft Structures
Kosík The CFD simulation of the flow around the aircraft using OpenFOAM and ANSA
Forsman et al. Numerical simulation of novel cooling fan concept for BEVs. Investigation of aerodynamic performance of novel cooling fan concept and comparison with conventional fan concept using Computational Fluid Dynamics
Chernukhin Global optimization algorithms for aerodynamic design

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20190207

A871 Explanation of circumstances concerning accelerated examination

Free format text: JAPANESE INTERMEDIATE CODE: A871

Effective date: 20190207

A975 Report on accelerated examination

Free format text: JAPANESE INTERMEDIATE CODE: A971005

Effective date: 20190220

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20190226

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20190311

R150 Certificate of patent or registration of utility model

Ref document number: 6504333

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250