ホーム > 特許ランキング > 住友重機械工業株式会社
知財求人 - 知財ポータルサイト「IP Force」
特開2023-102969シミュレーション方法、シミュレーション装置、及びプログラム
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】公開特許公報(A)
(11)【公開番号】P2023102969
(43)【公開日】2023-07-26
(54)【発明の名称】シミュレーション方法、シミュレーション装置、及びプログラム
(51)【国際特許分類】
G16Z 99/00 20190101AFI20230719BHJP
【ＦＩ】
G16Z99/00
【審査請求】未請求
【請求項の数】7
【出願形態】ＯＬ
(21)【出願番号】P 2022003755
(22)【出願日】2022-01-13
(71)【出願人】
【識別番号】000002107
【氏名又は名称】住友重機械工業株式会社
(74)【代理人】
【識別番号】100105887
【弁理士】
【氏名又は名称】来山 幹雄
(72)【発明者】
【氏名】廣田 真人
【テーマコード（参考）】
5L049
【Ｆターム（参考）】
5L049DD02
(57)【要約】
【課題】液膜流れの数値解析において、計算負荷を低減させることが可能なシミュレーション方法を提供する。
【解決手段】ＳＰＨ法を用いて浅水方程式を離散化して液膜流れを数値的に解析するシミュレーション方法において、まず、物理空間上で液膜が流れる曲面の形状を定義する。次に、曲面上の位置を二次元の解析空間上の位置に対応付ける。物理空間に定義される物理座標系と、解析空間に定義される一般座標系との間で物理量を変換する一般座標変換を、離散化された浅水方程式に施す。一般座標変換された浅水方程式を数値的に解く。
【選択図】図６
【解決手段】ＳＰＨ法を用いて浅水方程式を離散化して液膜流れを数値的に解析するシミュレーション方法において、まず、物理空間上で液膜が流れる曲面の形状を定義する。次に、曲面上の位置を二次元の解析空間上の位置に対応付ける。物理空間に定義される物理座標系と、解析空間に定義される一般座標系との間で物理量を変換する一般座標変換を、離散化された浅水方程式に施す。一般座標変換された浅水方程式を数値的に解く。
【選択図】図６
【特許請求の範囲】
【請求項1】
ＳＰＨ法を用いて浅水方程式を離散化して液膜流れを数値的に解析するシミュレーション方法であって、
物理空間上で液膜が流れる曲面の形状を定義し、
前記曲面上の位置を二次元の解析空間上の位置に対応付け、
前記物理空間に定義される物理座標系と、前記解析空間に定義される一般座標系との間で物理量を変換する一般座標変換を、離散化された前記浅水方程式に施し、
前記一般座標変換された前記浅水方程式を数値的に解くシミュレーション方法。
物理空間上で液膜が流れる曲面の形状を定義し、
前記曲面上の位置を二次元の解析空間上の位置に対応付け、
前記物理空間に定義される物理座標系と、前記解析空間に定義される一般座標系との間で物理量を変換する一般座標変換を、離散化された前記浅水方程式に施し、
前記一般座標変換された前記浅水方程式を数値的に解くシミュレーション方法。
【請求項2】
前記曲面の形状を、複数の頂点の座標を含むＣＡＤデータで定義し、前記頂点の座標を前記物理空間における座標として用い、
前記ＣＡＤデータに含まれるテクスチャ座標系を前記解析空間の座標系として用い、
前記ＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、前記物理空間上の物理座標と前記解析空間上の一般座標との間の座標変換係数を算出し、
算出された座標変換係数を用いて、前記一般座標変換を行う請求項１に記載のシミュレーション方法。
前記ＣＡＤデータに含まれるテクスチャ座標系を前記解析空間の座標系として用い、
前記ＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、前記物理空間上の物理座標と前記解析空間上の一般座標との間の座標変換係数を算出し、
算出された座標変換係数を用いて、前記一般座標変換を行う請求項１に記載のシミュレーション方法。
【請求項3】
物理空間における曲面の形状を、複数の頂点の座標を含むＣＡＤデータで定義し、前記頂点の座標を前記物理空間における座標として用い、
前記ＣＡＤデータに含まれるテクスチャ座標を解析空間の座標として用い、
支配方程式を数値的に解くシミュレーション方法。
前記ＣＡＤデータに含まれるテクスチャ座標を解析空間の座標として用い、
支配方程式を数値的に解くシミュレーション方法。
【請求項4】
前記支配方程式を、ＳＰＨ法を用いて数値的に解き、
さらに、前記ＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを利用して、前記物理空間と前記解析空間との間の座標変換係数を算出し、
算出された座標変換係数を用いて、前記物理空間上の物理座標と前記解析空間上の一般座標とを座標変換を行う請求項３に記載のシミュレーション方法。
さらに、前記ＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを利用して、前記物理空間と前記解析空間との間の座標変換係数を算出し、
算出された座標変換係数を用いて、前記物理空間上の物理座標と前記解析空間上の一般座標とを座標変換を行う請求項３に記載のシミュレーション方法。
【請求項5】
物理空間の曲面上の位置を二次元の解析空間上の位置に対応付けて、ＳＰＨ法を用いて解析を行い、得られた物理量を前記物理空間の物理量に変換するシミュレーション方法であって、
前記曲面の形状を複数の頂点の座標として定義するＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、前記物理空間上の物理座標系と前記解析空間上の一般座標系との間の座標変換係数を算出するシミュレーション方法。
前記曲面の形状を複数の頂点の座標として定義するＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、前記物理空間上の物理座標系と前記解析空間上の一般座標系との間の座標変換係数を算出するシミュレーション方法。
【請求項6】
ＳＰＨ法を用いて液膜の流れを解析するシミュレーション装置であって、
シミュレーション条件が入力される入力部と、
前記入力部に入力されたシミュレーション条件に基づいて、ＳＰＨ法を用いて複数の粒子のそれぞれの位置を時間発展させる処理部と
を備え、
前記シミュレーション条件は、液膜が沿う曲面の形状を物理空間に定義するＣＡＤデータ、液膜の物性値、及び初期条件を含み、
前記ＣＡＤデータは、前記物理空間に定義される物理座標系において、前記曲面の形状を複数の頂点の座標として定義する頂点の一覧、及び前記曲面に張られる面に定義されるテクスチャ座標系におけるテクスチャ座標の一覧、及び前記物理座標系と前記テクスチャ座標系との対応関係を定義する情報を含み、
前記処理部は、
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出し、
前記解析空間に複数の計算粒子を配置し、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返し、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換するシミュレーション装置。
シミュレーション条件が入力される入力部と、
前記入力部に入力されたシミュレーション条件に基づいて、ＳＰＨ法を用いて複数の粒子のそれぞれの位置を時間発展させる処理部と
を備え、
前記シミュレーション条件は、液膜が沿う曲面の形状を物理空間に定義するＣＡＤデータ、液膜の物性値、及び初期条件を含み、
前記ＣＡＤデータは、前記物理空間に定義される物理座標系において、前記曲面の形状を複数の頂点の座標として定義する頂点の一覧、及び前記曲面に張られる面に定義されるテクスチャ座標系におけるテクスチャ座標の一覧、及び前記物理座標系と前記テクスチャ座標系との対応関係を定義する情報を含み、
前記処理部は、
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出し、
前記解析空間に複数の計算粒子を配置し、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返し、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換するシミュレーション装置。
【請求項7】
ＣＡＤデータ、液膜の物性値、及び初期条件を含むシミュレーション条件を取得する手順であって、前記ＣＡＤデータは、物理空間に定義される物理座標系において、曲面の形状を複数の頂点の座標として定義する頂点の一覧、及び前記曲面に張られる面に定義されるテクスチャ座標系におけるテクスチャ座標の一覧、及び前記物理座標系と前記テクスチャ座標系との対応関係を定義する情報を含む手順と、
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出する手順と、
前記解析空間に複数の計算粒子を配置する手順と、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返す手順と、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換する手順と
をコンピュータに実行させるプログラム。
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出する手順と、
前記解析空間に複数の計算粒子を配置する手順と、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返す手順と、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換する手順と
をコンピュータに実行させるプログラム。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、シミュレーション方法、シミュレーション装置、及びプログラムに関する。
【背景技術】
【0002】
蒸気タービン内部において、非平衡凝縮によって蒸気の流れの中に水滴が発生し、発生した水滴が翼の表面に付着する。付着した水滴は液膜となって翼の表面を流れ、後縁に溜まる。最終的には、翼の後縁に溜まった液膜が周囲の気流によって翼からちぎれて飛散する。このような蒸気の流れをはじめとする湿り蒸気流れ等の気液二相流体の流れ場のシミュレーションにおいて、圧縮性流れの解法と、古典凝縮論とを組み合わせた数値解析法が用いられている（非特許文献１参照）。従来技術では、解析空間を格子状に分割し、圧縮性ナビエストークス方程式と凝集モデルからなる支配方程式を解くことにより、流れ場の時間発展を計算する。
【先行技術文献】
【非特許文献】
【0003】
【非特許文献1】F. Bakhtar and M.T. Mohammadi Tochai,"An investigation of two-dimensional flows of nucleating and wet steam bythe time-marching method", International Journal of Heat and Fluid Flow,vol. 2 (1), pp.5-18 (1980)
【発明の概要】
【発明が解決しようとする課題】
【0004】
液膜の流れを三次元のＳＰＨ法を用いて解くためには、液膜の厚さに対して十分小さい直径の計算粒子の集合で液膜を再現する必要がある。このため、計算粒子の個数が多くなり、計算負荷が高くなってしまう。
【0005】
また、液膜が沿う任意の形状の曲面を離散化するための格子を、曲面に沿って張ることが困難である。さらに、濡れている領域と濡れていない領域との境界で、液膜高さの勾配が不連続になるため、濡れている領域と濡れていない領域とを区別し、両者の境界において特別な処理を行うことが必要になる。ところが、濡れている領域と濡れていない領域とを区別することは困難である。
【0006】
本発明の目的は、液膜流れの数値解析において、計算負荷を低減させることが可能なシミュレーション方法、シミュレーション装置、及びプログラムを提供することである。本発明の他の目的は、液膜流れが沿う曲面に格子を張ることなく、解析を行うことが可能なシミュレーション方法、シミュレーション装置、及びプログラムを提供することである。本発明のさらに他の目的は、濡れている領域と濡れていない領域とを区別する必要がないシミュレーション方法、シミュレーション装置、及びプログラムを提供することである。
【課題を解決するための手段】
【0007】
本発明の一観点によると、
ＳＰＨ法を用いて浅水方程式を離散化して液膜流れを数値的に解析するシミュレーション方法であって、
物理空間上で液膜が流れる曲面の形状を定義し、
前記曲面上の位置を二次元の解析空間上の位置に対応付け、
前記物理空間に定義される物理座標系と、前記解析空間に定義される一般座標系との間で物理量を変換する一般座標変換を、離散化された前記浅水方程式に施し、
前記一般座標変換された前記浅水方程式を数値的に解くシミュレーション方法が提供される。
ＳＰＨ法を用いて浅水方程式を離散化して液膜流れを数値的に解析するシミュレーション方法であって、
物理空間上で液膜が流れる曲面の形状を定義し、
前記曲面上の位置を二次元の解析空間上の位置に対応付け、
前記物理空間に定義される物理座標系と、前記解析空間に定義される一般座標系との間で物理量を変換する一般座標変換を、離散化された前記浅水方程式に施し、
前記一般座標変換された前記浅水方程式を数値的に解くシミュレーション方法が提供される。
【0008】
本発明の他の観点によると、
物理空間における曲面の形状を、複数の頂点の座標を含むＣＡＤデータで定義し、前記頂点の座標を前記物理空間における座標として用い、
前記ＣＡＤデータに含まれるテクスチャ座標を解析空間の座標として用い、
支配方程式を数値的に解くシミュレーション方法が提供される。
物理空間における曲面の形状を、複数の頂点の座標を含むＣＡＤデータで定義し、前記頂点の座標を前記物理空間における座標として用い、
前記ＣＡＤデータに含まれるテクスチャ座標を解析空間の座標として用い、
支配方程式を数値的に解くシミュレーション方法が提供される。
【0009】
本発明のさらに他の観点によると、
物理空間の曲面上の位置を二次元の解析空間上の位置に対応付けて、ＳＰＨ法を用いて解析を行い、得られた物理量を前記物理空間の物理量に変換するシミュレーション方法であって、
前記曲面の形状を複数の頂点の座標として定義するＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、前記物理空間上の物理座標系と前記解析空間上の一般座標系との間の座標変換係数を算出するシミュレーション方法が提供される。
物理空間の曲面上の位置を二次元の解析空間上の位置に対応付けて、ＳＰＨ法を用いて解析を行い、得られた物理量を前記物理空間の物理量に変換するシミュレーション方法であって、
前記曲面の形状を複数の頂点の座標として定義するＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、前記物理空間上の物理座標系と前記解析空間上の一般座標系との間の座標変換係数を算出するシミュレーション方法が提供される。
【0010】
本発明のさらに他の観点によると、
ＳＰＨ法を用いて液膜の流れを解析するシミュレーション装置であって、
シミュレーション条件が入力される入力部と、
前記入力部に入力されたシミュレーション条件に基づいて、ＳＰＨ法を用いて複数の粒子のそれぞれの位置を時間発展させる処理部と
を備え、
前記シミュレーション条件は、液膜が沿う曲面の形状を物理空間に定義するＣＡＤデータ、液膜の物性値、及び初期条件を含み、
前記ＣＡＤデータは、前記物理空間に定義される物理座標系において、前記曲面の形状を複数の頂点の座標として定義する頂点の一覧、及び前記曲面に張られる面に定義されるテクスチャ座標系におけるテクスチャ座標の一覧、及び前記物理座標系と前記テクスチャ座標系との対応関係を定義する情報を含み、
前記処理部は、
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出し、
前記解析空間に複数の計算粒子を配置し、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返し、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換するシミュレーション装置が提供される。
ＳＰＨ法を用いて液膜の流れを解析するシミュレーション装置であって、
シミュレーション条件が入力される入力部と、
前記入力部に入力されたシミュレーション条件に基づいて、ＳＰＨ法を用いて複数の粒子のそれぞれの位置を時間発展させる処理部と
を備え、
前記シミュレーション条件は、液膜が沿う曲面の形状を物理空間に定義するＣＡＤデータ、液膜の物性値、及び初期条件を含み、
前記ＣＡＤデータは、前記物理空間に定義される物理座標系において、前記曲面の形状を複数の頂点の座標として定義する頂点の一覧、及び前記曲面に張られる面に定義されるテクスチャ座標系におけるテクスチャ座標の一覧、及び前記物理座標系と前記テクスチャ座標系との対応関係を定義する情報を含み、
前記処理部は、
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出し、
前記解析空間に複数の計算粒子を配置し、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返し、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換するシミュレーション装置が提供される。
【0011】
本発明のさらに他の観点によると、
ＣＡＤデータ、液膜の物性値、及び初期条件を含むシミュレーション条件を取得する手順であって、前記ＣＡＤデータは、物理空間に定義される物理座標系において、曲面の形状を複数の頂点の座標として定義する頂点の一覧、及び前記曲面に張られる面に定義されるテクスチャ座標系におけるテクスチャ座標の一覧、及び前記物理座標系と前記テクスチャ座標系との対応関係を定義する情報を含む手順と、
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出する手順と、
前記解析空間に複数の計算粒子を配置する手順と、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返す手順と、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換する手順と
をコンピュータに実行させるプログラムが提供される。
ＣＡＤデータ、液膜の物性値、及び初期条件を含むシミュレーション条件を取得する手順であって、前記ＣＡＤデータは、物理空間に定義される物理座標系において、曲面の形状を複数の頂点の座標として定義する頂点の一覧、及び前記曲面に張られる面に定義されるテクスチャ座標系におけるテクスチャ座標の一覧、及び前記物理座標系と前記テクスチャ座標系との対応関係を定義する情報を含む手順と、
前記テクスチャ座標系が定義された二次元平面を解析空間とし、前記テクスチャ座標系を前記解析空間に定義された一般座標系として用い、前記ＣＡＤデータに基づいて前記物理座標系と前記一般座標系との座標変換係数を算出する手順と、
前記解析空間に複数の計算粒子を配置する手順と、
前記物理空間上で支配方程式を数値的に解いて前記計算粒子の各々の速度を計算し、前記計算粒子の各々の速度を、前記座標変換係数を用いて座標変換して前記解析空間上の速度を求め、前記解析空間上で前記複数の計算粒子の各々の位置を更新する処理を繰り返す手順と、
求められた前記複数の計算粒子の前記解析空間上の位置を、前記座標変換係数を用いて座標変換し、前記物理空間上における位置に変換する手順と
をコンピュータに実行させるプログラムが提供される。
【発明の効果】
【0012】
液膜流れの解析に二次元の浅水方程式を用いることにより、三次元で解析する場合と比べて計算負荷を軽減させることができる。ＣＡＤデータのテクスチャ座標を解析空間の一般座標として用いることにより、液膜流れが沿う曲面に格子を張ることなく、かつ物理空間に定義される物理座標系と解析空間に定義される一般座標系との間の座標変換係数を容易に算出することができる。ＳＰＨ法を適用することにより、曲面の濡れている領域と濡れていない領域とを区別する必要がなくなる。
【図面の簡単な説明】
【0013】
【発明を実施するための形態】
【0014】
本願発明の実施例について説明する前に、浅水方程式について説明する。
［浅水方程式］
本願発明の一実施例においては、三次元のナビエストークス方程式に代えて、浅水方程式を使用する。浅水方程式を用いることにより、液膜流れを二次元的に解析することができる。このため、計算負荷が低減される。
［浅水方程式］
本願発明の一実施例においては、三次元のナビエストークス方程式に代えて、浅水方程式を使用する。浅水方程式を用いることにより、液膜流れを二次元的に解析することができる。このため、計算負荷が低減される。
【0015】
浅水方程式は、以下の式で与えられる。
ここで、Ｈは液膜の厚さ、ｔは時間、Ｄ／Ｄｔは時間での全微分、ρは密度、ｕは速度ベクトル、ｐは圧力、μは粘性係数、ｇは重力ベクトル、ｎは液面が流れる壁面（曲面）の法線ベクトル、σは表面張力係数、κは液膜の曲率、ｇｔは重力ベクトルの曲面に対する接線方向成分、θは液膜の接触角、ｎｃｌは液膜の縁の法線ベクトルを表す。式（１）の３番目の式の右辺第１項は重力による静水圧力を表し、第２項は液膜の表面張力によるラプラス圧を表す。式（１）の４番目の式の右辺第１項は、液膜に作用する重力加速度、第２項は粘性摩擦力、第３項は液膜の縁に作用する表面張力を表す。
【数1】
【0016】
液面の曲率κは、以下の式で表すことができる。
ここで、下付き文字ｘ、ｙは、対応する空間方向の導関数を表す。例えば、Ｈｘ、Ｈｘｘは、以下の式で表される。
【数2】
【数3】
【0017】
液膜の厚さＨが十分薄い場合、その１次導関数Ｈｘ、Ｈｙは１より十分小さくなるため、Ｈｘの２乗の項及びＨｙの２乗の項は無視することができる。したがって、式（２）は、以下のように簡略化することができる。
【数4】
【0018】
［ＳＰＨ法の概要］
本願発明の一実施例において、ＳＰＨ法が用いられる。次に、一般的なＳＰＨ法について簡単に説明する。ＳＰＨ法は、連続体を複数の粒子の集合を用いて離散化し、流れ場の時間発展を計算するラグランジュ型の数値計算手法である。粒子１個の質量分布として、連続かつ微分可能なカーネル関数を用いるとともに、カーネル関数の微分係数の重ね合わせを空間微分とする離散化モデルである。
本願発明の一実施例において、ＳＰＨ法が用いられる。次に、一般的なＳＰＨ法について簡単に説明する。ＳＰＨ法は、連続体を複数の粒子の集合を用いて離散化し、流れ場の時間発展を計算するラグランジュ型の数値計算手法である。粒子１個の質量分布として、連続かつ微分可能なカーネル関数を用いるとともに、カーネル関数の微分係数の重ね合わせを空間微分とする離散化モデルである。
【0019】
ＳＰＨ法では、粒子の質量ｍはカーネル関数Ｗによってその周囲に分配され、連続的に分布しており、密度場は、その重ね合わせとして与えられる。すなわち、ｉ番目の粒子の位置ｘｉにおける密度ρｉは、以下の式で表される。
カーネル関数Ｗとして、例えば、以下の関数を用いることができる。以下の関数は、Ｃ２Ｗｅｎｄｌａｎｄカーネルと呼ばれる。
【数5】
【数6】
【0020】
ここで、ｘｉｊは位置ベクトルｘｉと位置ベクトルｘｊとの差の絶対値であり、ｑ＝ｘｉｊ／ｈである。ｈは、平滑化長さであり、カーネル幅ともいわれる。カーネル幅ｈは、複数の粒子を同一の質量及び物性密度の球体に置き換えた等価球の直径に等しい。本明細書において、カーネル幅ｈを粒子の直径という場合がある。
【0021】
任意の位置ｘｉにおける任意のスカラー量Ａ及びベクトル量Ａの値は、カーネル関数の重ね合わせによって以下の式で表すことができる。
ここで、Ｖｊ＝ｍｊ／ρｊである。下付き文字ｊは、ｊ番目の粒子の持つ物理量の値であることを意味する。
【数7】
【0022】
任意のスカラー量Ａの勾配は、カーネル関数の微分係数の重ね合わせによって以下の式で表される。
【数8】
【0023】
ここで、∇Ｗ（ｘｉｊ，ｈ）はカーネル関数の導関数であり、カーネル関数としてＣ２Ｗｅｎｄｌａｎｄカーネルを用いる場合、以下の式で表される。
【数9】
【0024】
任意のベクトル量Ａに対する発散は、式（８）と同様に、以下の式で表すことができる。
【数10】
【0025】
任意のスカラー量Ａに対する２次導関数は、例えば以下の公式を用いて算出することができる。
【数11】
【0026】
［ＳＰＨ法を用いた浅水方程式の離散化］
浅水方程式（式（１））を、ＳＰＨ法を用いて離散化すると、以下のように表すことができる。
浅水方程式（式（１））を、ＳＰＨ法を用いて離散化すると、以下のように表すことができる。
【数12】
【0027】
［参考例による解析方法］
本願発明の一実施例について説明する前に、参考例による解析方法について説明する。
一参考例では、物理空間、すなわち液面流れが沿う解析対象の曲面を格子状に分割し、ナビエストークス方程式または浅水方程式（式（１））を数値的に解くことにより、流れ場の時間発展を計算する。その際、支配方程式は物理空間に張られた格子を直交等間隔格子に射影した空間を解析空間として定義し、物理空間と解析空間との間の座標変換係数を算出する。解析空間上で支配方程式を解き、その結果を物理空間上の物理量に変換する手法が用いられる。
本願発明の一実施例について説明する前に、参考例による解析方法について説明する。
一参考例では、物理空間、すなわち液面流れが沿う解析対象の曲面を格子状に分割し、ナビエストークス方程式または浅水方程式（式（１））を数値的に解くことにより、流れ場の時間発展を計算する。その際、支配方程式は物理空間に張られた格子を直交等間隔格子に射影した空間を解析空間として定義し、物理空間と解析空間との間の座標変換係数を算出する。解析空間上で支配方程式を解き、その結果を物理空間上の物理量に変換する手法が用いられる。
【0028】
本明細書において、物理空間上で用いられる座標系を物理座標系といい、解析空間上で用いられる座標系を一般座標系という。物理座標系と一般座標系との間で物理量を変換することを一般座標変換という。
【0029】
工学上問題となる曲面の形状は、ＣＡＤデータで表現される場合が多い。このとき、一般座標変換の座標変換係数を算出するためには、物理空間に張られた格子と解析空間に張られた格子との間の整合性が保たれていなければならない。ＣＡＤデータを作成する段階でその整合性を考慮しながらＣＡＤのモデリングを行うことは困難である。両者の整合性が保たれていない場合は、液膜流れのシミュレーションにおいて、一般座標変換の座標変換係数を算出する段階でＣＡＤデータを修正したり、物理空間と解析空間との間の不整合性を解消したりするといった余分な手順が必要になる。以下に説明する実施例では、これらの種々の課題が解決される。
【0030】
［汎用的なＣＡＤデータの構成］
次に、汎用的なＣＡＤデータの構成について説明する。ＣＡＤデータは、物体の形状を表現するための点列データ、物体表面の材質（質感）を表現するための点列データ、及びそれらを関係付けるための一覧表によって構成される。一例として、ＣＡＤデータとして用いられているＷａｖｅｆｒｏｎｔ形式のデータ構造について説明する。
次に、汎用的なＣＡＤデータの構成について説明する。ＣＡＤデータは、物体の形状を表現するための点列データ、物体表面の材質（質感）を表現するための点列データ、及びそれらを関係付けるための一覧表によって構成される。一例として、ＣＡＤデータとして用いられているＷａｖｅｆｒｏｎｔ形式のデータ構造について説明する。
【0031】
図１は、Ｗａｖｅｆｒｏｎｔ形式のＣＡＤデータのデータ構造の一例を示す図表である。図２は、物理空間と解析空間との関係を示す模式図である。物理空間Ｓｒ（図２）に曲面Ｆｃが定義されている。図１において、符号＃は、以下に続く文字列がコメントであることを表す。ＣＡＤデータは、物理空間の曲面上の複数の点（以下、「頂点」という。）の一覧、テクスチャ座標の一覧、物理空間の曲面の法線ベクトル（以下、「頂点法線」という。）の一覧、及び面要素の一覧を含む。
【0032】
図１に示した符号ｖ以下の３つの実数の組み合わせは、物理空間Ｓｒ（図２）の曲面Ｆｃ上の複数の点の座標（Ｘ，Ｙ，Ｚ）を表す。以下、曲面Ｆｃ上の点の座標を「頂点座標」という。符号ｖ以下の３つの実数の組み合わせのそれぞれに要素番号が付与されており、要素番号によって一つの頂点座標が特定される。
【0033】
図１に示した符号ｖｔ以下の２つの実数の組み合わせは、テクスチャ座標（ξ，η）（図２）の一覧を表す。符号ｖｔ以下の２つの実数の組み合わせのそれぞれに要素番号が付与されており、要素番号によって一つのテクスチャ座標（ξ，η）が特定される。各要素番号によって特定されるテクスチャ座標系における位置を「要素点」という。
【0034】
図１に示した符号ｎｖ以下の３つの実数の組み合わせは、曲面Ｆｃ上の点の法線ベクトルｎ（図２）を表す。符号ｎｖ以下の３つの実数の組み合わせのそれぞれに要素番号が付与されており、要素番号によって一つの法線ベクトルが特定される。
【0035】
図１に示した符号ｆ以下では、頂点座標（Ｘ，Ｙ，Ｚ）、テクスチャ座標（ξ，η）、頂点法線ｎの要素番号を表す整数が斜線記号を隔てて記述されている。斜線記号を隔てて記述された３つまたは４つの要素番号の組み合わせによって、三角形または四角形の物理空間Ｓｒ上の要素Ｅｒ、及びそれに対応付けられた解析空間Ｓａ上の要素Ｅａが定義される。
【0036】
図１に示した例では、要素番号が１の頂点座標と要素番号が６のテクスチャ座標と要素番号が２の頂点法線、要素番号が２の頂点座標と要素番号が３のテクスチャ座標と要素番号が１の頂点法線、要素番号が４の頂点座標と要素番号が５のテクスチャ座標と要素番号が６の頂点法線によって、一つの三角形の要素が定義される。
【0037】
物理空間Ｓｒの頂点座標（Ｘ，Ｙ，Ｚ）は、物体の輪郭を記述するための代表点であり、複数の代表点を接続することによって曲面Ｆｃを形成することができる。テクスチャ座標（ξ，η）は、物体の表面に張られる画像データの座標点を記述したものであり、二次元座標上に定義される画像データの各座標の値が｛０，１｝に正規化されたものである。頂点法線ｎは、頂点座標における曲面Ｆｃの法線方向ベクトルを記述したものである。
【0038】
上述の通り、汎用的なＣＡＤデータは、物理空間Ｓｒの座標（Ｘ，Ｙ，Ｚ）とテクスチャ座標（ξ，η）との対応関係が指定されている。テクスチャ座標は、｛０，１｝に正規化されているため、解析空間Ｓａの座標系として用いるのに適している。以下に説明する一実施例では、テクスチャ座標系が解析空間Ｓａの一般座標系として用いられる。
【0039】
［ＣＡＤデータからの座標変換係数の導出］
物理空間Ｓｒに定義された物理座標系における曲面Ｆｃ（図２）の基底ベクトルを（ｅｘ，ｅｙ）と標記し、一般座標系（テクスチャ座標系）の基底ベクトルを（ｅξ，ｅη）と標記する。物理空間Ｓｒの曲面Ｆｃ上の座標を（ｘ，ｙ）と標記する。解析空間Ｓａに定義される一般座標系の座標を（ξ，η）と標記する。物理空間Ｓｒと解析空間Ｓａとの間の座標変換係数は、カーネル関数の導関数の公式を用いて以下の式で算出することができる。
ここで、下付き文字ｉ、ｊは、ＣＡＤデータに格納されているテクスチャ座標の要素番号を表し、一般座標系における位置座標に対応づけられる。ベクトルｘｉは、テクスチャ座標の要素番号ｉと結びついている物理空間上の位置座標である。ベクトルｅｘｉは、物理空間上の位置座標ｘｉにおける曲面の基底ベクトルである。ｍ０、ρ０は、カーネル関数の公式を適用するための仮想的な質量及び密度であり、例えばそれぞれ１に設定することができる。
物理空間Ｓｒに定義された物理座標系における曲面Ｆｃ（図２）の基底ベクトルを（ｅｘ，ｅｙ）と標記し、一般座標系（テクスチャ座標系）の基底ベクトルを（ｅξ，ｅη）と標記する。物理空間Ｓｒの曲面Ｆｃ上の座標を（ｘ，ｙ）と標記する。解析空間Ｓａに定義される一般座標系の座標を（ξ，η）と標記する。物理空間Ｓｒと解析空間Ｓａとの間の座標変換係数は、カーネル関数の導関数の公式を用いて以下の式で算出することができる。
【数13】
【0040】
ｈ０は仮想的なカーネル幅であり、この値を大きくすると精度が向上する一方で、空間解像度が低下する。実用上は、仮想的なカーネル幅ｈ０は、例えば隣接するテクスチャ座標要素同士の最大距離に設定することができる。テクスチャ座標で定義される二次元平面を解析空間として用いることで、ＣＡＤデータに格納されている情報のみで、一般座標変換の座標変換係数を算出することができる。
【0041】
２次導関数の補間公式も同様に算出することができ、導関数の一般座標変換係数は、以下のように求められる。
【数14】
【0042】
［離散化された浅水方程式に対する一般座標変換］
次に、離散化された浅水方程式に一般座標変換を施す方法について説明する。
物理空間において、ＳＰＨ法を用いて離散化された浅水方程式は以下のように表される。
ここで、ｕｘ及びｕｙは、それぞれ速度ベクトルｕのｘ方向成分及びｙ方向成分を表す。また、Ｊｊは座標変換のヤコビアンを表す。ヤコビアンＪｊは以下の式で定義される。
次に、離散化された浅水方程式に一般座標変換を施す方法について説明する。
物理空間において、ＳＰＨ法を用いて離散化された浅水方程式は以下のように表される。
【数15】
【数16】
【数17】
【数18】
【数19】
【0043】
また、導関数を含む項は、解析空間上での導関数を求めた後に、物理空間上の値へと変換がなされている。式（１５）～式（１８）で得られる物理量は、物理空間上での値であるため、計算粒子の位置座標を更新する際には、それらを解析空間上における値に変換する必要がある。すなわち、計算粒子の位置座標は、以下の式を用いて時間発展させる。
【数20】
【0044】
解析空間上の任意の座標ξに対応する物理空間上のスカラー量Ａ及びベクトル量Ａは、以下の式を用いて算出することができる。
ここで、スカラー量Ａには、座標変換係数∂ξ／∂ｘ、∂ｘ／∂ξ等が含まれる。ベクトル量Ａには、計算粒子の物理空間上の位置ベクトルｘ、物理空間上の曲面の法線ベクトルｎ等が含まれる。式（１３）、式（１４）を用いて算出される座標変換係数は、ＣＡＤデータ（図１）に格納されている要素点の位置において定義される。式（１３）、式（１４）を用いて算出された要素点ξｊの座標変換係数から、式（２１）を用いて任意の位置ξの座標変換係数を算出することができる。また、ＣＡＤデータ（図１）で定義された要素点ξｊにおける物理空間上の位置ベクトルｘ及び法線ベクトルｎから、それぞれ式（２１）を用いて任意の位置ξにおける位置ベクトル及び法線ベクトルを算出することができる。
【数21】
【0045】
【0046】
図３は、一実施例によるシミュレーション装置のブロック図である。本実施例によるシミュレーション装置は、入力部４０、処理部４１、出力部４２、及び記憶部４３を含む。入力部４０から処理部４１にシミュレーション条件等が入力される。さらに、ユーザから入力部４０に各種指令（コマンド）等が入力される。入力部４０は、例えば通信装置、リムーバブルメディア読取装置、キーボード、ポインティングデバイス等で構成される。
【0047】
処理部４１は、入力されたシミュレーション条件及び指令に基づいてシミュレーションを実行する。さらに、シミュレーション結果を出力部４２に出力する。シミュレーション結果には、例えば、シミュレーション対象物を構成する粒子の位置の時間変化等を表す情報が含まれる。処理部４１は、例えばコンピュータの中央処理ユニット（ＣＰＵ）を含む。本実施例によるシミュレーション方法の各手順をコンピュータに実行させるためのプログラムが、記憶部４３に記憶されている。出力部４２は、通信装置、リムーバブルメディア書込み装置、ディスプレイ等を含む。
【0048】
図４は、本実施例によるシミュレーション方法の手順を示すフローチャートである。処理部４１（図１）が記憶部４３（図１）に格納されているプログラムを実行することにより、図４に示した各ステップの処理が実行される。
【0049】
まず、処理部４１が、ユーザから入力されたシミュレーション条件を取得する（ステップＳ０１）。シミュレーション条件には、シミュレーション対象の流体を定義する情報、初期条件、境界条件、解析終了条件等が含まれる。シミュレーション対象の流体を定義する情報には、流体の物性値、例えば密度、粘度等を指定する情報が含まれる。初期条件には、複数の計算粒子を解析空間に配置するための情報、計算粒子の初期速度を指定する方法等が含まれる。境界条件には、解析対象の空間の形状及び大きさを指定する情報、曲面の形状を指定するＣＡＤデータ等が含まれる。ＣＡＤデータによって、図２に示した物理空間Ｓｒ内の曲面Ｆｃの形状、曲面Ｆｃ上の複数の要素Ｅｒと、解析空間上の複数の要素Ｅａとの対応関係が定義される。
【0050】
次に、取得したＣＡＤデータに基づいて、物理空間Ｓｒと解析空間Ｓａとの間の座標変換係数を算出する（ステップＳ０２）。例えば、式（１３）及び式（１４）に示した変換係数を計算する。ステップＳ０２では、ＣＡＤデータで定義された複数の要素点のそれぞれの位置における座標変換係数が計算される。
【0051】
次に、入力された初期条件に基づいて、解析空間Ｓａに複数の計算粒子を配置する（ステップＳ０３）。ステップＳ０２で求められた各要素点の位置における座標変換係数から、式（２１）を用いて、計算粒子の各々の位置における座標変換係数を算出する（ステップＳ０４）。解析空間Ｓａ上における計算粒子の物理量の導関数の値を計算する（ステップＳ０５）。例えば、式（１６）、式（１７）のそれぞれの右辺のΣの値、すなわち計算粒子の圧力の導関数及び速度の導関数の値を計算する。ステップＳ０５で求められた解析空間Ｓａ上における計算粒子の物理量の導関数の値を、物理空間Ｓｒにおける値に変換する（ステップＳ０６）。例えば、式（１６）、式（１７）の右辺に示されているように、各項のΣの値のそれぞれに変換係数を乗じることにより、物理空間Ｓｒにおける計算粒子の圧力の導関数及び速度の導関数の値を求める。
【0052】
次に、物理空間Ｓｒで計算粒子に付与される速度を計算する（ステップＳ０７）。例えば、式（１６）、式（１７）を計算し、積分することにより、各計算粒子の物理空間Ｓｒにおける速度ｕｘ、ｕｙを求める。次に、物理空間Ｓｒにおける速度の値を、解析空間Ｓａにおける値に変換し、計算粒子の物理量及び位置を更新する（ステップＳ０８）。例えば、式（２０）の右辺に示されているように、物理空間Ｓｒにおける速度ｄｘ／ｄｔ、ｄｙ／ｄｔに座標変換係数を乗じることにより、解析空間Ｓａにおける速度の値を計算する。さらに、解析空間Ｓａにおける速度に時間刻み幅ｄｔを乗じることにより、解析空間Ｓａにおける計算粒子の変位量ｄξｉ、ｄηｉを計算する。この変位量に基づいて、解析空間Ｓａにおける計算粒子の位置を更新する。
【0053】
解析の終了条件が満たされるまで、ステップＳ０４からステップＳ０８までの手順を繰り返す。（ステップＳ０９）。
【0054】
次に、解析空間Ｓａにおける計算粒子の位置を、物理空間Ｓｒにおける位置に変換する（ステップＳ１０）。物理空間Ｓｒにおける位置は、式（２１）のベクトル量Ａに物理空間における位置ベクトルを与えることにより算出することができる。
【0055】
解析空間Ｓａにおける位置が物理空間Ｓｒにおける位置に変換されると、解析結果を出力する（ステップＳ１１）。ｉ番目の計算粒子の位置における液膜の厚さＨｉは、式（１５）を用いて計算することができる。
【0056】
次に、上位実施例の優れた効果について説明する。
上記実施例では、二次元の解析空間で解析を行っているため、三次元の物理空間で解析を行う場合と比べて、計算負荷を低減させることができる。
上記実施例では、二次元の解析空間で解析を行っているため、三次元の物理空間で解析を行う場合と比べて、計算負荷を低減させることができる。
【0057】
上記実施例では、ＣＡＤデータに格納されている情報（図１）のみで解析空間Ｓａ（図２）を生成し、座標変換係数（式（１３）、式（１４））を算出することができる。このため、本シミュレーション方法は汎用性が高く、シミュレーションに必要な工数も少ない。
【0058】
ＳＰＨ法のカーネル関数の導関数の公式（式（８）～式（１１））を用いて座標変換係数（式（１３）、式（１４））を算出するため、テクスチャ座標系（一般座標系）の点群（ＣＡＤデータの要素点）は、等間隔に整列している必要がなく、まばらであってもよい。
【0059】
格子を用いた解析では、濡れている領域と濡れていない領域との境界での物理量の不連続を考慮して解析を行う必要があるが、上記実施例ではＳＰＨ法を用いているため、濡れている領域と濡れていない領域とを区別する必要がない。
【0060】
【0061】
図５Ａは、解析対象の模式図である。液膜が流れる球面１０の直径を２ｍとした。この球面の中心Ｏを原点とするｘｙｚ直交座標系を定義する。重力ｇの向きをｙ軸の負の向きとする。重力加速度を９．８１ｍ／ｓ２とした。球面１０とｙ軸の正の部分とが交わる箇所に小孔１１を設置し、小孔１１から球面１０の表面に流量７３ｍｌ／ｓで水を供給した。球面１０における水の接触角を９０°とした。解析空間における計算粒子のカーネル幅を０．０４２（無次元）とした。球面１０の形状は汎用的なＣＡＤソフトウェアを用いて作成し、Ｗａｖｅｆｒｏｎｔ形式に変換して出力したものを用いて解析空間を生成した。
【0062】
図５Ｂは、解析結果を示すための球面１０、小孔１１、及び計算粒子１２の分布を示す模式図である。図５Ｂの左端、中央、及び右端の図は、それぞれ計算開始時点、計算開始時点からの経過時間ｔが０．９９秒、２．０１秒の状態を示す。小孔１１の近傍の淡い色の領域は、液膜の厚さが相対的に厚いことを表している。小孔１１から供給された水が球面１０に沿って広がり、下方に流れていく様子が再現された。上記実施例による方法を用いて、液膜の流れのシミュレーションが可能であることが確認された。
【0063】
【0064】
図６は、本実施例によるシミュレーション方法の手順を示すフローチャートである。
まず、物理空間Ｓｒ上で曲面Ｆｃ（図２）の形状を定義する（ステップＳＡ１）。これは、例えば、図１に示したＣＡＤデータの物理空間の頂点の一覧及び物理空間の頂点法線の一覧を作成することに相当する。次に、曲面上の位置を解析空間上の位置に対応付ける（ステップＳＡ２）。これは、図１に示したＣＡＤデータの面要素の一覧を作成することに相当する。
まず、物理空間Ｓｒ上で曲面Ｆｃ（図２）の形状を定義する（ステップＳＡ１）。これは、例えば、図１に示したＣＡＤデータの物理空間の頂点の一覧及び物理空間の頂点法線の一覧を作成することに相当する。次に、曲面上の位置を解析空間上の位置に対応付ける（ステップＳＡ２）。これは、図１に示したＣＡＤデータの面要素の一覧を作成することに相当する。
【0065】
次に、物理空間Ｓｒに定義される物理座標系と解析空間Ｓａに定義される一般座標系との間で物理量を変換する一般座標変換を、離散化された浅水方程式に施す（ステップＳＡ３）。これは、物理空間における浅水方程式（式（１２））から、式（１５）～式（１８）を導出する手順に相当する。例えば、式（１２）の物理座標系におけるカーネル関数の導関数が、式（１６）、式（１７）において一般座標系におけるカーネル関数の導関数に変換されるとともに、座標変換係数が付加されている。次に、一般座標変換された浅水方程式を数値的に解く（ステップＳＡ４）。これは、式（１５）～式（１８）を実際に計算することに相当する。
【0066】
次に、本実施例の優れた効果について説明する。
一般座標変換された浅水方程式を用いることにより、解析空間上で浅水方程式の値を計算することができる。
一般座標変換された浅水方程式を用いることにより、解析空間上で浅水方程式の値を計算することができる。
【0067】
【0068】
図７は、本実施例によるシミュレーション方法の手順を示すフローチャートである。
まず、物理空間Ｓｒにおける曲面Ｆｃ（図２）の形状を、複数の頂点の座標、及びテクスチャ座標を含むＣＡＤデータで定義する（ステップＳＢ１）。これは、例えば図１に示したＣＡＤデータを作成することに相当する。
まず、物理空間Ｓｒにおける曲面Ｆｃ（図２）の形状を、複数の頂点の座標、及びテクスチャ座標を含むＣＡＤデータで定義する（ステップＳＢ１）。これは、例えば図１に示したＣＡＤデータを作成することに相当する。
【0069】
ＣＡＤデータで定義した頂点の座標（図１の物理空間の頂点の一覧）を、物理空間の座標（Ｘ，Ｙ，Ｚ）として用いる（ステップＳＢ２）。ＣＡＤデータのテクスチャ座標を解析空間Ｓａの一般座標（ξ，η）として用いる（ステップＳＢ３）。
【0070】
次に、支配方程式を数値的に解き、得られた値を物理座標系と一般座標系との間で座標変換する（ステップＳＢ４）。支配方程式は、例えば、浅水方程式（（式１５）～式（１８））に相当する。なお、支配方程式として、離散化されたナビエストークス方程式を用いてもよい。得られた値を座標変換することは、例えば式（２０）を用いて物理座標系における速度（ｕｘ，ｕｙ）を、一般座標系における速度（ξｘｕｘ＋ξｙｕｙ）等に変換し、解析空間における移動距離ｄξｉ等を求めることに相当する。
【0071】
次に、本実施例の優れた効果について説明する。
ＣＡＤデータに基づいて、物理空間上の物理座標と、解析空間上の一般座標とを生成することができる。この物理座標と解析座標とは、種々のシミュレーションで利用することができる。
ＣＡＤデータに基づいて、物理空間上の物理座標と、解析空間上の一般座標とを生成することができる。この物理座標と解析座標とは、種々のシミュレーションで利用することができる。
【0072】
【0073】
図８は、本実施例によるシミュレーション方法の手順を示すフローチャートである。
まず、物理空間Ｓｒにおける曲面Ｆｃ（図２）の形状を、複数の頂点の座標、及びテクスチャ座標を含むＣＡＤデータで定義する（ステップＳＣ１）。この手順は、図７に示したステップＳＢ１の手順と同一である。次に、ＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、物理空間Ｓｒ上の物理座標系と解析空間Ｓａ上の一般座標系との間の座標変換係数を算出する（ステップＳＣ２）。ＳＰＨ法で用いられるカーネル関数の導関数の公式は、例えば式（８）～式（１１）で与えられる。物理空間Ｓｒ上の物理座標系と解析空間Ｓａ上の一般座標系との間の座標変換係数は、例えば、式（１３）及び式（１４）を用いて算出される。
まず、物理空間Ｓｒにおける曲面Ｆｃ（図２）の形状を、複数の頂点の座標、及びテクスチャ座標を含むＣＡＤデータで定義する（ステップＳＣ１）。この手順は、図７に示したステップＳＢ１の手順と同一である。次に、ＣＡＤデータの構造と、ＳＰＨ法で用いられるカーネル関数の導関数の公式とを用いて、物理空間Ｓｒ上の物理座標系と解析空間Ｓａ上の一般座標系との間の座標変換係数を算出する（ステップＳＣ２）。ＳＰＨ法で用いられるカーネル関数の導関数の公式は、例えば式（８）～式（１１）で与えられる。物理空間Ｓｒ上の物理座標系と解析空間Ｓａ上の一般座標系との間の座標変換係数は、例えば、式（１３）及び式（１４）を用いて算出される。
【0074】
計算によって得られた物理量を、ステップＳＣ２で算出された座標変換係数を用いて座標変換する（ステップＳＣ３）。例えば、一般座標系で計算された式（１５）のΣの値に座標変換係数を乗じることにより、物理座標系における値を求める。例えば、式（３０）の物理座標系における速度ｕｘ，ｕｙに座標変換係数を乗じることにより、一般座標系における速度を算出する。
【0075】
次に、本実施例の優れた効果について説明する。
本実施例においては、ＳＰＨ法で用いるカーネル関数の導関数の公式を用いて座標変換係数を算出するため、テクスチャ座標系（一般座標系）における点群が等間隔に整列している必要がなく、まばらであってもよい。
本実施例においては、ＳＰＨ法で用いるカーネル関数の導関数の公式を用いて座標変換係数を算出するため、テクスチャ座標系（一般座標系）における点群が等間隔に整列している必要がなく、まばらであってもよい。
【0076】
上述の各実施例は例示であり、異なる実施例で示した構成の部分的な置換または組み合わせが可能であることは言うまでもない。複数の実施例の同様の構成による同様の作用効果については実施例ごとには逐次言及しない。さらに、本発明は上述の実施例に制限されるものではない。例えば、種々の変更、改良、組み合わせ等が可能なことは当業者に自明であろう。
【符号の説明】
【0077】
１０ 液膜流れが沿う球面
１１ 小孔
１２ 計算粒子
４０ 入力部
４１ 処理部
４２ 出力部
４３ 記憶部
１１ 小孔
１２ 計算粒子
４０ 入力部
４１ 処理部
４２ 出力部
４３ 記憶部
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】